| Class | dynamics_hspl_vas83 |
| In: |
dynamics/dynamics_hspl_vas83.F90
|
Note that Japanese and English are described in parallel.
��å������æ¼�ç®������¢ã�¸ã�¥ã�¼ã���§ã��. æ°´å¹³�¢æ�£å�����¹ã��������� (Bourke, 1988) ��, ���´é�¢æ�£å������ Arakawa and Suarez (1983) �����������¾ã��.
����ç©���æ³��������¼ã�����ã���°ã�¹ã�ã�¼ã�������������¾ã��. �����������§ã��, $ Delta t $ ��大ã�������������, ����æ³¢é���� �»ã���¤ã�³ã�����·ã����æ³��������������¾ã��. NAMELIST#dynamics_hspl_vas83_nml �� TimeIntegScheme ��å¤��´ã����������, ����æ³¢é���������¹ã�����·ã����æ³������£ã��§£����������½ã�§ã��.
This is a dynamical core module. Spectral method (Bourke, 1988) (for horizontal) and Arakawa and Suarez (1983) method (for vertical) are used.
Leap-frog scheme is used as time integration. By default, semi-implicit scheme is applied to gravitational terms in order to enlarge the value of $ Delta t $ . Explicit scheme can be applied to gravitational terms by changing "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".
| DynamicsHSplVAS83 : | ���� |
| DynamicsHSplVAS83Init : | ������ |
| DynamicsHSplVAS83Finalize : | çµ�äº����� (�¢ã�¸ã�¥ã�¼ã����������°ã���²ã��ä»���解é��) |
| ——————— : | ———— |
| DynamicsHSplVAS83 : | Calculate dynamics |
| DynamicsHSplVAS83Init : | Initialization |
| DynamicsHSplVAS83Finalize : | Termination (deallocate variables in this module) |
NAMELIST#dynamics_hspl_vas83_nml
| Subroutine : | |||
| xyz_UB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_VB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_TempB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| xy_PsB(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xyz_UN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_VN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_TempN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| xy_PsN(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xyz_DUDtPhy(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_DVDtPhy(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_DTempDtPhy(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| xy_SurfHeight(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xyz_UA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_VA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_TempA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(out)
| ||
| xy_PsA(0:imax-1, 1:jmax) : | real(DP), intent(out)
| ||
| xyz_ArgOMG(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out), optional
|
��å�������ç®���è¡���, ä¸��������� $ t-\Delta t $ ������ $ t $ �� �±è¥¿é¢��� (xyz_UB, xyz_UN), ����é¢��� (xyz_VB, xyz_VN), æ¸�º¦ (xyz_TempB, xyz_TempN), è³��闋·å��æ¯� (xyzf_QMixB, xyzf_QMixN), �°è¡¨�¢æ��� (xyz_PsB, xyz_PsN), �����³å�°è¡¨�¢é�åº� (xy_SurfHeight) ����, $ t+Delta t $ �� �±è¥¿é¢��� (xyz_UA), ����é¢��� (xyz_VA), æ¸�º¦ (xyz_TempA), è³��闋·å��æ¯� (xyzf_QMixA), �°è¡¨�¢æ��� (xyz_PsA) ��è¿����¾ã��.
�¥ã���������ã�»ã�¹ã������æ¸�º¦, �ºæ��, æ¸�º¦, æ¯�湿ã�������, ��å�����������¶³��������¹ã��������è¨�ç®������´å������, ������������� xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy ��������������.
����ç©���æ³��������¼ã�����ã���°ã�¹ã�ã�¼ã�������������¾ã��. �����������§ã��, $ Delta t $ ��大ã�������������, ����æ³¢é���� �»ã���¤ã�³ã�����·ã����æ³��������������¾ã��. NAMELIST#dynamics_hspl_vas83_nml �� TimeIntegScheme ��å¤��´ã����������, ����æ³¢é���������¹ã�����·ã����æ³������£ã��§£����������½ã�§ã��.
Calculate dynamical processes. Input data are Eastward wind (xyz_UB, xyz_UN), northward wind (xyz_VB, xyz_VN), temperature (xyz_TempB, xyz_TempN), mass mixing ratio (xyzf_QMixB, xyzf_QMixN), surface pressure (xyz_PsB, xyz_PsN) at $ t-\Delta t $ and $ t $, and surface height (xy_SurfHeight). Eastward wind (xyz_UA), northward wind (xyz_VA), temperature (xyz_TempA), mass mixing ratio (xyzf_QMixA), surface pressure (xyz_PsA) are returned.
In order to add tendencies of vorticity, divergence, temperature and specific humidity calculated by physical processes other than those by this dynamical process, these tendencies should be given to "xyz_DUDtPhy", "xyz_DVDtPhy", "xyz_DTempDtPhy", "xyzf_DQMixDtPhy"
Leap-frog scheme is used as time integration. By default, semi-implicit scheme is applied to gravitational terms for extension of $ Delta t $ . Explicit scheme can be applied to gravitational terms by changing "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".
subroutine DynamicsHSplVAS83( xyz_UB, xyz_VB, xyz_TempB, xyzf_QMixB, xy_PsB, xyz_UN, xyz_VN, xyz_TempN, xyzf_QMixN, xy_PsN, xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy, xy_SurfHeight, xyz_UA, xyz_VA, xyz_TempA, xyzf_QMixA, xy_PsA, xyz_ArgOMG )
!
! ��������������, ��������� $ t-\Delta t $ ������ $ t $ ��
! �±è¥¿é¢��� (xyz_UB, xyz_UN), ����é¢��� (xyz_VB, xyz_VN),
! æ¸�º¦ (xyz_TempB, xyz_TempN), è³��闋·å��æ¯� (xyzf_QMixB, xyzf_QMixN),
! �°è¡¨�¢æ��� (xyz_PsB, xyz_PsN), �����³å�°è¡¨�¢é�åº� (xy_SurfHeight) ����,
! $ t+\Delta t $ ��
! �±è¥¿é¢��� (xyz_UA), ����é¢��� (xyz_VA), æ¸�º¦ (xyz_TempA),
! è³��闋·å��æ¯� (xyzf_QMixA), �°è¡¨�¢æ��� (xyz_PsA) ��è¿����¾ã��.
!
! �¥ã���������ã�»ã�¹ã������æ¸�º¦, �ºæ��, æ¸�º¦, æ¯�湿ã�������,
! ��å�����������¶³��������¹ã��������è¨�ç®������´å������,
! ������������� xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy
! ��������������.
!
! ����ç©���æ³��������¼ã�����ã���°ã�¹ã�ã�¼ã�������������¾ã��.
! �����������§ã��, $ \Delta t $ ��大ã�������������, ����æ³¢é����
! �»ã���¤ã�³ã�����·ã����æ³��������������¾ã��.
! NAMELIST#dynamics_hspl_vas83_nml �� TimeIntegScheme ��å¤��´ã����������,
! ����æ³¢é���������¹ã�����·ã����æ³������£ã��§£����������½ã�§ã��.
!
! Calculate dynamical processes.
! Input data are
! Eastward wind (xyz_UB, xyz_UN),
! northward wind (xyz_VB, xyz_VN),
! temperature (xyz_TempB, xyz_TempN),
! mass mixing ratio (xyzf_QMixB, xyzf_QMixN),
! surface pressure (xyz_PsB, xyz_PsN) at $ t-\Delta t $ and $ t $,
! and surface height (xy_SurfHeight).
! Eastward wind (xyz_UA), northward wind (xyz_VA),
! temperature (xyz_TempA),
! mass mixing ratio (xyzf_QMixA), surface pressure (xyz_PsA)
! are returned.
!
! In order to add tendencies of vorticity, divergence,
! temperature and specific humidity calculated by
! physical processes other than those by
! this dynamical process,
! these tendencies should be given to
! "xyz_DUDtPhy", "xyz_DVDtPhy", "xyz_DTempDtPhy", "xyzf_DQMixDtPhy"
!
! Leap-frog scheme is used as time integration.
! By default, semi-implicit scheme is applied to gravitational terms
! for extension of $ \Delta t $ .
! Explicit scheme can be applied to gravitational terms by changing
! "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".
!
! �¢ã�¸ã�¥ã�¼ã����� ; USE statements
!
use constants, only: RPlanet, CpDry, Grav
! $ g $ [m s-2].
! ���������.
! Gravitational acceleration
! �¼å��¹è¨å®�
! Grid points settings
!
use gridset, only: lmax, imax, jmax, kmax ! ���´å±¤��.
! Number of vertical level
! çµ������¢ã���������¨å®�
! Settings of array for atmospheric composition
!
use composition, only: ncmax, IndexH2OVap, CompositionInqFlagAdv
! ���»ç���
! Time control
!
use timeset, only: DelTime, TimeN, TimesetClockStart, TimesetClockStop
! 座æ����¼ã�¿è¨å®�
! Axes data settings
!
use axesset, only: r_Sigma, z_Sigma ! $ \sigma $ ������ (�´æ��).
! Full $ \sigma $ level
! SPMODEL ���¤ã������, ���¢ä����馹������¢è����½æ�°å���������解ã��(å¤�層å�å¿�)
! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
!
#ifdef LIB_MPI
#ifdef SJPACK
use wa_mpi_module_sjpack, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya => wa_DivLambda_xva, wa_DivMu_xya => wa_DivMu_xva, xy_GradLambda_w => xv_GradLambda_w, xya_GradLambda_wa => xva_GradLambda_wa, xy_GradMu_w => xv_GradMu_w, xya_GradMu_wa => xva_GradMu_wa, w_xy => w_xv, xy_w => xv_w, wa_xya => wa_xva, xya_wa => xva_wa
#else
use wa_mpi_module, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya => wa_DivLambda_xva, wa_DivMu_xya => wa_DivMu_xva, xy_GradLambda_w => xv_GradLambda_w, xya_GradLambda_wa => xva_GradLambda_wa, xy_GradMu_w => xv_GradMu_w, xya_GradMu_wa => xva_GradMu_wa, w_xy => w_xv, xy_w => xv_w, wa_xya => wa_xva, xya_wa => xva_wa
#endif
#elif AXISYMMETRY
use wa_zonal_module, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya, wa_DivMu_xya, xy_GradLambda_w, xya_GradLambda_wa, xy_GradMu_w, xya_GradMu_wa, w_xy, xy_w, wa_xya, xya_wa
#elif SJPACK
use wa_module_sjpack, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya, wa_DivMu_xya, xy_GradLambda_w, xya_GradLambda_wa, xy_GradMu_w, xya_GradMu_wa, w_xy, xy_w, wa_xya, xya_wa
#elif AXISYMMETRY_SJPACK
use wa_zonal_module_sjpack, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya, wa_DivMu_xya, xy_GradLambda_w, xya_GradLambda_wa, xy_GradMu_w, xya_GradMu_wa, w_xy, xy_w, wa_xya, xya_wa
#else
use wa_module, only: wa_Lapla_wa, w_Lapla_w, wa_LaplaInv_wa, wa_DivLambda_xya, wa_DivMu_xya, xy_GradLambda_w, xya_GradLambda_wa, xy_GradMu_w, xya_GradMu_wa, w_xy, xy_w, wa_xya, xya_wa
#endif
! ���¹ã�������¼ã�¿å�ºå��
! History data output
!
use gtool_historyauto, only: HistoryAutoPut
! �������
! Character handling
!
use dc_string, only: LChar
! ç¨��¥å�������¡ã��
! Kind type parameter
!
use dc_types, only: DP ! ��精度å®��°å��. Double precision.
! ç©�����¹³������ä½�
! Operation for integral and average
!
use intavr_operate, only: IntLonLat_xy
! è³�����æ�
! Mass fixer
!
use mass_fixer, only: MassFixer, MassFixerColumn
! �»ã�����°ã���³ã�¸ã�¥æ���������³ªç§»æ�è¨�ç®�
! Semi-Lagrangian method for tracer transport
use sltt, only: SLTTMain
! 宣�� ; Declaration statements
!
implicit none
real(DP), intent(in):: xyz_UB (0:imax-1, 1:jmax, 1:kmax)
! $ u (t-\Delta t) $ . �±è¥¿é¢���. Eastward wind
real(DP), intent(in):: xyz_VB (0:imax-1, 1:jmax, 1:kmax)
! $ v (t-\Delta t) $ . �������. Northward wind
real(DP), intent(in):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax)
! $ T (t-\Delta t) $ . æ¸�º¦. Temperature
real(DP), intent(in):: xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ q (t-\Delta t) $ . ��. Specific humidity
real(DP), intent(in):: xy_PsB (0:imax-1, 1:jmax)
! $ p_s (t-\Delta t) $ . �°è¡¨�¢æ���. Surface pressure
real(DP), intent(in):: xyz_UN (0:imax-1, 1:jmax, 1:kmax)
! $ u (t) $ . �±è¥¿é¢���. Eastward wind
real(DP), intent(in):: xyz_VN (0:imax-1, 1:jmax, 1:kmax)
! $ v (t) $ . �������. Northward wind
real(DP), intent(in):: xyz_TempN (0:imax-1, 1:jmax, 1:kmax)
! $ T (t) $ . æ¸�º¦. Temperature
real(DP), intent(in):: xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ q (t) $ . ��. Specific humidity
real(DP), intent(in):: xy_PsN (0:imax-1, 1:jmax)
! $ p_s (t) $ . �°è¡¨�¢æ���. Surface pressure
real(DP), intent(in):: xyz_DUDtPhy (0:imax-1, 1:jmax, 1:kmax)
! $ \left(\DP{u}{t}\right)^{phy} $ .
! å¤����� (�������) �������±è¥¿é¢���å¤���.
! Eastward wind tendency by external force terms (physical processes)
real(DP), intent(in):: xyz_DVDtPhy (0:imax-1, 1:jmax, 1:kmax)
! $ \left(\DP{v}{t}\right)^{phy} $ .
! ����� (�������) ����������������.
! Northward wind tendency by external force terms (physical processes)
real(DP), intent(in):: xyz_DTempDtPhy (0:imax-1, 1:jmax, 1:kmax)
! $ \left(\DP{T}{t}\right)^{phy} $ .
! å¤����� (�������) ������æ¸�º¦å¤���.
! Temperature tendency by external force terms (physical processes)
real(DP), intent(in):: xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ \left(\DP{q}{t}\right)^{phy} $ .
! ����� (�������) �������湿��.
! Temperature tendency by external force terms (physical processes)
real(DP), intent(in):: xy_SurfHeight (0:imax-1, 1:jmax)
! $ z_s $ . �°è¡¨�¢é�åº�.
! Surface height.
real(DP), intent(out):: xyz_UA (0:imax-1, 1:jmax, 1:kmax)
! $ u (t+\Delta t) $ . �±è¥¿é¢���. Eastward wind
real(DP), intent(out):: xyz_VA (0:imax-1, 1:jmax, 1:kmax)
! $ v (t+\Delta t) $ . �������. Northward wind
real(DP), intent(out):: xyz_TempA (0:imax-1, 1:jmax, 1:kmax)
! $ T (t+\Delta t) $ . æ¸�º¦. Temperature
real(DP), intent(out):: xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ q (t+\Delta t) $ . ��. Specific humidity
real(DP), intent(out):: xy_PsA (0:imax-1, 1:jmax)
! $ p_s (t+\Delta t) $ . �°è¡¨�¢æ���. Surface pressure
real(DP), intent(out), optional :: xyz_ArgOMG(0:imax-1, 1:jmax, 1:kmax)
!
! OMEGA, DP/Dt
! ä½�æ¥å���
! Work variables
!
real(DP):: xyz_UCosLatN (0:imax-1, 1:jmax, 1:kmax)
! $ U (t) = u (t) \cos \varphi $ .
real(DP):: xyz_VCosLatN (0:imax-1, 1:jmax, 1:kmax)
! $ V (t) = v (t) \cos \varphi $ .
real(DP):: xyz_VorN (0:imax-1, 1:jmax, 1:kmax)
! $ \zeta (t) $ . æ¸�º¦. Vorticity
real(DP):: xyz_DivN (0:imax-1, 1:jmax, 1:kmax)
! $ D (t) $ . �ºæ��. Divergence
real(DP):: xy_PiN (0:imax-1, 1:jmax)
! $ \pi = \ln p_s $
real(DP):: wz_DivN (lmax, 1:kmax)
! $ D (t) $ . �ºæ�� (�¹ã��������).
! Divergence (spectral)
real(DP):: wz_TempN (lmax, 1:kmax)
! $ T (t) $ . æ¸�º¦ (�¹ã��������).
! Temperature (spectral)
real(DP):: w_PiN (lmax)
! $ \pi $ �¹ã��������
real(DP):: xy_GradLambdaPiN (0:imax-1, 1:jmax)
! $ \DP{\pi}{\lambda} $
real(DP):: xy_GradMuPiN (0:imax-1, 1:jmax)
! $ (1-\mu^2) \DP{\pi}{\mu} $
real(DP):: xyz_PiAdv (0:imax-1, 1:jmax, 1:kmax)
! $ \Dvect{v} \cdot \nabla \pi $ .
! $ \pi $ ��§»æµ�. Advection of $ \pi $
real(DP):: xyz_UAdvN (0:imax-1, 1:jmax, 1:kmax)
! $ U_A (t) $ . �±è¥¿�����霡»æ���.
! Eastward advection of momentum
real(DP):: xyz_VAdvN (0:imax-1, 1:jmax, 1:kmax)
! $ V_A (t) $ . ���������霡»æ���.
! Northward advection of momentum
real(DP):: xyz_TempNonLinearN (0:imax-1, 1:jmax, 1:kmax)
! $ H (t) $ . æ¸�º¦����å¤�����.
! Temperature tendency
real(DP):: xyzf_QMixNonLinearN (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ R (t) $ . �湿��������.
! Specific humidity tendency
real(DP):: xyz_KinEngyN (0:imax-1, 1:jmax, 1:kmax)
! $ KE (t) $ . �������������¼é��.
! Kinetic energy
real(DP):: xyz_TempUAdvN (0:imax-1, 1:jmax, 1:kmax)
! $ UT' (t) $ . æ¸�º¦�±è¥¿ç§»æ���.
! Eastward advection of temperature
real(DP):: xyz_TempVAdvN (0:imax-1, 1:jmax, 1:kmax)
! $ VT' (t) $ . æ¸�º¦����移æ���.
! Northward advection of temperature
real(DP):: xyr_SigDotN (0:imax-1, 1:jmax, 0:kmax)
! $ \dot{\sigma} $ .
! ���´æ�. Vertical flow
real(DP):: xy_DPiDtNG (0:imax-1, 1:jmax)
! $ Z $ . �°è¡¨�¢æ��§æ����å¤���������æ³¢é��.
! Non-gravity wave component of surface pressure tendency
real(DP):: xyzf_QMixUAdvN (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ Uq (t) $ . �湿�西移��.
! Eastward advection of specific humidity
real(DP):: xyzf_QMixVAdvN (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ Vq (t) $ . �湿���移��.
! Northward advection of specific humidity
real(DP):: wz_DVorDtNG (lmax, 1:kmax)
! $ \left( \DD{\zeta}{t} (t) \right)^{NG}$ . æ¸�º¦å¤�����������æ³¢æ���� (�¹ã��������).
! Non-gravity wave component of vorticity tendency (spectral)
real(DP):: wz_DDivDtNG (lmax, 1:kmax)
! $ \left( \DD{D}{t} (t) \right)^{NG}$ . �ºæ�£å�����������æ³¢æ���� (�¹ã��������).
! Non-gravity wave component of divergence tendency (spectral)
real(DP):: wz_DTempDtNG (lmax, 1:kmax)
! $ \left( \DD{T}{t} (t) \right)^{NG}$ . æ¸�º¦å¤�����������æ³¢æ���� (�¹ã��������).
! Non-gravity wave component of temperature tendency (spectral)
real(DP):: wzf_DQMixDtN(lmax, 1:kmax, 1:ncmax)
! $ \DD{q}{t} (t) $ . æ¯�湿å��� (�¹ã��������).
! Specific humidity tendency (spectral)
real(DP):: w_DPiDtNG(lmax)
! $ \left( \DD{p_s}{t} (t) \right)^{NG}$ . �°è¡¨�¢æ��§å���è²»é����æ³¢é�� (�¹ã��������).
! Non-gravity wave component of surface pressure tendency (spectral)
real(DP):: xyz_UCosLatB (0:imax-1, 1:jmax, 1:kmax)
! $ U (t-\Delta t) = u (t-\Delta t) \cos \varphi $ .
real(DP):: xyz_VCosLatB (0:imax-1, 1:jmax, 1:kmax)
! $ V (t-\Delta t) = v (t-\Delta t) \cos \varphi $ .
real(DP):: wz_VorB (lmax, 1:kmax)
! $ \zeta (t-\Delta t) $ . æ¸�º¦ (�¹ã��������).
! Vorticity (spectral)
real(DP):: wz_DivB (lmax, 1:kmax)
! $ D (t-\Delta t) $ . �ºæ�� (�¹ã��������).
! Divergence (spectral)
real(DP):: wz_TempB (lmax, 1:kmax)
! $ T (t-\Delta t) $ . æ¸�º¦ (�¹ã��������).
! Temperature (spectral)
real(DP):: w_PiB (lmax)
! $ \pi = \ln p_s (t-\Delta t) $ . �°è¡¨�¢æ��� (�¹ã��������).
! Surface pressure (spectral)
real(DP):: wzf_QMixB(lmax, 1:kmax, 1:ncmax)
! $ q (t-\Delta t) $ . æ¯�æ¹� (�¹ã��������).
! Specific humidity (spectral)
real(DP):: xyz_DUDtPhyCosLat (0:imax-1, 1:jmax, 1:kmax)
! $ \left(\DP{u}{t}\right)^{phy} \cos \varphi $ .
real(DP):: xyz_DVDtPhyCosLat (0:imax-1, 1:jmax, 1:kmax)
! $ \left(\DP{v}{t}\right)^{phy} \cos \varphi $ .
real(DP):: wz_VorA (lmax, 1:kmax)
! $ \zeta (t+\Delta t) $ . æ¸�º¦ (�¹ã��������).
! Vorticity (spectral)
real(DP):: wz_DivA (lmax, 1:kmax)
! $ D (t+\Delta t) $ . �ºæ�� (�¹ã��������).
! Divergence (spectral)
real(DP):: wz_TempA (lmax, 1:kmax)
! $ T (t+\Delta t) $ . æ¸�º¦ (�¹ã��������).
! Temperature (spectral)
real(DP):: w_PiA (lmax)
! $ \pi = \ln p_s (t+\Delta t) $ . �°è¡¨�¢æ��� (�¹ã��������).
! Surface pressure (spectral)
real(DP):: wzf_QMixA(lmax, 1:kmax, 1:ncmax)
! $ q (t+\Delta t) $ . æ¯�æ¹� (�¹ã��������).
! Specific humidity (spectral)
real(DP):: wz_Psi (lmax, 1:kmax)
! $ \psi $ . æµ�ç·��¢æ��. Streamline function
real(DP):: wz_Chi (lmax, 1:kmax)
! $ \chi $ . �����³ã�·ã�£ã��. Potential
real(DP):: xyz_UCosLatA (0:imax-1, 1:jmax, 1:kmax)
! $ U (t+\Delta t) = u (t+\Delta t) \cos \varphi $ .
real(DP):: xyz_VCosLatA (0:imax-1, 1:jmax, 1:kmax)
! $ V (t+\Delta t) = v (t+\Delta t) \cos \varphi $ .
real(DP):: wz_VorDiffA (lmax, 1:kmax)
! $ \mathscr{D}(\zeta) $ .
! �����閶´å¹³�¡æ�£ã������æ¸�º¦å¤��� (�¹ã��������).
! Vorticity tendency by
! horizontal momentum diffusion (spectral)
real(DP):: wz_DivDiffA (lmax, 1:kmax)
! $ \mathscr{D}(D) $ .
! �����閶´å¹³�¡æ�£ã�������ºæ�£å��� (�¹ã��������).
! Divergence tendency by
! horizontal momentum diffusion (spectral)
real(DP):: wz_PsiDiff (lmax, 1:kmax)
! �����閶´å¹³�¡æ�£ã������æµ�ç·��¢æ�� $ \psi $ å¤���
! Streamline function tendency by
! horizontal momentum diffusion
real(DP):: wz_ChiDiff (lmax, 1:kmax)
! �����閶´å¹³�¡æ�£ã�����������³ã�·ã�£ã�� $ \chi $ å¤���
! Potential tendency by
! horizontal momentum diffusion
real(DP):: xyz_UDiff (0:imax-1, 1:jmax, 1:kmax)
! �����閶´å¹³�¡æ�£ã�������±è¥¿é¢����.
! Eastward wind tendency by
! horizontal momentum diffusion
real(DP):: xyz_VDiff (0:imax-1, 1:jmax, 1:kmax)
! �����閶´å¹³�¡æ�£ã����������é¢����.
! Northward wind tendency by
! horizontal momentum diffusion
! �°è¡¨�¢é�åº��¢é��
! Surface height etc.
!
real(DP):: w_SurfGeoPot (lmax)
! $ \Phi_s $ . �°è¡¨�¸ã�������³ã�·ã�£ã��.
! Surface geo-potential
! Variables for mass fixing
!
real(DP):: MassB
real(DP):: MassA
real(DP):: xyr_PressA(0:imax-1, 1:jmax, 0:kmax)
real(DP):: xyr_PressB(0:imax-1, 1:jmax, 0:kmax)
real(DP) :: xyz_PressA (0:imax-1, 1:jmax, 1:kmax)
real(DP) :: xy_HDifPsA (0:imax-1, 1:jmax)
real(DP) :: xyz_HDifPressA (0:imax-1, 1:jmax, 1:kmax)
real(DP) :: xyzf_DQMixDtHDCor(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
integer :: kp, kn
real(DP) :: xyz_OMG (0:imax-1, 1:jmax, 1:kmax)
!
! OMEGA, DP/Dt
real(DP) :: xyzf_QMixTentative(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
logical :: FlagExistNonAdvComp
integer:: k ! ���´æ�¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in vertical direction
integer:: n ! çµ����¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in dimension of constituents
! ���� ; Executable statement
!
! ������確è�
! Initialization check
!
if ( .not. dynamics_hspl_vas83_inited ) then
call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
end if
! �����������
! Start measurement of computation time
!
call TimesetClockStart( module_name )
! TimeIntegration �§ä½¿������ä¿��°ã��¨å®�
! Configure coefficients for "TimeIntegration"
!
call SemiImplMatrix
! �¼å��¹å�¤ã���¹ã���������¤ã�� ( $ t-\Delta t$ )
! Exchange grid values to spectral values ( $ t-\Delta t$ )
!
! æ¸�º¦�ºæ�£ã���ç®�
! Calculate vorticity and divergence
do k = 1, kmax
xyz_UCosLatB(:,:,k) = xyz_UB(:,:,k) * xy_CosLat
xyz_VCosLatB(:,:,k) = xyz_VB(:,:,k) * xy_CosLat
end do
wz_VorB = ( wa_DivLambda_xya( xyz_VCosLatB ) - wa_DivMu_xya ( xyz_UCosLatB ) ) / RPlanet
wz_DivB = ( wa_DivLambda_xya( xyz_UCosLatB ) + wa_DivMu_xya ( xyz_VCosLatB ) ) / RPlanet
!
wz_TempB = wa_xya( xyz_TempB )
!
if (.not. FlagSLTT) then
do n = 1, ncmax
wzf_QMixB(:,:,n) = wa_xya( xyzf_QMixB(:,:,:,n) )
end do
endif
!
w_PiB = w_xy( log( xy_PsB ) )
! �¼å��¹å�¤ã���¹ã���������¤ã�� ( $ t $ )
! Exchange grid values to spectral values ( $ t $ )
!
! æ¸�º¦�ºæ�£ã���ç®�
! Calculate vorticity and divergence
!
do k = 1, kmax
xyz_UCosLatN(:,:,k) = xyz_UN(:,:,k) * xy_CosLat
xyz_VCosLatN(:,:,k) = xyz_VN(:,:,k) * xy_CosLat
end do
xyz_VorN = xya_wa( wa_DivLambda_xya( xyz_VCosLatN ) - wa_DivMu_xya ( xyz_UCosLatN ) ) / RPlanet
wz_DivN = ( wa_DivLambda_xya( xyz_UCosLatN ) + wa_DivMu_xya ( xyz_VCosLatN ) ) / RPlanet
xyz_DivN = xya_wa( wz_DivN )
!
select case ( IDTimeIntegScheme )
case ( IDTimeIntegSchemeExplicit )
wz_TempN = wa_xya( xyz_TempN )
end select
!
xy_PiN = log( xy_PsN )
w_PiN = w_xy( xy_PiN )
xy_GradLambdaPiN = xy_GradLambda_w( w_PiN ) / RPlanet
xy_GradMuPiN = xy_GradMu_w ( w_PiN ) / RPlanet
! �°è¡¨�¸ã�������³ã�·ã�£ã�����ç®�
! Calculate surface geo-potential
!
w_SurfGeoPot = w_xy( Grav * xy_SurfHeight )
! �¼å��¹ä��§ã����ç·�å½¢å��å������ç®�
! Calculate non-linear dynamical terms on grid points
!
call NonLinearOnGrid( xyz_UCosLatN, xyz_VCosLatN, xyz_VorN, xyz_DivN, xyz_TempN, xyzf_QMixN, xy_GradLambdaPiN, xy_GradMuPiN, xyz_PiAdv, xyz_UAdvN, xyz_VAdvN, xyz_TempNonLinearN, xyzf_QMixNonLinearN, xyz_KinEngyN, xyz_TempUAdvN, xyz_TempVAdvN, xyr_SigDotN, xy_DPiDtNG, xyzf_QMixUAdvN, xyzf_QMixVAdvN )
! ��ç·�å½¢é���������������������å¤������������¹ã�����������
! Spectral transformation of sum of non-linear term and tendencies by physical
! processes
!
do k = 1, kmax
xyz_DUDtPhyCosLat(:,:,k) = xyz_DUDtPhy(:,:,k) * xy_CosLat
xyz_DVDtPhyCosLat(:,:,k) = xyz_DVDtPhy(:,:,k) * xy_CosLat
end do
wz_DVorDtNG = ( wa_DivLambda_xya( xyz_VAdvN + xyz_DVDtPhyCosLat ) - wa_DivMu_xya ( xyz_UAdvN + xyz_DUDtPhyCosLat ) ) / RPlanet
wz_DDivDtNG = ( wa_DivLambda_xya( xyz_UAdvN + xyz_DUDtPhyCosLat ) + wa_DivMu_xya ( xyz_VAdvN + xyz_DVDtPhyCosLat ) ) / RPlanet - wa_Lapla_wa( wa_xya(xyz_KinEngyN) ) / RPlanet**2
wz_DTempDtNG = - ( wa_DivLambda_xya( xyz_TempUAdvN ) + wa_DivMu_xya ( xyz_TempVAdvN ) ) / RPlanet + wa_xya( xyz_TempNonLinearN + xyz_DTempDtPhy )
w_DPiDtNG = w_xy( xy_DPiDtNG )
if (.not. FlagSLTT) then
do n = 1, ncmax
wzf_DQMixDtN(:,:,n) = - ( wa_DivLambda_xya( xyzf_QMixUAdvN(:,:,:,n) ) + wa_DivMu_xya ( xyzf_QMixVAdvN(:,:,:,n) ) ) / RPlanet + wa_xya( xyzf_QMixNonLinearN(:,:,:,n) + xyzf_DQMixDtPhy (:,:,:,n) )
end do
endif
! �������
! Time integration
!
call TimeIntegration( w_SurfGeoPot, wz_DVorDtNG, wz_DDivDtNG, wz_DTempDtNG, wzf_DQMixDtN, w_DPiDtNG, wz_DivN, wz_TempN, w_PiN, wz_VorB, wz_DivB, wz_TempB, wzf_QMixB, w_PiB, wz_VorA, wz_DivA, wz_TempA, wzf_QMixA, w_PiA )
! Divergence damping
!
call DivergenceDamping( wz_DivA )
! �¹ã���������¤ã���¼å��¹å�¤ã�� ( $ t+\Delta t$ )
! Exchange spectral values to grid values ( $ t+\Delta t$ )
!
wz_Psi = wa_LaplaInv_wa( wz_VorA ) * RPlanet**2
wz_Chi = wa_LaplaInv_wa( wz_DivA ) * RPlanet**2
xyz_UCosLatA = ( xya_GradLambda_wa( wz_Chi ) - xya_GradMu_wa ( wz_Psi ) ) / RPlanet
xyz_VCosLatA = ( xya_GradLambda_wa( wz_Psi ) + xya_GradMu_wa ( wz_Chi ) ) / RPlanet
do k = 1, kmax
xyz_UA(:,:,k) = xyz_UCosLatA(:,:,k) / xy_CosLat
xyz_VA(:,:,k) = xyz_VCosLatA(:,:,k) / xy_CosLat
end do
xyz_TempA = xya_wa( wz_TempA )
if (.not. FlagSLTT) then
do n = 1, ncmax
xyzf_QMixA(:,:,:,n) = xya_wa( wzf_QMixA(:,:,n) )
end do
endif
xy_PsA = exp( xy_w( w_PiA ) )
! æ°´å¹³�¡æ�£ã���¹ã���³ã�¸å±¤�������£é�¸ã���æ�
! Correction for dissipation of kinetic enery by horizontal diffusion and sponge
! layer
!
! æ°´å¹³�¡æ�£ã���¹ã���³ã�¸å±¤������æ¸�º¦�ºæ�£ã������å¤���
! Vorticity and divergence tendency by horizontal diffusion and damping in a sponge
! layer
!
wz_VorDiffA = wz_VorA * wz_DisCoefM
wz_DivDiffA = wz_DivA * wz_DisCoefM
! æ°´å¹³�¡æ�£ã���¹ã���³ã�¸å±¤�������£é�¸ã�������±è�æ�
! Correction for internal energy by horizontal diffusion and damping in a sponge
! layer
!
wz_PsiDiff = wa_LaplaInv_wa( wz_VorDiffA ) * RPlanet**2
wz_ChiDiff = wa_LaplaInv_wa( wz_DivDiffA ) * RPlanet**2
xyz_UDiff = ( xya_GradLambda_wa( wz_ChiDiff ) - xya_GradMu_wa ( wz_PsiDiff ) ) / RPlanet
xyz_VDiff = ( xya_GradLambda_wa( wz_PsiDiff ) + xya_GradMu_wa ( wz_ChiDiff ) ) / RPlanet
do k = 1, kmax
xyz_TempA(:,:,k) = xyz_TempA(:,:,k) - ( xyz_UA(:,:,k) * xyz_UDiff(:,:,k) + xyz_VA(:,:,k) * xyz_VDiff(:,:,k) ) / xy_CosLat / CpDry * 2.0_DP * DelTime
end do
if ( FlagMassHorDifCor ) then
! Correction for horizontal diffusion of constituents.
! This is performed for the aim of correcting the diffusion in the viscinity of
! steep terrain slope.
!
do k = 1, kmax
xyz_PressA(:,:,k) = xy_PsA * z_Sigma(k)
end do
xy_HDifPsA = xy_w( w_DisCoefQ * w_xy( xy_PsA ) )
do k = 1, kmax
xyz_HDifPressA(:,:,k) = z_Sigma(k) * xy_HDifPsA(:,:)
end do
if (.not. FlagSLTT) then
do n = 1, ncmax
do k = 1, kmax
kp = k-1
kn = k+1
if( k == 1 ) kp = 1
if( k == kmax ) kn = kmax
xyzf_DQMixDtHDCor(:,:,k,n) = - ( xyzf_QMixA(:,:,kn,n) - xyzf_QMixA(:,:,kp,n) ) / ( xyz_PressA(:,:,kn) - xyz_PressA(:,:,kp) ) * xyz_HDifPressA(:,:,k)
end do
end do
xyzf_QMixA = xyzf_QMixA + xyzf_DQMixDtHDCor * 2.0_DP * DelTime
endif
end if
! Mass fixer
! Total Mass
! NOTICE: It is assumed that the total atmospheric mass does not change.
!
if ( FlagTotMassFix ) then
MassB = IntLonLat_xy( xy_PsB )
MassA = IntLonLat_xy( xy_PsA )
if ( MassA /= 0.0_DP ) then
xy_PsA = MassB / MassA * xy_PsA
end if
end if
! 主ã��§»æµ����¹ã����������½¿�� if ��
if ( .not. FlagCalcUVTPs ) then
xyz_UA = xyz_UB
xyz_VA = xyz_VB
xyz_TempA = xyz_TempB
xy_PsA = xy_PsB
xyr_SigDotN = 0.0_DP
end if
if (FlagSLTT) then
do k = 0, kmax
xyr_PressA(:,:,k) = xy_PsA * r_Sigma(k)
xyr_PressB(:,:,k) = xy_PsB * r_Sigma(k)
end do
! Mass fixer is included in SLTTMain
call SLTTMain( xyr_PressB, xyr_PressA, xyz_UN, xyz_VN, xyr_SigDotN, xyzf_DQMixDtPhy, xyzf_QMixB, xyzf_QMixA )
else
! Mass fixer
! Constituents
!
do k = 0, kmax
xyr_PressA(:,:,k) = xy_PsA * r_Sigma(k)
xyr_PressB(:,:,k) = xy_PsB * r_Sigma(k)
end do
call MassFixer( xyr_PressA, xyzf_QMixA, xyr_PressRef = xyr_PressB, xyzf_QMixRef = ( xyzf_QMixB+xyzf_DQMixDtPhy*2.0_DP*DelTime ) )
endif
! Calculation for non-advection case
! ���¤ã���¼ç§»æµ�è¨�ç®������������¼ã�ã�¼ã����å¸�����³ª��移æ�����������è¨�ç®�å®�å®���
! ������
!
FlagExistNonAdvComp = .false.
do n = 1, ncmax
if ( .not. CompositionInqFlagAdv( n ) ) then
xyzf_QMixA(:,:,:,n) = xyzf_QMixB(:,:,:,n) + xyzf_DQMixDtPhy(:,:,:,n) * 2.0_DP * DelTime
xyzf_QMixTentative(:,:,:,n) = xyzf_QMixA(:,:,:,n)
FlagExistNonAdvComp = .true.
else
xyzf_QMixTentative(:,:,:,n) = 0.0_DP
end if
end do
if ( FlagExistNonAdvComp ) then
! Mass fixer
! Constituents
call MassFixerColumn( xyr_PressB, xyzf_QMixTentative, xyr_PressRef = xyr_PressB, xyzf_QMixRef = ( xyzf_QMixB+xyzf_DQMixDtPhy*2.0_DP*DelTime ) )
do n = 1, ncmax
if ( .not. CompositionInqFlagAdv( n ) ) then
xyzf_QMixA(:,:,:,n) = xyzf_QMixTentative(:,:,:,n)
end if
end do
end if
! ���¹ã�������¼ã�¿å�ºå��
! History data output
!
call OutputDiagnosedVariables( xyz_UB, xyz_VB, xyz_TempB, xyzf_QMixB, xy_PsB, xyz_UN, xyz_VN, xyz_TempN, xyzf_QMixN, xy_PsN, xyz_UA, xyz_VA, xyz_TempA, xyzf_QMixA, xy_PsA, xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy, xyr_SigDotN, xy_DPiDtNG, xyz_PiAdv, xyz_OMG )
if ( present( xyz_ArgOMG ) ) then
xyz_ArgOMG = xyz_OMG
end if
! è¨�ç®�����è¨�æ¸������æ�
! Pause measurement of computation time
!
call TimesetClockStop( module_name )
end subroutine DynamicsHSplVAS83
| Subroutine : |
�¢ã�¸ã�¥ã�¼ã����������°ã���²ã��ä»���解é�¤ã��è¡����¾ã��.
Deallocate variables in this module.
subroutine DynamicsHSplVAS83Finalize
!
! �¢ã�¸ã�¥ã�¼ã����������°ã���²ã��ä»���解é�¤ã��è¡����¾ã��.
!
! Deallocate variables in this module.
!
! 宣�� ; Declaration statements
!
implicit none
! ���� ; Executable statement
!
if ( .not. dynamics_hspl_vas83_inited ) return
! �����������¤ã�¸æ�»ã��
! Return to default values
!
DelTimeSave = - 1.0_DP
! �²ã��ä»���解é��
! Deallocation
!
if ( allocated( xy_SinLat ) ) deallocate( xy_SinLat )
if ( allocated( xy_CosLat ) ) deallocate( xy_CosLat )
if ( allocated( xy_Cori ) ) deallocate( xy_Cori )
if ( allocated( z_HydroAlpha ) ) deallocate( z_HydroAlpha )
if ( allocated( z_HydroBeta ) ) deallocate( z_HydroBeta )
if ( allocated( z_TInpCoefA ) ) deallocate( z_TInpCoefA )
if ( allocated( z_TInpCoefB ) ) deallocate( z_TInpCoefB )
if ( allocated( z_TInpCoefK ) ) deallocate( z_TInpCoefK )
if ( allocated( z_RefTemp ) ) deallocate( z_RefTemp )
if ( allocated( r_RefTemp ) ) deallocate( r_RefTemp )
if ( allocated( w_LaplaEigVal ) ) deallocate( w_LaplaEigVal )
if ( allocated( wz_DisCoefM ) ) deallocate( wz_DisCoefM )
if ( allocated( wz_DisCoefH ) ) deallocate( wz_DisCoefH )
if ( allocated( w_DisCoefQ ) ) deallocate( w_DisCoefQ )
if ( allocated( z_siMtxC ) ) deallocate( z_siMtxC )
if ( allocated( z_siMtxG ) ) deallocate( z_siMtxG )
if ( allocated( zz_siMtxH ) ) deallocate( zz_siMtxH )
if ( allocated( zz_siMtxDiH ) ) deallocate( zz_siMtxDiH )
if ( allocated( wzz_siMtxWDiH ) ) deallocate( wzz_siMtxWDiH )
if ( allocated( zz_siMtxGCt ) ) deallocate( zz_siMtxGCt )
if ( allocated( zz_siMtxW ) ) deallocate( zz_siMtxW )
if ( allocated( zz_siMtxQ ) ) deallocate( zz_siMtxQ )
if ( allocated( zz_siMtxS ) ) deallocate( zz_siMtxS )
if ( allocated( zz_siMtxR ) ) deallocate( zz_siMtxR )
!!$ if ( allocated(nmo ) ) deallocate( nmo )
!!$ if ( allocated(wzz_siMtxM ) ) deallocate( wzz_siMtxM )
!!$ if ( allocated(z_siMtxPivWork) ) deallocate( z_siMtxPivWork )
if ( allocated(wzz_siMtxLU ) ) deallocate( wzz_siMtxLU )
if ( allocated(wz_siMtxPiv ) ) deallocate( wz_siMtxPiv )
dynamics_hspl_vas83_inited = .false.
end subroutine DynamicsHSplVAS83Finalize
| Subroutine : |
è¨�ç®����è¦��������¡ã�¿ã��¨å®��� NAMELIST#dynamics_hspl_vas83_nml ����¿è¾¼�¿ã��è¡����¾ã��.
Configure parameters for calculation, and load "NAMELIST#dynamics_hspl_vas83_nml"
This procedure input/output NAMELIST#dynamics_hspl_vas83_nml .
subroutine DynamicsHSplVAS83Init
!
! è¨�ç®����è¦��������¡ã�¿ã��¨å®��� NAMELIST#dynamics_hspl_vas83_nml
! ����¿è¾¼�¿ã��è¡����¾ã��.
!
! Configure parameters for calculation,
! and load "NAMELIST#dynamics_hspl_vas83_nml"
!
! �¢ã�¸ã�¥ã�¼ã����� ; USE statements
!
! ����å®��°è¨å®�
! Physical constants settings
!
use constants, only: RPlanet, Omega, GasRDry, CpDry
! $ C_p $ [J kg-1 K-1].
! ä¹¾ç�¥å¤§æ°�����§æ���.
! Specific heat of air at constant pressure
! �¼å��¹è¨å®�
! Grid points settings
!
use gridset, only: nmax, lmax, imax, jmax, kmax ! ���´å±¤��.
! Number of vertical level
! 座æ����¼ã�¿è¨å®�
! Axes data settings
!
use axesset, only: z_Sigma, r_Sigma, z_DelSigma ! $ \Delta \sigma $ (�´æ��).
! $ \Delta \sigma $ (Full)
! SPMODEL ���¤ã������, ���¢ä����馹������¢è����½æ�°å���������解ã��(å¤�層å�å¿�)
! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
!
#ifdef LIB_MPI
#ifdef SJPACK
use wa_mpi_module_sjpack, only: xy_Lat => xv_Lat, w_xy => w_xv, rn, nm_l
#else
use wa_mpi_module, only: xy_Lat => xv_Lat, w_xy => w_xv, rn, nm_l
#endif
#elif AXISYMMETRY
use wa_zonal_module, only: xy_Lat, w_xy, rn, nm_l
#elif SJPACK
use wa_module_sjpack, only: xy_Lat, w_xy, rn, nm_l
#elif AXISYMMETRY_SJPACK
use wa_zonal_module_sjpack, only: xy_Lat, w_xy, rn, nm_l
#else
use wa_module, only: xy_Lat, w_xy, rn, nm_l
#endif
! NAMELIST ���¡ã�¤ã���¥å�����¢ã�������¼ã���£ã������
! Utilities for NAMELIST file input
!
use namelist_util, only: namelist_filename, NmlutilMsg
! ���¹ã�������¼ã�¿å�ºå��
! History data output
!
use gtool_historyauto, only: HistoryAutoAddVariable
! gtool4 ���¼ã�¿å�¥å��
! Gtool4 data input
!
use gtool_history, only: HistoryGet
! ���¡ã�¤ã���¥å�ºå��è£���
! File I/O support
!
use dc_iounit, only: FileOpen
! ç¨��¥å�������¡ã��
! Kind type parameter
!
use dc_types, only: STDOUT, STRING ! ����. Strings.
! �¥ä������³æ���»ã������±ã��
! Date and time handler
!
use dc_calendar, only: DCCalConvertByUnit
! çµ��¿è¾¼�¿é�¢æ�� PRESENT ���¡å¼µ���¢æ��
! Extended functions of intrinsic function "PRESENT"
!
use dc_present, only: present_and_not_empty
! �������
! Character handling
!
use dc_string, only: LChar
! è³�����æ�
! Mass fixer
!
use mass_fixer, only : MassFixerInit
! �»ã�����°ã���³ã�¸ã�¥æ���������³ªç§»æ�è¨�ç®�
! Semi-Lagrangian method for tracer transport
use sltt, only: SLTTInit
! 宣�� ; Declaration statements
!
implicit none
! �ºæ�æ¸�º¦��¨å®����������æ¥å���
! Work variable for reference temperature
!
character(TOKEN) :: TimeIntegScheme
! ��������.
! 以ä����¹æ����¸æ������.
!
! Time integration scheme.
! Available schemes are as follows.
!
! * "Semi-implicit"
! * "Explicit"
real(DP):: RefTemp
! $ \overline{T} $ . �ºæ�æ¸�º¦.
! Reference temperature
! æ°´å¹³�¡æ��, �¹ã���³ã�¸å±¤���������æ¥å���
! Work variable for horizontal diffusion and sponge layer
!
real(DP) :: w_HDifCoefM (lmax)
! $ {\cal D}_M = -K_{HD} [(-1)^{N_D/2}\nabla^{N_D}- (2/a^2)^{N_D/2}] $ .
! �����閶´å¹³�¡æ�£ä���.
! Coefficient of horizontal momentum diffusion
real(DP) :: w_HDifCoefH (lmax)
! $ {\cal D}_H = -(-1)^{N_D/2}K_{HD}\nabla^{N_D} $ .
! ��, 水水平����.
! Coefficient of horizontal thermal and water diffusion
real(DP) :: wz_SpongeLayerCoefM(lmax, 1:kmax)
! �¹ã���³ã�¸å±¤����������������è¡°ä���
! Damping coefficient for momentum in a sponge layer
real(DP) :: wz_SpongeLayerCoefH(lmax, 1:kmax)
! �¹ã���³ã�¸å±¤���������±ã���è¡°ä���
! Damping coefficient for temperature zonal perturbation
! in a sponge layer
integer:: HDOrder ! è¶�ç²��§ã�����. Order of hyper-viscosity
real(DP):: VisCoef ! è¶�ç²��§ä���. Hyper-viscosity coefficient
real(DP):: HDEFoldTimeValue ! ��大波�°ã������� e-folding time.
! è²����¤ã��ä¸�������, æ°´å¹³�¡æ�£ä��°ã���¼ã�ã�����¾ã��.
!
! E-folding time for maximum wavenumber.
! If negative value is given,
! coefficients of horizontal diffusion become zero.
character(TOKEN):: HDEFoldTimeUnit
! ��大波�°ã������� e-folding time ����ä½�.
! Unit of e-folding time for maximum wavenumber
real(DP):: HDEFoldTime ! ��大波�°ã������� e-folding time [��ä½�: ç§�].
! E-folding time for maximum wavenumber [Unit: sec]
logical :: FlagSpongeLayer
logical :: FlagSpongeLayerforZonalMean
logical :: FlagSpongeLayerforHeat
real(DP) :: SLEFoldTimeValue
! �¹ã���³ã�¸å±¤����ä¸�層ã��������æ¸�è¡°æ��å®���
! Damping time constant at the top model layer
! in a sponge layer
character(TOKEN) :: SLEFoldTimeUnit
! SLEFoldTimeValue �����
! Unit of SLEFoldTimeValue
real(DP) :: SLEFoldTime
! �¹ã���³ã�¸å±¤����ä¸�層ã��������æ¸�è¡°æ��å®��� (ç§�)
! Damping time constant at the top model layer
! in a sponge layer (sec)
integer :: SLOrder
! �¹ã���³ã�¸å±¤���è¡°ä��°ã�� sigma ä¾�å��§ã�����¼ã��
! Sigma dependence of damping coefficient
! in a sponge layer
integer :: SLNumLayer
! �¹ã���³ã�¸å±¤������������¢ã�����ç«�������±¤����
! Number of layers which the sponge layer is applied.
integer :: a_DegOrd(2)
real(DP) :: DivDampPeriodValue
! Period for divergence damping application
character(TOKEN) :: DivDampPeriodUnit
! Unit of DivDampPeriodValue
! NonLinearOnGrid ç��§ä½¿������ä¿��°ã��¨å®����������æ¥å���
! Work variable for coefficients for "NonLinearOnGrid", etc.
!
real(DP):: Kappa ! $ \kappa = R/C_p $ .
! æ°�ä½�å®��°ã����§æ��±ã�������æ¯�.
! Ratio of gas constant to specific heat
integer:: k ! ���´æ�¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in vertical direction
integer:: l ! æ³¢æ�°æ�¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in wavenumber
integer:: unit_nml ! NAMELIST ���¡ã�¤ã�����¼ã���³ç���ç½����.
! Unit number for NAMELIST file open
integer:: iostat_nml ! NAMELIST èªã�¿è¾¼�¿æ���� IOSTAT.
! IOSTAT of NAMELIST read
! NAMELIST å¤��°ç¾¤
! NAMELIST group name
!
namelist /dynamics_hspl_vas83_nml/ TimeIntegScheme, HDOrder, HDEFoldTimeValue, HDEFoldTimeUnit, FlagSpongeLayer, FlagSpongeLayerforZonalMean, FlagSpongeLayerforHeat, SLEFoldTimeValue, SLEFoldTimeUnit, SLOrder, SLNumLayer, RefTemp, FlagDivDamp, DivDampPeriodValue, DivDampPeriodUnit, FlagTotMassFix, FlagMassHorDifCor, FlagSLTT, FlagCalcUVTPs
!
! �����������¤ã���¤ã��������������ç¶� "dynamics_hspl_vas83#DynamicsInit"
! ���½ã�¼ã�¹ã�³ã�¼ã�������§ã������.
!
! Refer to source codes in the initialization procedure
! "dynamics_hspl_vas83#DynamicsInit" for the default values.
!
! ���� ; Executable statement
!
if ( dynamics_hspl_vas83_inited ) return
! �����������¤ã��¨å®�
! Default values settings
!
TimeIntegScheme = 'Semi-implicit'
HDOrder = 8
HDEFoldTimeValue = 8640.0_DP
HDEFoldTimeUnit = 'sec'
FlagSpongeLayer = .false.
FlagSpongeLayerforZonalMean = .false.
FlagSpongeLayerforHeat = .false.
SLEFoldTimeValue = 86400.0_DP
SLEFoldTimeUnit = 'sec'
SLOrder = 1
SLNumLayer = kmax
RefTemp = 300.
FlagDivDamp = .false.
DivDampPeriodValue = 2.0_DP
DivDampPeriodUnit = 'day'
FlagTotMassFix = .true.
FlagMassHorDifCor = .false.
FlagSLTT = .false.
FlagCalcUVTPs = .true.
! NAMELIST ����¿è¾¼��
! NAMELIST is input
!
if ( trim(namelist_filename) /= '' ) then
call FileOpen( unit_nml, namelist_filename, mode = 'r' ) ! (in)
rewind( unit_nml )
read( unit_nml, nml = dynamics_hspl_vas83_nml, iostat = iostat_nml ) ! (out)
close( unit_nml )
call NmlutilMsg( iostat_nml, module_name ) ! (in)
if ( iostat_nml == 0 ) write( STDOUT, nml = dynamics_hspl_vas83_nml )
end if
! ����ç©���æ³������§ã����
! Check time integration scheme
!
select case ( LChar( trim( TimeIntegScheme ) ) )
case ('semi-implicit')
IDTimeIntegScheme = IDTimeIntegSchemeSemiImplicit
case ('explicit')
IDTimeIntegScheme = IDTimeIntegSchemeExplicit
case default
call MessageNotify( 'E', module_name, '"TimeIntegScheme" must be "Semi-Implicit" or "Explicit".' )
end select
! SemiImplMatrix �µã�����¼ã���³ç���� $ \Delta t $ ���å��¤ã�������¤ã��è¨å�.
! Configure initial value to saved value of $ \Delta t $ for a subroutine "SemiImplMatrix"
!
DelTimeSave = - 1.0_DP
! NonLinearOnGrid ç��§ä½¿������ä¿��°ã��¨å®�
! Configure coefficients for "NonLinearOnGrid", etc.
!
! $ \sin \varphi $ �� $ \cos \varphi $ ����
! Calculate $ \sin \varphi $ and $ \cos \varphi $
!
allocate( xy_SinLat (0:imax-1, 1:jmax) )
allocate( xy_CosLat (0:imax-1, 1:jmax) )
xy_SinLat = sin( xy_Lat )
xy_CosLat = cos( xy_Lat )
! �³ã�����������¡ã�¼ã�¿ã���ç®�
! Calculate Coriolis parameter
!
allocate( xy_Cori (0:imax-1, 1:jmax) )
xy_Cori = 2.0_DP * Omega * xy_SinLat
! ��水�������� $ \alpha $ , $ \beta $ ����
! Calculate coefficients $ \alpha $ and $ \beta $ in hydrostatic equation.
!
allocate( z_HydroAlpha(1:kmax) )
allocate( z_HydroBeta (1:kmax) )
Kappa = GasRDry / CpDry
do k = 1, kmax
z_HydroAlpha(k) = ( r_Sigma(k-1) / z_Sigma(k) ) ** Kappa - 1.0_DP
z_HydroBeta(k) = 1.0_DP - ( r_Sigma(k) / z_Sigma(k) ) ** Kappa
enddo
! æ¸�º¦���´è�������� $ \kappa $, $ a $ , $ b $ ���ç®�
! Calculate coefficients $ \kappa $, $ a $ , $ b $
! for interpolation of temperature
!
allocate( z_TInpCoefA (1:kmax) )
allocate( z_TInpCoefB (1:kmax) )
allocate( z_TInpCoefK (1:kmax) )
do k = 1, kmax
z_TInpCoefK(k) = ( r_Sigma(k-1) * z_HydroAlpha(k) + r_Sigma(k ) * z_HydroBeta (k) ) / z_DelSigma(k)
enddo
z_TInpCoefA = 0.
do k = 2, kmax
z_TInpCoefA(k) = z_HydroAlpha(k) * ( 1.0_DP - (z_Sigma(k) / z_Sigma(k-1)) ** Kappa ) ** ( -1 )
end do
z_TInpCoefB = 0.
do k = 1, kmax - 1
z_TInpCoefB(k) = z_HydroBeta(k) * ( (z_Sigma(k) / z_Sigma(k+1) ) ** Kappa - 1.0_DP ) ** ( -1 )
end do
! �ºæ�æ¸�º¦ (���´æ�°ã������) ���ç®�
! Calculate reference temperature on half levels
!
allocate( z_RefTemp(1:kmax) )
allocate( r_RefTemp(0:kmax) )
z_RefTemp = RefTemp
r_RefTemp(0) = 0.
r_RefTemp(kmax) = 0.
do k = 1, kmax - 1
r_RefTemp(k) = z_TInpCoefA(k+1) * z_RefTemp(k+1) + z_TInpCoefB( k ) * z_RefTemp( k )
enddo
! ���¸ã�£ã�³ã�������¢æ�°ã���ºæ���¤ã��¨å®�
! Set eigen value of Associated Legendre Function
!
allocate( w_LaplaEigVal(lmax) )
w_LaplaEigVal(:) = rn(:,1) / RPlanet**2
! æ°´å¹³�¡æ�£ä��°ã��¨å®�
! Configure coefficient of horizontal diffusion
!
! ç²��§ä��°ã���ç®� (��大波�°ã�� e-folding time �� HDEFoldTime ������������)
! Calculate viscosity coefficient
!
HDEFoldTime = DCCalConvertByUnit( HDEFoldTimeValue, HDEFoldTimeUnit, 'sec' ) ! (in)
if ( HDEFoldTimeValue > 0.0_DP ) then
VisCoef = ( (nmax*(nmax+1)) / RPlanet**2 )**(-HDOrder / 2) / HDEFoldTime
else
VisCoef = 0.0_DP
end if
w_HDifCoefH = - VisCoef * ( ( - w_LaplaEigVal )**( HDOrder / 2 ) )
w_HDifCoefM = w_HDifCoefH - VisCoef * ( - ( 2.0_DP / RPlanet**2 )**( HDOrder / 2 ) )
! �¹ã���³ã�¸å±¤���è¡°ä��°ã��¨å®�
! Set damping coefficient in a sponge layer
!
if ( SLEFoldTimeValue <= 0.0_DP ) then
call MessageNotify( 'E', module_name, 'SLEFoldTimeValue must be greater than zero, SLEFoldTimeValue=<%f>.', d = (/ SLEFoldTimeValue /) )
end if
if ( ( SLNumLayer <= 0 ) .or. ( SLNumLayer > kmax ) ) then
call MessageNotify( 'E', module_name, 'SLNumLayer must be greater than zero and less than/equal to kmax, SLNumLayer=<%d>.', i = (/ SLNumLayer /) )
end if
! �¹ã���³ã�¸å±¤���è¡°ä��°ã���ç®�
! Calculate damping coefficient for momentum in a sponge layer
!
SLEFoldTime = DCCalConvertByUnit( SLEFoldTimeValue, SLEFoldTimeUnit, 'sec' ) ! (in)
if ( FlagSpongeLayer ) then
! Sponge layer is not applied for model layers, k = 1, kmax-SLNumLayers.
!
do k = 1, kmax-SLNumLayer
wz_SpongeLayerCoefM(:,k) = 0.0_DP
end do
do k = kmax-SLNumLayer+1, kmax
wz_SpongeLayerCoefM(:,k) = 1.0_DP / ( SLEFoldTime * ( z_Sigma(k) / z_Sigma(kmax) )**SLOrder )
end do
if ( .not. FlagSpongeLayerforZonalMean ) then
! Sponge layer is not applied for zonal mean component.
! i.e., sponge layer is applied for zonal wave component only.
!
do k = kmax-SLNumLayer+1, kmax
do l = 1, lmax
a_DegOrd = nm_l( l )
if ( a_DegOrd(2) == 0 ) then
wz_SpongeLayerCoefM(l,k) = 0.0_DP
end if
end do
end do
end if
if ( FlagSpongeLayerforHeat ) then
wz_SpongeLayerCoefH = wz_SpongeLayerCoefM
! Sponge layer is not applied for zonal mean component.
! i.e., sponge layer is applied for zonal wave component only.
!
do k = kmax-SLNumLayer+1, kmax
do l = 1, lmax
a_DegOrd = nm_l( l )
if ( a_DegOrd(2) == 0 ) then
wz_SpongeLayerCoefH(l,k) = 0.0_DP
end if
end do
end do
else
wz_SpongeLayerCoefH = 0.0_DP
end if
else
! Sponge layer is not applied.
!
wz_SpongeLayerCoefM = 0.0_DP
wz_SpongeLayerCoefH = 0.0_DP
end if
! �������°´å¹³æ�¡æ�£ä��°ã���¹ã���³ã�¸å±¤æ¸�è¡°ä��°ã����
! Damping coefficients by horizonatl diffusion and a sponge layer
!
allocate( wz_DisCoefM(lmax, 1:kmax) )
allocate( wz_DisCoefH(lmax, 1:kmax) )
allocate( w_DisCoefQ (lmax) )
do k = 1, kmax
wz_DisCoefM(:,k) = w_HDifCoefM - wz_SpongeLayerCoefM(:,k)
wz_DisCoefH(:,k) = w_HDifCoefH - wz_SpongeLayerCoefH(:,k)
end do
w_DisCoefQ = w_HDifCoefH
! Calculation of divergence damping period
!
DivDampPeriod = DCCalConvertByUnit( DivDampPeriodValue, DivDampPeriodUnit, 'sec' )
! ���¹ã�������¼ã�¿å�ºå�����������¸ã����°ç�»é��
! Register of variables for history data output
!
call HistoryAutoAddVariable( 'DUDtDyn', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'dynamical tendency of zonal wind', 'm s-2' )
call HistoryAutoAddVariable( 'DVDtDyn', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'dynamical tendency of meridional wind', 'm s-2' )
call HistoryAutoAddVariable( 'DTempDtDyn', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'dynamical tendency of temperature', 'K s-1' )
call HistoryAutoAddVariable( 'DQVapDtDyn', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'dynamical tendency of water vapor', 's-1' )
call HistoryAutoAddVariable( 'DPsDtDyn', (/ 'lon ', 'lat ', 'time' /), 'dynamical tendency of surface pressure', 'Ps s-1' )
call HistoryAutoAddVariable( 'OMG', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'vertical velocity in pressure coordinate (omega, DP/Dt)', 'Pa s-1' )
call HistoryAutoAddVariable( 'SigDot', (/ 'lon ', 'lat ', 'sigm', 'time' /), 'sigma-vertical velocity', '1 s-1' )
call HistoryAutoAddVariable( 'DPiDt', (/ 'lon ', 'lat ', 'time' /), 'Pi (log Ps) tendency)', 'Pa s-1' )
call HistoryAutoAddVariable( 'Vor', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'vorticity', 's-1' )
call HistoryAutoAddVariable( 'Div', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'divergence', 's-1' )
call HistoryAutoAddVariable( 'Mass', (/ 'lon ', 'lat ', 'time' /), 'mass', 'kg' )
call HistoryAutoAddVariable( 'KinEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'kinetic energy', 'J' )
call HistoryAutoAddVariable( 'IntEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'internal energy', 'J' )
call HistoryAutoAddVariable( 'PotEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'potential energy', 'J' )
call HistoryAutoAddVariable( 'LatEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'latent energy', 'J' )
call HistoryAutoAddVariable( 'TotEngy', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'total energy', 'J' )
call HistoryAutoAddVariable( 'Enstro', (/ 'lon ', 'lat ', 'sig ', 'time' /), 'enstrophy', 'kg' )
! Initialization of modules used in this module
!
! è³�����æ�
! Mass fixer
!
call MassFixerInit
if (FlagSLTT) then
! �»ã�����°ã���³ã�¸ã�¥æ���������³ªç§»æ�
! Semi-Lagrangian method for tracer transport
call SLTTInit
endif
! �°å� ; Print
!
call MessageNotify( 'M', module_name, '----- Initialization Messages -----' )
call MessageNotify( 'M', module_name, ' TimeIntegScheme = %c', c1 = trim( TimeIntegScheme ) )
call MessageNotify( 'M', module_name, ' HDEFoldTime = %f [%c]', d = (/ HDEFoldTimeValue /), c1 = trim(HDEFoldTimeUnit) )
call MessageNotify( 'M', module_name, ' HDOrder = %d', i = (/ HDOrder /) )
call MessageNotify( 'M', module_name, ' VisCoef = %f', d = (/ VisCoef /) )
call MessageNotify( 'M', module_name, ' FlagSpongeLayer = %b', l = (/ FlagSpongeLayer /) )
call MessageNotify( 'M', module_name, ' FlagSpongeLayerforZonalMean = %b', l = (/ FlagSpongeLayerforZonalMean /) )
call MessageNotify( 'M', module_name, ' FlagSpongeLayerforHeat = %b', l = (/ FlagSpongeLayerforHeat /) )
call MessageNotify( 'M', module_name, ' SLEFoldTime = %f [%c]', d = (/ SLEFoldTimeValue /), c1 = trim(SLEFoldTimeUnit) )
call MessageNotify( 'M', module_name, ' SLOrder = %d', i = (/ SLOrder /) )
call MessageNotify( 'M', module_name, ' SLNumLayer = %d', i = (/ SLNumLayer /) )
call MessageNotify( 'M', module_name, ' RefTemp = %f', d = (/ RefTemp /) )
call MessageNotify( 'M', module_name, ' FlagDivDamp = %b', l = (/ FlagDivDamp /) )
call MessageNotify( 'M', module_name, ' DivDampPeriod = %f', d = (/ DivDampPeriod /) )
call MessageNotify( 'M', module_name, ' FlagTotMassFix = %b', l = (/ FlagTotMassFix /) )
call MessageNotify( 'M', module_name, ' FlagMassHorDifCor = %b', l = (/ FlagMassHorDifCor /) )
call MessageNotify( 'M', module_name, ' FlagSLTT = %b', l = (/ FlagSLTT /) )
call MessageNotify( 'M', module_name, ' FlagCalcUVTPs = %b', l = (/ FlagCalcUVTPs /) )
call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )
dynamics_hspl_vas83_inited = .true.
end subroutine DynamicsHSplVAS83Init
| Variable : | |||
| DelTimeSave : | real(DP), save
|
| Variable : | |||
| DivDampPeriod : | real(DP), save
|
| Subroutine : | |||
| wz_DivA(lmax, 1:kmax) : | real(DP), intent(inout)
|
�ºæ�£ã���è¡°é������� (�´ã���¤ã�����£ã��������è¨�ç®��������¾ä¹±���è¡°ã���³å�)
Addition of divergence damping term
subroutine DivergenceDamping( wz_DivA )
!
! �ºæ�£ã���è¡°é������� (�´ã���¤ã�����£ã��������è¨�ç®��������¾ä¹±���è¡°ã���³å�)
!
! Addition of divergence damping term
!
! �¢ã�¸ã�¥ã�¼ã����� ; USE statements
!
! ���»ç���
! Time control
!
use timeset, only: DelTime, TimeN ! �¹ã������ $ t $ ������. Time of step $ t $.
! �¼å��¹è¨å®�
! Grid points settings
!
use gridset, only: lmax, imax, jmax, kmax ! ���´å±¤��.
! Number of vertical level
implicit none
real(DP), intent(inout):: wz_DivA (lmax, 1:kmax)
! $ D (t+\Delta t) $ . �ºæ�� (�¹ã��������).
! Divergence (spectral)
! ä½�æ¥å���
! Work variables
!
real(DP) :: DivDampCoef
! ���� ; Executable statement
!
if ( FlagDivDamp ) then
DivDampCoef = 1.0_DP / DelTime * ( DivDampPeriod - TimeN ) / DivDampPeriod
if ( DivDampCoef < 0.0_DP ) DivDampCoef = 0.0_DP
wz_DivA = wz_DivA / ( 1.0_DP + 2.0_DP * DelTime * DivDampCoef )
end if
end subroutine DivergenceDamping
| Variable : | |||
| FlagCalcUVTPs : | logical , save
|
| Variable : | |||
| FlagMassHorDifCor : | logical , save
|
| Variable : | |||
| FlagSLTT : | logical , save
|
| Subroutine : | |||
| xyz_Temp(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_Phi(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
|
�¼å��¹ã���¼ã�¿ã�§ã����æ¸�º¦ $ T $ ����, ��æ°´å�§ã����������� �¼å��¹ã���¼ã�¿ã���¸ã�������³ã�·ã�£ã���åº� $ Phi $ ��æ±����¾ã��.
subroutine HydroGrid( xyz_Temp, xyz_Phi )
!
! �¼å��¹ã���¼ã�¿ã�§ã����æ¸�º¦ $ T $ ����, ��æ°´å�§ã�����������
! �¼å��¹ã���¼ã�¿ã���¸ã�������³ã�·ã�£ã���åº� $ \Phi $ ��æ±����¾ã��.
!
! �¢ã�¸ã�¥ã�¼ã����� ; USE statements
!
use constants, only: CpDry
! $ C_p $ [J kg-1 K-1].
! ä¹¾ç�¥å¤§æ°�����§æ���.
! Specific heat of air at constant pressure
! �¼å��¹è¨å®�
! Grid points settings
!
use gridset, only: imax, jmax, kmax ! ���´å±¤��.
! Number of vertical level
! 宣�� ; Declaration statements
!
implicit none
real(DP), intent(in):: xyz_Temp (0:imax-1, 1:jmax, 1:kmax)
! $ T $ . æ¸�º¦. Temperature
real(DP), intent(out):: xyz_Phi (0:imax-1, 1:jmax, 1:kmax)
! $ \Phi $ . �¸ã�������³ã�·ã�£ã���åº�.
! Getpotential height
! ä½�æ¥å���
! Work variables
!
integer:: k ! ���´æ�¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in vertical direction
! ���� ; Executable statement
!
xyz_Phi(:,:,1) = CpDry * z_HydroAlpha(1) * xyz_Temp(:,:,1)
do k = 2, kmax
xyz_Phi(:,:,k) = xyz_Phi(:,:,k-1) + CpDry * z_HydroAlpha(k ) * xyz_Temp(:,:,k ) + CpDry * z_HydroBeta (k-1) * xyz_Temp(:,:,k-1)
enddo
end subroutine HydroGrid
| Variable : | |||
| IDTimeIntegScheme : | integer , save
|
| Subroutine : | |||
| xyz_UCosLatN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_VCosLatN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_VorN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_DivN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_TempN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| xy_GradLambdaPiN(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xy_GradMuPiN(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xyz_PiAdv(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_UAdvN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_VAdvN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_TempNonLinearN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyzf_QMixNonLinearN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(out)
| ||
| xyz_KinEngyN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_TempUAdvN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyz_TempVAdvN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
| ||
| xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(out)
| ||
| xy_DPiDtNG(0:imax-1, 1:jmax) : | real(DP), intent(out)
| ||
| xyzf_QMixUAdvN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(out)
| ||
| xyzf_QMixVAdvN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(out)
|
��ç·�å½¢é�� (������æ³¢é��) ���¼å��¹ä��§è�ç®����¾ã��.
Non-linear terms (non gravitational terms) are calculated on grid points
subroutine NonLinearOnGrid( xyz_UCosLatN, xyz_VCosLatN, xyz_VorN, xyz_DivN, xyz_TempN, xyzf_QMixN, xy_GradLambdaPiN, xy_GradMuPiN, xyz_PiAdv, xyz_UAdvN, xyz_VAdvN, xyz_TempNonLinearN, xyzf_QMixNonLinearN, xyz_KinEngyN, xyz_TempUAdvN, xyz_TempVAdvN, xyr_SigDotN, xy_DPiDtNG, xyzf_QMixUAdvN, xyzf_QMixVAdvN )
!
! ��ç·�å½¢é�� (������æ³¢é��) ���¼å��¹ä��§è�ç®����¾ã��.
!
! Non-linear terms (non gravitational terms) are calculated on
! grid points
!
! �¢ã�¸ã�¥ã�¼ã����� ; USE statements
!
! 座æ����¼ã�¿è¨å®�
! Axes data settings
!
use axesset, only: r_Sigma, z_DelSigma ! $ \Delta \sigma $ (�´æ��).
! $ \Delta \sigma $ (Full)
! ����å®��°è¨å®�
! Physical constants settings
!
use constants, only: CpDry, EpsV
! $ \epsilon_v $ .
! æ°´è�¸æ���å��閴�.
! Molecular weight of water vapor
! �������
! Character handling
!
use dc_string, only: LChar
! �¼å��¹è¨å®�
! Grid points settings
!
use gridset, only: imax, jmax, kmax ! ���´å±¤��.
! Number of vertical level
! çµ������¢ã���������¨å®�
! Settings of array for atmospheric composition
!
use composition, only: ncmax, IndexH2OVap
implicit none
real(DP), intent(in):: xyz_UCosLatN (0:imax-1, 1:jmax, 1:kmax)
! $ U (t) = u (t) \cos \varphi $ .
real(DP), intent(in):: xyz_VCosLatN (0:imax-1, 1:jmax, 1:kmax)
! $ V (t) = v (t) \cos \varphi $ .
real(DP), intent(in):: xyz_VorN (0:imax-1, 1:jmax, 1:kmax)
! $ \zeta (t) $ . æ¸�º¦. Vorticity
real(DP), intent(in):: xyz_DivN (0:imax-1, 1:jmax, 1:kmax)
! $ D (t) $ . �ºæ��. Divergence
real(DP), intent(in):: xyz_TempN (0:imax-1, 1:jmax, 1:kmax)
! $ T (t) $ . æ¸�º¦. Temperature
real(DP), intent(in):: xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ q (t) $ . ��. Specific humidity
real(DP), intent(in):: xy_GradLambdaPiN (0:imax-1, 1:jmax)
! $ \DP{\pi}{\lambda} $
real(DP), intent(in):: xy_GradMuPiN (0:imax-1, 1:jmax)
! $ (1-\mu^2) \DP{\pi}{\mu} $
real(DP), intent(out):: xyz_PiAdv (0:imax-1, 1:jmax, 1:kmax)
! $ \Dvect{v} \cdot \nabla \pi $ .
! $ \pi $ ��§»æµ�. Advection of $ \pi $
real(DP), intent(out):: xyz_UAdvN (0:imax-1, 1:jmax, 1:kmax)
! $ U_A (t) $ . �±è¥¿�����霡»æ���.
! Eastward advection of momentum
real(DP), intent(out):: xyz_VAdvN (0:imax-1, 1:jmax, 1:kmax)
! $ V_A (t) $ . ���������霡»æ���.
! Northward advection of momentum
real(DP), intent(out):: xyz_TempNonLinearN (0:imax-1, 1:jmax, 1:kmax)
! $ \hat{H} (t) $ . æ¸�º¦����å¤�����.
! Temperature tendency
real(DP), intent(out):: xyzf_QMixNonLinearN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ R (t) $ . �湿��������.
! Specific humidity tendency
real(DP), intent(out):: xyz_KinEngyN (0:imax-1, 1:jmax, 1:kmax)
! $ KE (t) $ . �������������¼é��.
! Kinetic energy
real(DP), intent(out):: xyz_TempUAdvN (0:imax-1, 1:jmax, 1:kmax)
! $ UT' (t) $ . æ¸�º¦�±è¥¿ç§»æ���.
! Eastward advection of temperature
real(DP), intent(out):: xyz_TempVAdvN (0:imax-1, 1:jmax, 1:kmax)
! $ VT' (t) $ . æ¸�º¦����移æ���.
! Northward advection of temperature
real(DP), intent(out):: xyr_SigDotN (0:imax-1, 1:jmax, 0:kmax)
! $ \dot{\sigma} (t) $ .
! ���´æ�. Vertical flow
real(DP), intent(out):: xy_DPiDtNG (0:imax-1, 1:jmax)
! $ Z (t) $ . �°è¡¨�¢æ��§æ����å¤���������æ³¢é��.
! Non-gravity wave component of surface pressure tendency
real(DP), intent(out):: xyzf_QMixUAdvN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ Uq (t) $ . �湿�西移��.
! Eastward advection of specific humidity
real(DP), intent(out):: xyzf_QMixVAdvN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ Vq (t) $ . �湿���移��.
! Northward advection of specific humidity
!-----------------------------------
! ä½�æ¥å���
! Work variables
!
real(DP):: xyz_PiAdvSum (0:imax-1, 1:jmax, 1:kmax)
! $ \sum_k^K(\Dvect{v}\cdot\nabla\pi)\Delta\sigma $ .
! $ \pi $ 移�������. Integral downward of advection of $ \pi $
real(DP):: xyz_DivSum (0:imax-1, 1:jmax, 1:kmax)
! $ \sum_k^K D\Delta\sigma $ .
! �ºæ�£ã�������. Integral downward of divergence
real(DP):: xyr_SigDotNonG (0:imax-1, 1:jmax, 0:kmax)
! $ \dot{\sigma} $ .
! ���´æ� (������æ³�). Vertical flow (non gravitational)
real(DP):: xyz_TempEdd (0:imax-1, 1:jmax, 1:kmax)
! $ T' = T - \bar{T} $ .
! æ¸�º¦���¾ä¹± (�´æ�°ã������). Temperature eddy (full level)
real(DP):: xyr_TempEdd (0:imax-1, 1:jmax, 0:kmax)
! $ \hat{T}' $ .
! æ¸�º¦���¾ä¹± (���´æ�°ã������). Temperature eddy (half level)
real(DP):: xyz_VirTemp (0:imax-1, 1:jmax, 1:kmax)
! $ T_v $ .
! ä»�¸©åº�. Virtual temperature
real(DP):: xyz_VirTempEdd (0:imax-1, 1:jmax, 1:kmax)
! $ {T_v}' = T_v - \bar{T} $ .
! ä»�¸©åº����¾ä¹±. Virtual temperature eddy
integer:: k ! ���´æ�¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in vertical direction
integer:: n ! çµ����¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in dimension of constituents
! ���� ; Executable statement
!
! $ \pi $ ��§»æµ�, $ \pi $ 移æ����ºæ�£ã��¤§æ°�ä¸�ç«�����ä¸����������
! Calculate advection of $ \pi $, integral advection of $ \pi $ and divergence
! from the top of the model downward
!
do k = 1, kmax
xyz_PiAdv(:,:,k) = ( xyz_UCosLatN(:,:,k) * xy_GradLambdaPiN + xyz_VCosLatN(:,:,k) * xy_GradMuPiN ) / ( 1.0_DP - xy_SinLat**2 )
enddo
xyz_PiAdvSum(:,:,kmax) = xyz_PiAdv(:,:,kmax) * z_DelSigma(kmax)
do k = kmax-1, 1, -1
xyz_PiAdvSum(:,:,k) = xyz_PiAdvSum(:,:,k+1) + xyz_PiAdv(:,:,k) * z_DelSigma(k)
enddo
xyz_DivSum(:,:,kmax) = xyz_DivN(:,:,kmax) * z_DelSigma(kmax)
do k = kmax-1, 1, -1
xyz_DivSum(:,:,k) = xyz_DivSum(:,:,k+1) + xyz_DivN(:,:,k) * z_DelSigma(k)
enddo
! �°è¡¨�¢æ��§æ����å¤�����������æ³¢é�� $ Z $ ���ç®�
! Calculate non-gravity wave component of surface pressure tendency $ Z $
!
xy_DPiDtNG = - xyz_PiAdvSum(:,:,1)
! $ \dot{\sigma} $ ����
! Calculate $ \dot{\sigma} $
!
do k = 1, kmax-1
xyr_SigDotN(:,:,k) = r_Sigma(k) * ( xyz_PiAdvSum(:,:,1) + xyz_DivSum(:,:,1) ) - ( xyz_PiAdvSum(:,:,k+1) + xyz_DivSum(:,:,k+1) )
xyr_SigDotNonG(:,:,k) = r_Sigma(k) * xyz_PiAdvSum(:,:,1) - xyz_PiAdvSum(:,:,k+1)
enddo
! $ \dot{\sigma} $ ���������
! $ \dot{\sigma} $ on upper and lower boundary
!
xyr_SigDotN (:,:,0 ) = 0.
xyr_SigDotN (:,:,kmax) = 0.
xyr_SigDotNonG(:,:,0 ) = 0.
xyr_SigDotNonG(:,:,kmax) = 0.
! æ¸�º¦���¾ä¹± (�´æ�°ã������), ä»�¸©åº�, ä»�¸©åº����¾ä¹±���ç®�
! Calculate temperature eddy (full level), virtual temperature,
! virtual temperature eddy
!
do k = 1, kmax
xyz_VirTemp(:,:,k) = xyz_TempN(:,:,k) * ( 1.0_DP + ((( 1.0_DP / EpsV ) - 1.0_DP ) * xyzf_QMixN(:,:,k,IndexH2OVap)) )
xyz_TempEdd(:,:,k) = xyz_TempN(:,:,k) - z_RefTemp(k)
xyz_VirTempEdd(:,:,k) = xyz_VirTemp(:,:,k) - z_RefTemp(k)
enddo
! æ¸�º¦���¾ä¹± (���´æ�°ã������) ���ç®�
! Calculate temperature eddy (half level)
!
xyr_TempEdd(:,:,0 ) = 0.0_DP
xyr_TempEdd(:,:,kmax) = 0.0_DP
do k = 1, kmax-1
xyr_TempEdd(:,:,k) = z_TInpCoefA(k+1) * xyz_TempN(:,:,k+1) + z_TInpCoefB(k ) * xyz_TempN(:,:,k ) - r_RefTemp(k)
enddo
! �±è¥¿�����霡»æ������ç®�
! Calculate advection of eastward momentum
!
xyz_UAdvN(:,:,1) = ( xyz_VorN(:,:,1) + xy_Cori ) * xyz_VCosLatN(:,:,1) - 1.0_DP / ( 2.0_DP * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyz_UCosLatN(:,:,1) - xyz_UCosLatN(:,:,2) ) - CpDry * z_TInpCoefK(1) * xyz_VirTempEdd(:,:,1) * xy_GradLambdaPiN
do k = 2, kmax-1
xyz_UAdvN(:,:,k) = ( xyz_VorN(:,:,k) + xy_Cori ) * xyz_VCosLatN(:,:,k) - 1.0_DP / ( 2.0_DP * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyz_UCosLatN(:,:,k-1) - xyz_UCosLatN(:,:,k) ) + xyr_SigDotN(:,:,k) * ( xyz_UCosLatN(:,:,k) - xyz_UCosLatN(:,:,k+1) ) ) - CpDry * z_TInpCoefK(k) * xyz_VirTempEdd(:,:,k) * xy_GradLambdaPiN
end do
xyz_UAdvN(:,:,kmax) = ( xyz_VorN(:,:,kmax) + xy_Cori ) * xyz_VCosLatN(:,:,kmax) - 1.0_DP / ( 2.0_DP * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyz_UCosLatN(:,:,kmax-1) - xyz_UCosLatN(:,:,kmax) ) - CpDry * z_TInpCoefK(kmax) * xyz_VirTempEdd(:,:,kmax) * xy_GradLambdaPiN
! ���������霡»æ������ç®�
! Calculate advection of northward momentum
!
xyz_VAdvN(:,:,1) = - ( xyz_VorN(:,:,1) + xy_Cori ) * xyz_UCosLatN(:,:,1) - 1.0_DP / ( 2.0_DP * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyz_VCosLatN(:,:,1) - xyz_VCosLatN(:,:,2) ) - CpDry * z_TInpCoefK(1) * xyz_VirTempEdd(:,:,1) * xy_GradMuPiN
do k = 2, kmax-1
xyz_VAdvN(:,:,k) = - ( xyz_VorN(:,:,k) + xy_Cori ) * xyz_UCosLatN(:,:,k) - 1.0_DP / ( 2.0_DP * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyz_VCosLatN(:,:,k-1) - xyz_VCosLatN(:,:,k) ) + xyr_SigDotN(:,:,k) * ( xyz_VCosLatN(:,:,k) - xyz_VCosLatN(:,:,k+1) ) ) - CpDry * z_TInpCoefK(k) * xyz_VirTempEdd(:,:,k) * xy_GradMuPiN
end do
xyz_VAdvN(:,:,kmax) = - ( xyz_VorN(:,:,kmax) + xy_Cori ) * xyz_UCosLatN(:,:,kmax) - 1.0_DP / ( 2.0_DP * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyz_VCosLatN(:,:,kmax-1) - xyz_VCosLatN(:,:,kmax) ) - CpDry * z_TInpCoefK(kmax) * xyz_VirTempEdd(:,:,kmax) * xy_GradMuPiN
! �������������¼é�� (ä»�¸©åº��æ£å����) ���ç®�
! Calculate kinematic energy term
! (including virtual temperature correction)
!
call HydroGrid( xyz_VirTemp - xyz_TempN, xyz_KinEngyN ) ! (out)
do k = 1, kmax
xyz_KinEngyN(:,:,k) = xyz_KinEngyN(:,:,k) + ( xyz_UCosLatN(:,:,k)**2 + xyz_VCosLatN(:,:,k)**2 ) / ( 2.0_DP * ( 1.0_DP - xy_SinLat**2 ) )
end do
! æ¸�º¦�±è¥¿ç§»æ���, æ¸�º¦����移æ������ç®�
! Calculate eastward and northward advection of temperature
!
do k = 1, kmax
xyz_TempUAdvN(:,:,k) = xyz_UCosLatN(:,:,k) * xyz_TempEdd(:,:,k)
xyz_TempVAdvN(:,:,k) = xyz_VCosLatN(:,:,k) * xyz_TempEdd(:,:,k)
end do
! æ¸�º¦������å¤����� $ \hat{H} $ ���ç®�
! Calculate temperature tendency term $ \hat{H} $
!
do k = 1, kmax-1
xyz_TempNonLinearN(:,:,k) = xyz_TempEdd(:,:,k) * xyz_DivN(:,:,k) - 1.0_DP / z_DelSigma(k) * ( xyr_SigDotN(:,:,k-1) * ( xyr_TempEdd(:,:,k-1) - xyz_TempEdd(:,:,k) ) + xyr_SigDotN(:,:,k) * ( xyz_TempEdd(:,:,k) - xyr_TempEdd(:,:,k) ) ) - 1.0_DP / z_DelSigma(k) * ( xyr_SigDotNonG(:,:,k-1) * ( r_RefTemp(k-1) - z_RefTemp(k) ) + xyr_SigDotNonG(:,:,k) * ( z_RefTemp(k) - r_RefTemp(k) ) ) + z_TInpCoefK(k) * xyz_VirTemp(:,:,k) * xyz_PiAdv(:,:,k) - z_HydroAlpha(k) / z_DelSigma(k) * ( xyz_VirTemp(:,:,k) * xyz_PiAdvSum(:,:,k) + xyz_VirTempEdd(:,:,k) * xyz_DivSum(:,:,k) ) - z_HydroBeta(k) / z_DelSigma(k) * ( xyz_VirTemp(:,:,k) * xyz_PiAdvSum(:,:,k+1) + xyz_VirTempEdd(:,:,k) * xyz_DivSum(:,:,k+1) )
enddo
xyz_TempNonLinearN(:,:,kmax) = xyz_TempEdd(:,:,kmax) * xyz_DivN(:,:,kmax) - 1.0_DP / z_DelSigma(kmax) * ( xyr_SigDotN(:,:,kmax-1) * ( xyr_TempEdd(:,:,kmax-1) - xyz_TempEdd(:,:,kmax) ) + xyr_SigDotN(:,:,kmax) * ( xyz_TempEdd(:,:,kmax) - xyr_TempEdd(:,:,kmax) ) ) - 1.0_DP / z_DelSigma(kmax) * ( xyr_SigDotNonG(:,:,kmax-1) * ( r_RefTemp(kmax-1) - z_RefTemp(kmax) ) + xyr_SigDotNonG(:,:,kmax) * ( z_RefTemp(kmax) - r_RefTemp(kmax) ) ) + z_TInpCoefK(kmax) * xyz_VirTemp(:,:,kmax) * xyz_PiAdv(:,:,kmax) - z_HydroAlpha(kmax) / z_DelSigma(kmax) * ( xyz_VirTemp(:,:,kmax) * xyz_PiAdvSum(:,:,kmax) + xyz_VirTempEdd(:,:,kmax) * xyz_DivSum(:,:,kmax) )
if (.not. FlagSLTT) then
! �湿�西移��, �湿���移������
! Calculate eastward and northward advection of specific humidity
!
do n = 1, ncmax
do k = 1, kmax
xyzf_QMixUAdvN(:,:,k,n) = xyz_UCosLatN(:,:,k) * xyzf_QMixN(:,:,k,n)
xyzf_QMixVAdvN(:,:,k,n) = xyz_VCosLatN(:,:,k) * xyzf_QMixN(:,:,k,n)
end do
end do
! �湿�������� $ R $ ����
! Calculate specific humidity tendency $ R $
!
do n = 1, ncmax
xyzf_QMixNonLinearN(:,:,1,n) = xyzf_QMixN(:,:,1,n) * xyz_DivN(:,:,1) - 1.0_DP / ( 2.0_DP * z_DelSigma(1) ) * xyr_SigDotN(:,:,1) * ( xyzf_QMixN(:,:,1,n) - xyzf_QMixN(:,:,2,n) )
do k = 2, kmax - 1
xyzf_QMixNonLinearN(:,:,k,n) = xyzf_QMixN(:,:,k,n) * xyz_DivN(:,:,k) - 1.0_DP / ( 2.0_DP * z_DelSigma(k) ) * ( xyr_SigDotN(:,:,k-1) * ( xyzf_QMixN(:,:,k-1,n) - xyzf_QMixN(:,:,k ,n) ) + xyr_SigDotN(:,:,k) * ( xyzf_QMixN(:,:,k ,n) - xyzf_QMixN(:,:,k+1,n) ) )
end do
xyzf_QMixNonLinearN(:,:,kmax,n) = xyzf_QMixN(:,:,kmax,n) * xyz_DivN(:,:,kmax) - 1.0_DP / ( 2.0_DP * z_DelSigma(kmax) ) * xyr_SigDotN(:,:,kmax-1) * ( xyzf_QMixN(:,:,kmax-1,n) - xyzf_QMixN(:,:,kmax,n) )
end do
endif
end subroutine NonLinearOnGrid
| Subroutine : | |||
| xyz_UB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_VB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_TempB(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| xy_PsB(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xyz_UN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_VN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_TempN(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| xy_PsN(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xyz_UA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_VA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_TempA(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| xy_PsA(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xyz_DUDtPhy(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_DVDtPhy(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_DTempDtPhy(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| xyr_SigDot(0:imax-1, 1:jmax, 0:kmax) : | real(DP), intent(in)
| ||
| xy_DPiDt(0:imax-1, 1:jmax) : | real(DP), intent(in)
| ||
| xyz_PiAdv(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in)
| ||
| xyz_OMG(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out)
|
診æ�é����ºå����è¡����¾ã��.
Diagnostic variables are output.
subroutine OutputDiagnosedVariables( xyz_UB, xyz_VB, xyz_TempB, xyzf_QMixB, xy_PsB, xyz_UN, xyz_VN, xyz_TempN, xyzf_QMixN, xy_PsN, xyz_UA, xyz_VA, xyz_TempA, xyzf_QMixA, xy_PsA, xyz_DUDtPhy, xyz_DVDtPhy, xyz_DTempDtPhy, xyzf_DQMixDtPhy, xyr_SigDot, xy_DPiDt, xyz_PiAdv, xyz_OMG )
!
! 診æ�é����ºå����è¡����¾ã��.
!
! Diagnostic variables are output.
!
! �¢ã�¸ã�¥ã�¼ã����� ; USE statements
!
! ����å®��°è¨å®�
! Physical constants settings
!
use constants, only: RPlanet, Grav, CpDry, LatentHeat
! $ L $ [J kg-1] .
! ��������.
! Latent heat of condensation
! ���»ç���
! Time control
!
use timeset, only: DelTime, TimeN ! �¹ã������ $ t $ ������. Time of step $ t $.
! �¼å��¹è¨å®�
! Grid points settings
!
use gridset, only: lmax, imax, jmax, kmax ! ���´å±¤��.
! Number of vertical level
! çµ������¢ã���������¨å®�
! Settings of array for atmospheric composition
!
use composition, only: ncmax, IndexH2OVap
! 座æ����¼ã�¿è¨å®�
! Axes data settings
!
use axesset, only: z_Sigma, r_Sigma, z_DelSigma ! $ \Delta \sigma $ (�´æ��).
! $ \Delta \sigma $ (Full)
! SPMODEL ���¤ã������, ���¢ä����馹������¢è����½æ�°å���������解ã��(å¤�層å�å¿�)
! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
!
#ifdef LIB_MPI
#ifdef SJPACK
use wa_mpi_module_sjpack, only: w_xy => w_xv, xy_GradLon_w => xv_GradLon_w, xy_GradLat_w => xv_GradLat_w, wa_DivLambda_xya => wa_DivLambda_xva, wa_DivMu_xya => wa_DivMu_xva, xya_wa => xva_wa
#else
use wa_mpi_module, only: w_xy => w_xv, xy_GradLon_w => xv_GradLon_w, xy_GradLat_w => xv_GradLat_w, wa_DivLambda_xya => wa_DivLambda_xva, wa_DivMu_xya => wa_DivMu_xva, xya_wa => xva_wa
#endif
#elif AXISYMMETRY
use wa_zonal_module, only: w_xy, xy_GradLon_w, xy_GradLat_w, wa_DivLambda_xya, wa_DivMu_xya, xya_wa
#elif SJPACK
use wa_module_sjpack, only: w_xy, xy_GradLon_w, xy_GradLat_w, wa_DivLambda_xya, wa_DivMu_xya, xya_wa
#elif AXISYMMETRY_SJPACK
use wa_zonal_module_sjpack, only: w_xy, xy_GradLon_w, xy_GradLat_w, wa_DivLambda_xya, wa_DivMu_xya, xya_wa
#else
use wa_module, only: w_xy, xy_GradLon_w, xy_GradLat_w, wa_DivLambda_xya, wa_DivMu_xya, xya_wa
#endif
! ���¹ã�������¼ã�¿å�ºå��
! History data output
!
use gtool_historyauto, only: HistoryAutoPut
! 宣�� ; Declaration statements
!
implicit none
real(DP), intent(in):: xyz_UB (0:imax-1, 1:jmax, 1:kmax)
! $ u (t-\Delta t) $ . �±è¥¿é¢���. Eastward wind
real(DP), intent(in):: xyz_VB (0:imax-1, 1:jmax, 1:kmax)
! $ v (t-\Delta t) $ . �������. Northward wind
real(DP), intent(in):: xyz_TempB (0:imax-1, 1:jmax, 1:kmax)
! $ T (t-\Delta t) $ . æ¸�º¦. Temperature
real(DP), intent(in):: xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ q (t-\Delta t) $ . ��. Specific humidity
real(DP), intent(in):: xy_PsB (0:imax-1, 1:jmax)
! $ p_s (t-\Delta t) $ . �°è¡¨�¢æ���. Surface pressure
real(DP), intent(in):: xyz_UN (0:imax-1, 1:jmax, 1:kmax)
! $ u (t) $ . �±è¥¿é¢���. Eastward wind
real(DP), intent(in):: xyz_VN (0:imax-1, 1:jmax, 1:kmax)
! $ v (t) $ . �������. Northward wind
real(DP), intent(in):: xyz_TempN (0:imax-1, 1:jmax, 1:kmax)
! $ T (t) $ . æ¸�º¦. Temperature
real(DP), intent(in):: xyzf_QMixN(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ q (t) $ . ��. Specific humidity
real(DP), intent(in):: xy_PsN (0:imax-1, 1:jmax)
! $ p_s (t) $ . �°è¡¨�¢æ���. Surface pressure
real(DP), intent(in):: xyz_UA (0:imax-1, 1:jmax, 1:kmax)
! $ u (t+\Delta t) $ . �±è¥¿é¢���. Eastward wind
real(DP), intent(in):: xyz_VA (0:imax-1, 1:jmax, 1:kmax)
! $ v (t+\Delta t) $ . �������. Northward wind
real(DP), intent(in):: xyz_TempA (0:imax-1, 1:jmax, 1:kmax)
! $ T (t+\Delta t) $ . æ¸�º¦. Temperature
real(DP), intent(in):: xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ q (t+\Delta t) $ . ��. Specific humidity
real(DP), intent(in):: xy_PsA (0:imax-1, 1:jmax)
! $ p_s (t+\Delta t) $ . �°è¡¨�¢æ���. Surface pressure
real(DP), intent(in):: xyz_DUDtPhy (0:imax-1, 1:jmax, 1:kmax)
! $ \left(\DP{u}{t}\right)^{phy} $ .
! å¤����� (�������) �������±è¥¿é¢���å¤���.
! Eastward wind tendency by external force terms (physical processes)
real(DP), intent(in):: xyz_DVDtPhy (0:imax-1, 1:jmax, 1:kmax)
! $ \left(\DP{v}{t}\right)^{phy} $ .
! ����� (�������) ����������������.
! Northward wind tendency by external force terms (physical processes)
real(DP), intent(in):: xyz_DTempDtPhy (0:imax-1, 1:jmax, 1:kmax)
! $ \left(\DP{T}{t}\right)^{phy} $ .
! å¤����� (�������) ������æ¸�º¦å¤���.
! Temperature tendency by external force terms (physical processes)
real(DP), intent(in):: xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
! $ \left(\DP{q}{t}\right)^{phy} $ .
! ����� (�������) �������湿��.
! Temperature tendency by external force terms (physical processes)
real(DP), intent(in):: xyr_SigDot (0:imax-1, 1:jmax, 0:kmax)
! $ \dot{\sigma} $ .
! ���´æ�. Vertical flow
real(DP), intent(in):: xy_DPiDt (0:imax-1, 1:jmax)
! $ Z $ . �°è¡¨�¢æ��§æ����å¤�����.
! Surface pressure tendency
real(DP), intent(in):: xyz_PiAdv (0:imax-1, 1:jmax, 1:kmax)
! $ \Dvect{v} \cdot \nabla \pi $ .
! $ \pi $ ��§»æµ�. Advection of $ \pi $
real(DP), intent(out):: xyz_OMG (0:imax-1, 1:jmax, 1:kmax)
!
! OMEGA, DP/Dt
! ä½�æ¥å���
! Work variables
!
real(DP):: xyz_DUDtDyn (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_DVDtDyn (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyz_DTempDtDyn (0:imax-1, 1:jmax, 1:kmax)
real(DP):: xyzf_DQMixDtDyn(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
real(DP):: xy_DPsDtDyn (0:imax-1, 1:jmax)
real(DP):: xyz_UCosLat (0:imax-1, 1:jmax, 1:kmax)
! $ U = u \cos \varphi $ .
real(DP):: xyz_VCosLat (0:imax-1, 1:jmax, 1:kmax)
! $ V = v \cos \varphi $ .
real(DP):: xyz_Vor (0:imax-1, 1:jmax, 1:kmax)
! $ \zeta $ . æ¸�º¦. Vorticity
real(DP):: xyz_Div (0:imax-1, 1:jmax, 1:kmax)
! $ D $ . �ºæ��. Divergence
real(DP):: xyr_PiAdvDivSum(0:imax-1, 1:jmax, 0:kmax)
real(DP):: xy_Mass (0:imax-1, 1:jmax)
! ���.
! Mass
real(DP):: xyz_KinEngy (0:imax-1, 1:jmax, 1:kmax)
! $ KE $ . ��������������.
! Kinetic energy
real(DP):: xyz_IntEngy (0:imax-1, 1:jmax, 1:kmax)
! $ IE $ . ��������������.
! Internal energy
real(DP):: xyz_PotEngy (0:imax-1, 1:jmax, 1:kmax)
! $ PE $ . �����³ã�·ã�£ã������������.
! Potential energy
real(DP):: xyz_LatEngy (0:imax-1, 1:jmax, 1:kmax)
! $ LE $ . æ½��±ã����������.
! Latent heat energy
real(DP):: xyz_TotEngy (0:imax-1, 1:jmax, 1:kmax)
! $ TE $ . ������������.
! Total energy
real(DP):: xyz_Enstro (0:imax-1, 1:jmax, 1:kmax)
! ���³ã�¹ã���ã���£ã��.
! Enstrophy
integer:: k ! ���´æ�¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in vertical direction
! ���� ; Executable statement
!
! Calculate tendencies
!
xyz_DUDtDyn = ( xyz_UA - xyz_UB ) / ( 2.0_DP * DelTime ) - xyz_DUDtPhy
xyz_DVDtDyn = ( xyz_VA - xyz_VB ) / ( 2.0_DP * DelTime ) - xyz_DVDtPhy
xyz_DTempDtDyn = ( xyz_TempA - xyz_TempB ) / ( 2.0_DP * DelTime ) - xyz_DTempDtPhy
xyzf_DQMixDtDyn = ( xyzf_QMixA - xyzf_QMixB ) / ( 2.0_DP * DelTime ) - xyzf_DQMixDtPhy
xy_DPsDtDyn = ( xy_PsA - xy_PsB ) / ( 2.0_DP * DelTime )
call HistoryAutoPut( TimeN, 'DUDtDyn' , xyz_DUDtDyn )
call HistoryAutoPut( TimeN, 'DVDtDyn' , xyz_DVDtDyn )
call HistoryAutoPut( TimeN, 'DTempDtDyn', xyz_DTempDtDyn )
call HistoryAutoPut( TimeN, 'DQVapDtDyn', xyzf_DQMixDtDyn(:,:,:,IndexH2OVap) )
call HistoryAutoPut( TimeN, 'DPsDtDyn' , xy_DPsDtDyn )
! ���´æ����°è¡¨�¢æ��§æ����å¤��������ºå��
! Output vertical flow and surface pressure tendency
!
call HistoryAutoPut( TimeN, 'SigDot', xyr_SigDot )
call HistoryAutoPut( TimeN, 'DPiDt', xy_DPiDt )
! é¢�������æ¸�º¦�ºæ�£ã���ç®�
! Calculate vorticity and divergence from wind velocity
!
do k = 1, kmax
xyz_UCosLat(:,:,k) = xyz_UN(:,:,k) * xy_CosLat
xyz_VCosLat(:,:,k) = xyz_VN(:,:,k) * xy_CosLat
end do
xyz_Vor = xya_wa( wa_DivLambda_xya( xyz_VCosLat ) - wa_DivMu_xya( xyz_UCosLat ) ) / RPlanet
xyz_Div = xya_wa( wa_DivLambda_xya( xyz_UCosLat ) + wa_DivMu_xya( xyz_VCosLat ) ) / RPlanet
call HistoryAutoPut( TimeN, 'Vor', xyz_Vor )
call HistoryAutoPut( TimeN, 'Div', xyz_Div )
! Calculation of Omega (DPressDt)
! It should be recognized that the value of xy_Ps here may have been modified
! by mass fixer.
!
! Integration from top of the model to k's layer upper interface
xyr_PiAdvDivSum(:,:,kmax) = 0.0_DP
do k = kmax-1, 0, -1
xyr_PiAdvDivSum(:,:,k) = xyr_PiAdvDivSum(:,:,k+1) + ( xyz_PiAdv(:,:,k+1) + xyz_Div(:,:,k+1) ) * z_DelSigma(k+1)
end do
do k = 1, kmax
xyz_OMG(:,:,k) = xy_PsN * ( z_Sigma(k) * xyz_PiAdv(:,:,k) - xyr_PiAdvDivSum(:,:,k) - ( xyz_PiAdv(:,:,k) + xyz_Div(:,:,k) ) * ( z_Sigma(k) - r_Sigma(k) ) )
end do
call HistoryAutoPut( TimeN, 'OMG', xyz_OMG )
! ������
! Calculate mass
!
xy_Mass = xy_PsN / Grav
! ����������, ���³ã�¹ã���ã���£ã�¼ã���ç®�
! Calculate energy and enstrophy
!
call HydroGrid( xyz_TempN, xyz_PotEngy )
do k = 1, kmax
xyz_KinEngy(:,:,k) = ( xyz_UN(:,:,k) ** 2 + xyz_VN(:,:,k) ** 2 ) / 2.0_DP * xy_Mass
xyz_IntEngy(:,:,k) = CpDry * xyz_TempN(:,:,k) * xy_Mass
xyz_PotEngy(:,:,k) = xyz_PotEngy(:,:,k) * xy_Mass
xyz_LatEngy(:,:,k) = LatentHeat * xyzf_QMixN(:,:,k,IndexH2OVap) * xy_Mass
end do
xyz_TotEngy = xyz_KinEngy + xyz_IntEngy + xyz_PotEngy + xyz_LatEngy
do k = 1, kmax
xyz_Enstro(:,:,k) = xyz_Vor(:,:,k) ** 2 * xy_Mass
end do
call HistoryAutoPut( TimeN, 'Mass', xy_Mass )
call HistoryAutoPut( TimeN, 'KinEngy', xyz_KinEngy )
call HistoryAutoPut( TimeN, 'IntEngy', xyz_IntEngy )
call HistoryAutoPut( TimeN, 'PotEngy', xyz_PotEngy )
call HistoryAutoPut( TimeN, 'LatEngy', xyz_LatEngy )
call HistoryAutoPut( TimeN, 'TotEngy', xyz_TotEngy )
call HistoryAutoPut( TimeN, 'Enstro', xyz_Enstro )
end subroutine OutputDiagnosedVariables
| Subroutine : |
TimeIntegration �§ä½¿������ä¿��°ã��¨å®���è¡����¾ã��. ���������� $ Delta t $ ��å¤��´ã�������´å��以å���, ������¨å®������¤ã�������¾ã�¾ç�����¾ã��.
Configure coefficients for "TimeIntegration". Setting values that are set last time are used, except when first time and $ Delta t $ are changed.
subroutine SemiImplMatrix
!
! TimeIntegration �§ä½¿������ä¿��°ã��¨å®���è¡����¾ã��.
! ���������� $ \Delta t $ ��å¤��´ã�������´å��以å���,
! ������¨å®������¤ã�������¾ã�¾ç�����¾ã��.
!
! Configure coefficients for "TimeIntegration".
! Setting values that are set last time are used,
! except when first time and $ \Delta t $ are changed.
!
! �¢ã�¸ã�¥ã�¼ã����� ; USE statements
!
! ����å®��°è¨å®�
! Physical constants settings
!
use constants, only: RPlanet, CpDry
! $ C_p $ [J kg-1 K-1].
! ä¹¾ç�¥å¤§æ°�����§æ���.
! Specific heat of air at constant pressure
! �¼å��¹è¨å®�
! Grid points settings
!
use gridset, only: lmax, imax, jmax, kmax ! ���´å±¤��.
! Number of vertical level
! 座æ����¼ã�¿è¨å®�
! Axes data settings
!
use axesset, only: r_Sigma, z_DelSigma ! $ \Delta \sigma $ (�´æ��).
! $ \Delta \sigma $ (Full)
! ���»ç���
! Time control
!
use timeset, only: DelTime ! $ \Delta t $ [s]
! SPMODEL ���¤ã������, ���¢ä����馹������¢è����½æ�°å���������解ã��(å¤�層å�å¿�)
! SPMODEL library, problems on sphere are solved with spherical harmonics (multi layer is supported)
!
#ifdef LIB_MPI
#ifdef SJPACK
use wa_mpi_module_sjpack, only: l_nm
#else
use wa_mpi_module, only: l_nm
#endif
#elif AXISYMMETRY
use wa_zonal_module, only: l_nm
#elif SJPACK
use wa_module_sjpack, only: l_nm
#elif AXISYMMETRY_SJPACK
use wa_zonal_module_sjpack, only: l_nm
#else
use wa_module, only: l_nm
#endif
! LU ��解æ��������£ç� 1 次æ�¹ç�å¼��解ã��������¢æ�� (spml ��梱ã�¢ã�¸ã�¥ã�¼ã��)
! Functions to solve linear simultaneous equation by LU decomposition
! (a module included in spml)
!
use lumatrix, only: LUDecomp
! �������°ç�����¼ã���£ã������
! Utilities for debug
!
use dc_trace, only: DbgMessage, BeginSub, EndSub, Debug
! 宣�� ; Declaration statements
!
implicit none
! TimeIntegration ç��§ä½¿������ä¿��°ã��¨å®����������æ¥å���
! Work variable for coefficients for "TimeIntegration", etc.
!
integer:: k, l, kk ! ���´æ�¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in vertical direction
integer:: n
! ���� ; Executable statement
!
! $ \Delta t $ [s] �����§ã�����»ä�å�
! Check and save $ \Delta t $ [s]
!
if ( DelTimeSave == DelTime ) return
DelTimeSave = DelTime
call DbgMessage( '%c: %c: (DelTime=%f [sec]) coefficients for "TimeIntegration" is generated. ', c1 = module_name, c2 = 'SemiImplMatrix', d = (/ DelTime /) )
! TimeIntegration �§ä½¿������ä¿��°ã��¨å®�
! Configure coefficients for "TimeIntegration"
!
if ( .not. allocated( z_siMtxC ) ) allocate( z_siMtxC (1:kmax) )
if ( .not. allocated( z_siMtxG ) ) allocate( z_siMtxG (1:kmax) )
if ( .not. allocated( zz_siMtxH ) ) allocate( zz_siMtxH (1:kmax, 1:kmax) )
if ( .not. allocated( wz_siMtxDi ) ) allocate( wz_siMtxDi (lmax,1:kmax) )
if ( .not. allocated( zz_siMtxDiH ) ) allocate( zz_siMtxDiH (1:kmax, 1:kmax) )
if ( .not. allocated( wzz_siMtxWDiH ) ) allocate( wzz_siMtxWDiH(lmax,1:kmax, 1:kmax) )
if ( .not. allocated( zz_siMtxGCt ) ) allocate( zz_siMtxGCt (1:kmax, 1:kmax) )
if ( .not. allocated( zz_siMtxW ) ) allocate( zz_siMtxW (1:kmax, 1:kmax) )
if ( .not. allocated( zz_siMtxQ ) ) allocate( zz_siMtxQ (1:kmax, 1:kmax) )
if ( .not. allocated( zz_siMtxS ) ) allocate( zz_siMtxS (1:kmax, 1:kmax) )
if ( .not. allocated( zz_siMtxR ) ) allocate( zz_siMtxR (1:kmax, 1:kmax) )
z_siMtxC = z_DelSigma
z_siMtxG = CpDry * z_TInpCoefK * z_RefTemp
do k = 1, kmax
do l = 1, kmax
zz_siMtxGCt(k,l) = z_siMtxG(k) * z_siMtxC(l)
end do
end do
zz_siMtxW = 0.
do k = 1, kmax
do l = 1, k
zz_siMtxW(k,l) = CpDry * z_HydroAlpha(l)
enddo
do l = 1, k-1
zz_siMtxW(k,l) = zz_siMtxW(k,l) + CpDry * z_HydroBeta(l)
enddo
enddo
zz_siMtxS = 0.
do k = 1, kmax
do l = 1, kmax
zz_siMtxS(k,l) = r_Sigma(k-1) * z_DelSigma(l)
enddo
do l = k, kmax
zz_siMtxS(k,l) = zz_siMtxS(k,l) - z_DelSigma(l)
enddo
enddo
zz_siMtxQ = 0.
do k = 1, kmax
zz_siMtxQ(k,k) = ( r_RefTemp(k-1) - z_RefTemp(k) ) / z_DelSigma(k)
enddo
do k = 1, kmax-1
zz_siMtxQ(k,k+1) = ( z_RefTemp(k) - r_RefTemp(k) ) / z_DelSigma(k)
enddo
zz_siMtxR = 0.
do k = 1, kmax
do l = k, kmax
zz_siMtxR(k,l) = - z_HydroAlpha(k) / z_DelSigma(k) * z_DelSigma(l) * z_RefTemp(k)
enddo
do l = k + 1, kmax
zz_siMtxR(k,l) = zz_siMtxR(k,l) - z_HydroBeta(k) / z_DelSigma(k) * z_DelSigma(l) * z_RefTemp(k)
enddo
enddo
zz_siMtxH = matmul(zz_siMtxQ, zz_siMtxS) - zz_siMtxR
! Check of coefficients for horizontal diffusion and sponge layer
! A threshold value used below is arbitrary.
!
do n = 1, lmax
do k = 1, kmax
if ( abs( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefM(n,k) ) < 1.0e-10_DP ) then
call MessageNotify( 'E', module_name, 'Dissipation coefficient for momentum is inappropriate.' )
end if
if ( abs( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefH(n,k) ) < 1.0e-10_DP ) then
call MessageNotify( 'E', module_name, 'Dissipation coefficient for heat is inappropriate.' )
end if
end do
if ( abs( 1.0_DP - 2.0_DP * DelTime * w_DisCoefQ(n) ) < 1.0e-10_DP ) then
call MessageNotify( 'E', module_name, 'Dissipation coefficient for composition is inappropriate.' )
end if
end do
wz_siMtxDi = 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefH )
do n = 1, lmax
! NOTE: wz_siMtiDi is a diagonal matrix.
!
do k = 1, kmax
do kk = 1, kmax
zz_siMtxDiH(k,kk) = wz_siMtxDi(n,k) * zz_siMtxH(k,kk)
end do
end do
zz_siMtxDiH = matmul( zz_siMtxW, zz_siMtxDiH )
do k = 1, kmax
do kk = 1, kmax
wzz_siMtxWDiH(n,k,kk) = zz_siMtxDiH(k,kk)
end do
end do
end do
if ( .not. allocated(wzz_siMtxLU) ) allocate( wzz_siMtxLU(lmax, 1:kmax, 1:kmax) )
if ( .not. allocated(wz_siMtxPiv) ) allocate( wz_siMtxPiv(lmax, 1:kmax) )
! ��� $ \underline{M} $ ����
! Calculate matrix $ \underline{M} $
!
do k = 1, kmax
do kk = 1, kmax
wzz_siMtxLU ( :,k,kk ) = - DelTime**2 * w_LaplaEigVal(:) * ( wzz_siMtxWDiH(:,k,kk) + zz_siMtxGCt(k,kk) )
if ( k == kk ) then
wzz_siMtxLU ( :,k,kk ) = wzz_siMtxLU ( :,k,kk ) + ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefM (:,k) )
endif
end do
end do
! LU �����
! LU matrix calculation
!
call LUDecomp( wzz_siMtxLU, wz_siMtxPiv ) ! (out)
end subroutine SemiImplMatrix
| Subroutine : | |||
| w_SurfGeoPot(lmax) : | real(DP), intent(in)
| ||
| wz_DVorDtNG(lmax, 1:kmax) : | real(DP), intent(in)
| ||
| wz_DDivDtNG(lmax, 1:kmax) : | real(DP), intent(in)
| ||
| wz_DTempDtNG(lmax, 1:kmax) : | real(DP), intent(in)
| ||
| wzf_DQMixDtN(lmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| w_DPiDtNG(lmax) : | real(DP), intent(in)
| ||
| wz_DivN(lmax, 1:kmax) : | real(DP), intent(in)
| ||
| wz_TempN(lmax, 1:kmax) : | real(DP), intent(in)
| ||
| w_PiN(lmax) : | real(DP), intent(in)
| ||
| wz_VorB(lmax, 1:kmax) : | real(DP), intent(in)
| ||
| wz_DivB(lmax, 1:kmax) : | real(DP), intent(in)
| ||
| wz_TempB(lmax, 1:kmax) : | real(DP), intent(in)
| ||
| wzf_QMixB(lmax, 1:kmax, 1:ncmax) : | real(DP), intent(in)
| ||
| w_PiB(lmax) : | real(DP), intent(in)
| ||
| wz_VorA(lmax, 1:kmax) : | real(DP), intent(out)
| ||
| wz_DivA(lmax, 1:kmax) : | real(DP), intent(out)
| ||
| wz_TempA(lmax, 1:kmax) : | real(DP), intent(out)
| ||
| wzf_QMixA(lmax, 1:kmax, 1:ncmax) : | real(DP), intent(out)
| ||
| w_PiA(lmax) : | real(DP), intent(out)
|
����ç©�����è¡���, ���� $ t $ �������������å¤����� $ t-\Delta t$ ����������� ���� $ t+Delta t $ ���������è¨�ç®����¾ã��.
����ç©���æ³��������¼ã�����ã���°ã�¹ã�ã�¼ã�������������¾ã��. �������¡æ�£é��������¹å·®�������������¾ã��. �����������§ã��, $ Delta t $ ��大ã�������������, ����æ³¢é���� �»ã���¤ã�³ã�����·ã����æ³��������������¾ã��. NAMELIST#dynamics_hspl_vas83_nml �� TimeIntegScheme ��å¤��´ã����������, ����æ³¢é���������¹ã�����·ã����æ³������£ã��§£����������½ã�§ã��.
With time integration, physical values at $ t+Delta t $ is calculated from tendency at $ t $ and physical values at $ t-\Delta t $ .
Leap-frog scheme is used as time integration scheme. And backward scheme is used for diffusion terms. As default setting, semi-implicit scheme is applied to gravitational terms in order to enlarge the value of $ Delta t $ . Explicit scheme can be applied to gravitational terms by changing "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".
subroutine TimeIntegration( w_SurfGeoPot, wz_DVorDtNG, wz_DDivDtNG, wz_DTempDtNG, wzf_DQMixDtN, w_DPiDtNG, wz_DivN, wz_TempN, w_PiN, wz_VorB, wz_DivB, wz_TempB, wzf_QMixB, w_PiB, wz_VorA, wz_DivA, wz_TempA, wzf_QMixA, w_PiA )
!
! ������������,
! ���� $ t $ ������������������ $ t-\Delta t$ �����������
! ���� $ t+\Delta t $ ���������è¨�ç®����¾ã��.
!
! ����ç©���æ³��������¼ã�����ã���°ã�¹ã�ã�¼ã�������������¾ã��.
! �������¡æ�£é��������¹å·®�������������¾ã��.
! �����������§ã��, $ \Delta t $ ��大ã�������������, ����æ³¢é����
! �»ã���¤ã�³ã�����·ã����æ³��������������¾ã��.
! NAMELIST#dynamics_hspl_vas83_nml �� TimeIntegScheme ��å¤��´ã����������,
! ����æ³¢é���������¹ã�����·ã����æ³������£ã��§£����������½ã�§ã��.
!
! With time integration, physical values at $ t+\Delta t $ is calculated
! from tendency at $ t $ and physical values at $ t-\Delta t $ .
!
! Leap-frog scheme is used as time integration scheme.
! And backward scheme is used for diffusion terms.
! As default setting, semi-implicit scheme is applied to gravitational terms
! in order to enlarge the value of $ \Delta t $ .
! Explicit scheme can be applied to gravitational terms by changing
! "TimeIntegScheme" in "NAMELIST#dynamics_hspl_vas83_nml".
!
! �¢ã�¸ã�¥ã�¼ã����� ; USE statements
!
! ����å®��°è¨å®�
! Physical constants settings
!
use constants, only: RPlanet, CpDry
! $ C_p $ [J kg-1 K-1].
! ä¹¾ç�¥å¤§æ°�����§æ���.
! Specific heat of air at constant pressure
! ���»ç���
! Time control
!
use timeset, only: DelTime ! $ \Delta t $ [s]
! 座æ����¼ã�¿è¨å®�
! Axes data settings
!
use axesset, only: z_DelSigma ! $ \Delta \sigma $ (�´æ��).
! $ \Delta \sigma $ (Full)
! LU ��解æ��������£ç� 1 次æ�¹ç�å¼��解ã��������¢æ�� (spml ��梱ã�¢ã�¸ã�¥ã�¼ã��)
! Functions to solve linear simultaneous equation by LU decomposition
! (a module included in spml)
!
use lumatrix, only: LUSolve
! �������
! Character handling
!
use dc_string, only: LChar
! �¼å��¹è¨å®�
! Grid points settings
!
use gridset, only: lmax, imax, jmax, kmax ! ���´å±¤��.
! Number of vertical level
! çµ������¢ã���������¨å®�
! Settings of array for atmospheric composition
!
use composition, only: ncmax, IndexH2OVap
implicit none
real(DP), intent(in):: w_SurfGeoPot (lmax)
! $ \Phi_s $ . �°è¡¨�¸ã�������³ã�·ã�£ã��.
! Surface geo-potential
real(DP), intent(in):: wz_DVorDtNG (lmax, 1:kmax)
! $ \left( \DD{\zeta}{t} (t) \right)^{NG}$ . æ¸�º¦å¤�����������æ³¢æ���� (�¹ã��������).
! Non-gravity wave component of vorticity tendency (spectral)
real(DP), intent(in):: wz_DDivDtNG (lmax, 1:kmax)
! $ \DD{D}{t} (t) $ . �ºæ�£å�����������æ³¢æ���� (�¹ã��������).
! Divergence tendency (spectral)
real(DP), intent(in):: wz_DTempDtNG(lmax, 1:kmax)
! $ \left( \DD{T}{t} (t) \right)^{NG}$ . æ¸�º¦å¤�����������æ³¢æ���� (�¹ã��������).
! Temperature tendency (spectral)
real(DP), intent(in):: wzf_DQMixDtN(lmax, 1:kmax, 1:ncmax)
! $ \DD{q}{t} (t) $ . æ¯�湿å��� (�¹ã��������).
! Specific humidity tendency (spectral)
real(DP), intent(in):: w_DPiDtNG(lmax)
! $ \left( \DD{p_s}{t} (t) \right) $ . �°è¡¨�¢æ��§å�����������æ³¢é�� (�¹ã��������).
! Non-gravity wave component of surface pressure tendency (spectral)
real(DP), intent(in):: wz_DivN (lmax, 1:kmax)
! $ D (t) $ . �ºæ�� (�¹ã��������).
! Divergence (spectral)
real(DP), intent(in):: wz_TempN (lmax, 1:kmax)
! $ T (t) $ . æ¸�º¦ (�¹ã��������).
! Temperature (spectral)
real(DP), intent(in):: w_PiN (lmax)
! $ \pi = \ln p_s (t) $ . �°è¡¨�¢æ��� (�¹ã��������).
real(DP), intent(in):: wz_VorB (lmax, 1:kmax)
! $ \zeta (t-\Delta t) $ . æ¸�º¦ (�¹ã��������).
! Vorticity (spectral)
real(DP), intent(in):: wz_DivB (lmax, 1:kmax)
! $ D (t-\Delta t) $ . �ºæ�� (�¹ã��������).
! Divergence (spectral)
real(DP), intent(in):: wz_TempB (lmax, 1:kmax)
! $ T (t-\Delta t) $ . æ¸�º¦ (�¹ã��������).
! Temperature (spectral)
real(DP), intent(in):: w_PiB (lmax)
! $ \pi = \ln p_s (t-\Delta t) $ . �°è¡¨�¢æ��� (�¹ã��������).
! Surface pressure (spectral)
real(DP), intent(in):: wzf_QMixB(lmax, 1:kmax, 1:ncmax)
! $ q (t-\Delta t) $ . æ¯�æ¹� (�¹ã��������).
! Specific humidity (spectral)
real(DP), intent(out):: wz_VorA (lmax, 1:kmax)
! $ \zeta (t+\Delta t) $ . æ¸�º¦ (�¹ã��������).
! Vorticity (spectral)
real(DP), intent(out):: wz_DivA (lmax, 1:kmax)
! $ D (t+\Delta t) $ . �ºæ�� (�¹ã��������).
! Divergence (spectral)
real(DP), intent(out):: wz_TempA (lmax, 1:kmax)
! $ T (t+\Delta t) $ . æ¸�º¦ (�¹ã��������).
! Temperature (spectral)
real(DP), intent(out):: w_PiA (lmax)
! $ \pi = \ln p_s (t+\Delta t) $ . �°è¡¨�¢æ��� (�¹ã��������).
! Surface pressure (spectral)
real(DP), intent(out):: wzf_QMixA(lmax, 1:kmax, 1:ncmax)
! $ q (t+\Delta t) $ . æ¯�æ¹� (�¹ã��������).
! Specific humidity (spectral)
! ä½�æ¥å���
! Work variables
!
real(DP):: wz_HDiv (lmax, 1:kmax)
real(DP):: w_CtDiv (lmax)
real(DP):: wz_WT (lmax, 1:kmax)
real(DP):: wz_siTemp (lmax, 1:kmax)
! æ¸�º¦ (�»ã���¤ã�³ã�����·ã����æ³����������æ¥å���).
! Temperature (work variable for semi-implicit scheme)
real(DP):: w_siPi (lmax)
! $ \pi $ (�»ã���¤ã�³ã�����·ã����æ³����������æ¥å���).
! $ \pi $ (work variable for semi-implicit scheme)
real(DP):: wz_siPhi (lmax, 1:kmax)
! $ \Phi = \underline{W} ( 1 - 2 \Delta t \underline{D_H} ) \overline{ \Dvect{T} }^{t}$ .
! (�»ã���¤ã�³ã�����·ã����æ³����������æ¥å���).
! (Work variable for semi-implicit scheme)
real(DP):: wz_siVectF (lmax, 1:kmax)
! $ \Dvect{f} $ .
! �ºæ�£é�����¢ã�����»ã���¤ã�³ã�����·ã�����¹ç�å¼��勈¾º.
! Right-hand side of a semi-implicit equation of a divergence term.
real(DP):: wz_siDivAvrTime (lmax, 1:kmax)
! $ \overline{\Dvect{D}}^{t} $ .
! ����å¹³å���� $ \Dvect{D} $ (�»ã���¤ã�³ã�����·ã����æ³����������æ¥å���).
! Time average $ \Dvect{D} $ (a work variable for semi-implicit scheme)
integer:: k, kk ! ���´æ�¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in vertical direction
integer:: n ! çµ����¹å�������� DO ���¼ã�����æ¥å���
! Work variables for DO loop in dimension of constituents
! ���� ; Executable statement
!
! ����ç©���. �¡æ�£ã����¹å·®��
! Time integration. Backward difference is applied to diffusion
!
! æ¸�º¦ ; Vorticity
!
wz_VorA = ( 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefM ) ) * ( wz_VorB + 2.0_DP * DelTime * wz_DVorDtNG )
! �ºæ�� ; Divergence
!
select case ( IDTimeIntegScheme )
case ( IDTimeIntegSchemeSemiImplicit )
! �»ã���¤ã�³ã�����·ã����æ³��§ç������è¡������ç®���½¿����������ç®�
! Calculate terms used in making a matrix for semi-implicit method
!
wz_siTemp = ( 1.0_DP - DelTime * wz_DisCoefH ) * wz_TempB + DelTime * wz_DTempDtNG
! NOTE:
! The matrix wz_siMtxDi = ( 1 - 2 \Delta t D_H )^{-1} is diagonal matrix.
!
wz_siPhi (:,1) = CpDry * z_HydroAlpha(1) * wz_siMtxDi(:,1) * wz_siTemp(:,1)
do k = 2, kmax
wz_siPhi (:,k) = wz_siPhi(:,k-1) + CpDry * z_HydroAlpha(k ) * wz_siMtxDi(:,k ) * wz_siTemp(:,k ) + CpDry * z_HydroBeta (k-1) * wz_siMtxDi(:,k-1) * wz_siTemp(:,k-1)
end do
w_siPi = w_PiB + DelTime * w_DPiDtNG
! �ºæ�£æ�¹ç�å¼��勈¾º���ç®�
! Calculate right side of divergence equation
!
do k = 1, kmax
wz_siVectF(:,k) = ( 1.0_DP - DelTime * wz_DisCoefM (:,k) ) * wz_DivB(:,k) + DelTime * wz_DDivDtNG(:,k) - DelTime * w_LaplaEigVal(:) * ( w_SurfGeoPot + wz_siPhi(:,k) + z_siMtxG(k) * w_siPi )
end do
! ����å¹³å���� $ \Dvect{D} $ �� LU è¡����§è§£��
! Solve time average $ \Dvect{D} $ with LU matrix
!
wz_siDivAvrTime = LUSolve( wzz_siMtxLU, wz_siMtxPiv, wz_siVectF )
wz_DivA = 2. * wz_siDivAvrTime - wz_DivB
case ( IDTimeIntegSchemeExplicit )
wz_WT = 0.0_DP
do k = 1, kmax
do kk = 1, kmax
wz_WT(:,k) = wz_WT(:,k) + zz_siMtxW(k,kk) * wz_TempN(:,kk)
end do
end do
do k = 1, kmax
wz_DivA(:,k) = ( 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefM(:,k) ) ) * ( wz_DivB(:,k) + 2.0_DP * DelTime * ( wz_DDivDtNG(:,k) - w_LaplaEigVal(:) * ( w_SurfGeoPot + wz_WT(:,k) + z_siMtxG(k) * w_PiN ) ) )
end do
end select
! æ¸�º¦ ; Temperature
!
select case ( IDTimeIntegScheme )
case ( IDTimeIntegSchemeSemiImplicit )
wz_HDiv = 0.
do k = 1, kmax
do kk = 1, kmax
wz_HDiv(:,k) = wz_HDiv(:,k) + zz_siMtxH(k,kk) * wz_siDivAvrTime(:,kk)
end do
end do
case ( IDTimeIntegSchemeExplicit )
wz_HDiv = 0.
do k = 1, kmax
do kk = 1, kmax
wz_HDiv(:,k) = wz_HDiv(:,k) + zz_siMtxH(k,kk) * wz_DivN(:,kk)
end do
end do
end select
wz_TempA = ( 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * wz_DisCoefH ) ) * ( wz_TempB + 2.0_DP * DelTime * ( wz_DTempDtNG - wz_HDiv ) )
! �°è¡¨�¢æ��� ; Surface pressure
!
select case ( IDTimeIntegScheme )
case ( IDTimeIntegSchemeSemiImplicit )
w_CtDiv = 0.0_DP
do k = 1, kmax
w_CtDiv = w_CtDiv + z_siMtxC(k) * wz_siDivAvrTime(:,k)
end do
case ( IDTimeIntegSchemeExplicit )
w_CtDiv = 0.0_DP
do k = 1, kmax
w_CtDiv = w_CtDiv + z_siMtxC(k) * wz_DivN(:,k)
end do
end select
w_PiA = w_PiB + 2.0_DP * DelTime * ( w_DPiDtNG - w_CtDiv )
if (.not. FlagSLTT) then
! �� ; Specific humidity
!
do n = 1, ncmax
do k = 1, kmax
wzf_QMixA(:,k,n) = ( 1.0_DP / ( 1.0_DP - 2.0_DP * DelTime * w_DisCoefQ ) ) * ( wzf_QMixB(:,k,n) + 2.0_DP * DelTime * wzf_DQMixDtN(:,k,n) )
end do
end do
endif
end subroutine TimeIntegration
| Variable : | |||
| dynamics_hspl_vas83_inited = .false. : | logical, save
|
| Constant : | |||
| module_name = ‘dynamics_hspl_vas83‘ : | character(*), parameter
|
| Variable : | |||
| r_RefTemp(:) : | real(DP), allocatable
|
| Constant : | |||
| version = ’$Name: $’ // ’$Id: dynamics_hspl_vas83.F90,v 1.83 2014/05/07 09:39:17 murashin Exp $’ : | character(*), parameter
|
| Variable : | |||
| w_DisCoefQ(:) : | real(DP), save, allocatable
|
| Variable : | |||
| w_LaplaEigVal(:) : | real(DP), save, allocatable
|
| Variable : | |||
| wz_DisCoefH(:,:) : | real(DP), save, allocatable
|
| Variable : | |||
| wz_DisCoefM(:,:) : | real(DP), save, allocatable
|
| Variable : | |||
| wz_siMtxDi(:,:) : | real(DP), save, allocatable
|
| Variable : | |||
| wz_siMtxPiv(:,:) : | integer , save, allocatable
|
| Variable : | |||
| wzz_siMtxLU(:,:,:) : | real(DP), save, allocatable
|
| Variable : | |||
| xy_Cori(:,:) : | real(DP), allocatable
|
| Variable : | |||
| z_HydroAlpha(:) : | real(DP), allocatable
|
| Variable : | |||
| z_HydroBeta(:) : | real(DP), allocatable
|
| Variable : | |||
| z_RefTemp(:) : | real(DP), allocatable
|
| Variable : | |||
| z_TInpCoefA(:) : | real(DP), allocatable
|
| Variable : | |||
| z_TInpCoefB(:) : | real(DP), allocatable
|
| Variable : | |||
| z_TInpCoefK(:) : | real(DP), allocatable
|
| Variable : | |||
| zz_siMtxGCt(:,:) : | real(DP), save, allocatable
|
| Variable : | |||||||||
| zz_siMtxR(:,:) : | real(DP), save, allocatable
|
| Variable : | |||
| zz_siMtxW(:,:) : | real(DP), save, allocatable
|