Class sltt
In: sltt/sltt.F90

��³ªç§»æ� (�»ã�����°ã���³ã�¸ã�¥æ�, Enomoto (2008) modified)

Tracer Transport (Semi-Lagrangian method, Enomoto (2008) modified)

Note that Japanese and English are described in parallel.

��³ªç§»æ�����ä¿�å­������»ã�����°ã���³ã�¸ã�¥æ��§æ�ç®������¢ã�¸ã�¥ã�¼ã���§ã��. ä¸�æµ��¹æŽ¢ç´¢ã���� Williamson and Rasch (1989, MWR) �� è£������� Enomoto (2008) ��å¿��������¹æ������������¾ã���� ������¡ã�¹ã������������æ±�����ï¼���¾®�����¤ã����������ï¼�次精åº�������������¼ã��è£����§ã���� ��è²���ä¿�証ã���������� arcsine å¤������£ã���¿ã�����������¾ã���� �¹ã������������»é�精度è£������±æ�¥ã����人工�����­æ³¢���¤åŽ»���������� Sun et al. (1996) �� ��調ã���£ã���¿ã��å¿���������������������������������

This is a tracer transport module. Semi-Lagrangian method (Enomoto 2008 modified) Arcsine transformation filter is used to avoid negative values. Monotonicity filter (Sun et al 1996) is partly used.

Procedures List

SLTTMain :移æ�è¨�ç®�
SLTTInit :������
SLTTTest :移æ����¹ã����
——————— :————
SLTTMain :Main subroutine for SLTT
SLTTInit :Initialization for SLTT
SLTTTest :Generate velocity for SLTT Test

NAMELIST

NAMELIST#

References

  • Kashimura, H., T. Enomoto, Y. O. Takahashi, 2013: Non-negative filter using arcsine transformation for tracer advection with semi-Lagrangian scheme. NCTAM, 62.
  • Enomoto, T., 2008: Bicubic Interpolation with Spectral Derivatives. SOLA, 4, 5-8. doi:10.2151/sola.2008-002
  • Williamson, D. L., and Rasch, P. J., 1989: Two-dimensional semi-Lagrangian transport with shape-preserving interpolation. Mon. Wea. Rev., 117, 102-129.
  • Sun, W.-Y., Yeh, K.-S., and Sun, R.-Y., 1996: A simple semi-Lagrangian scheme for advection equations. Quarterly Journal of the Royal Meteorological Society, 122(533), 1211-1226. doi:10.1002/qj.49712253310

Methods

FlagSLTTArcsine   SLTTHorAdv   SLTTInit   SLTTIntHor   SLTTIntVer   SLTTMain   SLTTTest   SLTTVerAdv   module_name   sltt_inited   version   x_CosLonN   x_CosLonS   x_ExtLonN   x_ExtLonS   x_LonN   x_LonS   x_SinLonN   x_SinLonS   y_CosLatN   y_CosLatS   y_ExtLatN   y_ExtLatS   y_LatN   y_LatS   y_SinLatN   y_SinLatS  

Included Modules

dc_types dc_message mpi_wrapper timeset gridset composition mass_fixer axesset sltt_const sltt_extarr sltt_dp sltt_lagint wa_mpi_module_sjpack wa_mpi_module wa_zonal_module wa_module_sjpack wa_zonal_module_sjpack wa_module namelist_util dc_iounit constants0 constants

Public Instance methods

SLTTInit( )
Subroutine :

�»ã�����°ã���³ã�¸ã�¥æ������������� Initialization for Semi-Lagrangian method

This procedure input/output NAMELIST#sltt_nml .

[Source]

  subroutine SLTTInit
    ! �»ã�����°ã���³ã�¸ã�¥æ�������������
    ! Initialization for Semi-Lagrangian method


    use axesset, only : x_Lon, y_Lat


    ! 座æ����¼ã�¿è¨­å®�
    ! Axes data settings
    !
    use axesset, only: r_Sigma, z_Sigma               ! $ \sigma $ ������ (�´æ��).
                              ! Full $ \sigma $ level

    use sltt_const , only : SLTTConstInit
    use sltt_extarr, only : SLTTExtArrInit


    ! NAMELIST ���¡ã�¤ã���¥å�����¢ã�������¼ã���£ã������
    ! Utilities for NAMELIST file input
    !
    use namelist_util, only: namelist_filename, NmlutilMsg

    ! ç¨��¥å�������¡ã��
    ! Kind type parameter
    !
    use dc_types, only: STDOUT, STRING                ! �����.       Strings. 
    ! ���¡ã�¤ã���¥å�ºå��è£���
    ! File I/O support
    !
    use dc_iounit, only: FileOpen

    use sltt_const , only : iexmin, iexmax, jexmins, jexmaxs, jexminn, jexmaxn

    !
    ! local variables
    !
    integer:: i               ! �±è¥¿�¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in zonal direction
    integer:: j               ! �����¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in meridional direction
    integer:: k               ! ���´æ�¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in vertical direction

    integer:: unit_nml        ! NAMELIST ���¡ã�¤ã�����¼ã���³ç���ç½����. 
                              ! Unit number for NAMELIST file open
    integer:: iostat_nml      ! NAMELIST 読ã�¿è¾¼�¿æ���� IOSTAT. 
                              ! IOSTAT of NAMELIST read
    ! NAMELIST å¤��°ç¾¤
    ! NAMELIST group name
    !
    namelist /sltt_nml/ FlagSLTTArcsine, SLTTIntHor, SLTTIntVer

    ! ���� ; Executable statement
    !

    if ( sltt_inited ) return

    if ( mod( jmax, 2 ) /= 0 ) then
      stop 'jmax cannot be divided by 2.'
    end if

    call SLTTConstInit


    ! �����������¤ã��¨­å®�
    ! Default values settings
    !
    FlagSLTTArcsine             = .true.
    SLTTIntHor                  = "HQ"
    SLTTIntVer                  = "HQ"


    ! 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 = sltt_nml, iostat = iostat_nml )        ! (out)
      close( unit_nml )

      call NmlutilMsg( iostat_nml, module_name ) ! (in)
      if ( iostat_nml == 0 ) write( STDOUT, nml = sltt_nml )
    end if



    allocate( x_LonS   (0:imax-1) )
    allocate( x_SinLonS(0:imax-1) )
    allocate( x_CosLonS(0:imax-1) )
    allocate( y_latS   (1:jmax/2) )
    allocate( y_SinLatS(1:jmax/2) )
    allocate( y_CosLatS(1:jmax/2) )
    do i = 0, imax-1
      x_LonS   (i) = x_Lon(i)
      x_SinLonS(i) = sin( x_LonS(i) )
      x_CosLonS(i) = cos( x_LonS(i) )
    end do
    do j = 1, jmax/2
      y_LatS   (j) = y_Lat(j)
      y_SinLatS(j) = sin( y_LatS(j) )
      y_CosLatS(j) = cos( y_LatS(j) )
    end do

    allocate( x_LonN   (0:imax-1) )
    allocate( x_SinLonN(0:imax-1) )
    allocate( x_CosLonN(0:imax-1) )
    allocate( y_latN   (1:jmax/2) )
    allocate( y_SinLatN(1:jmax/2) )
    allocate( y_CosLatN(1:jmax/2) )
    do i = 0, imax-1
      x_LonN   (i) = x_Lon(i)
      x_SinLonN(i) = sin( x_LonN(i) )
      x_CosLonN(i) = cos( x_LonN(i) )
    end do
    do j = 1, jmax/2
      y_LatN   (j) = y_Lat(j+jmax/2)
      y_SinLatN(j) = sin( y_LatN(j) )
      y_CosLatN(j) = cos( y_LatN(j) )
    end do

    allocate( x_ExtLonS( iexmin:iexmax ) )
    allocate( x_ExtLonN( iexmin:iexmax ) )

    allocate( y_ExtLatS( jexmins:jexmaxs ) )
    allocate( y_ExtLatN( jexminn:jexmaxn ) )


    call SLTTExtArrInit( x_LonS, y_LatS, x_LonN, y_LatN, x_ExtLonS, y_ExtLatS, x_ExtLonN, y_ExtLatN )


    sltt_inited = .true.

  end subroutine SLTTInit
SLTTMain( xyr_PressB, xyr_PressA, xyz_UN, xyz_VN, xyr_SigDotN, xyzf_DQMixDtPhy, xyzf_QMixB, xyzf_QMixA )
Subroutine :
xyr_PressB(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in )
: Pressure at current time step
xyr_PressA(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in )
: Pressure at next time step
xyz_UN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
: �±è¥¿é¢��� Zonal Wind
xyz_VN(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
: ������� Meridional Wind
xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in )
: ���´æ���ï¼�SigmaDotï¼�
xyzf_DQMixDtPhy(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in )
: $ left(DP{q}{t}right)^{phy} $ . ����� (�������) �������湿��. Temperature tendency by external force terms (physical processes)
xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in )
: ��³ªæ··å��æ¯� Mix ratio of the tracers
xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(out)
: ��³ªæ··å��æ¯� Mix ratio of the tracers

�»ã�����°ã���³ã�¸ã�¥æ���������³ªç§»æ�è¨�ç®���è¡����� Calculates tracer transports by Semi-Lagrangian method

[Source]

  subroutine SLTTMain( xyr_PressB, xyr_PressA, xyz_UN, xyz_VN, xyr_SigDotN, xyzf_DQMixDtPhy, xyzf_QMixB, xyzf_QMixA )
    ! �»ã�����°ã���³ã�¸ã�¥æ���������³ªç§»æ�è¨�ç®���è¡�����
    ! Calculates tracer transports by Semi-Lagrangian method


!!$    ! 座æ����¼ã�¿è¨­å®�
!!$    ! Axes data settings
!!$    !
!!$    use axesset, only: &
!!$      & z_DelSigma            ! $ \Delta \sigma $ (�´æ��).
!!$                              ! $ \Delta \sigma $ (Full)

    real(DP), intent(in ) :: xyr_PressB(0:imax-1, 1:jmax, 0:kmax)
                              !
                              ! Pressure at current time step
    real(DP), intent(in ) :: xyr_PressA(0:imax-1, 1:jmax, 0:kmax)
                              !
                              ! Pressure at next time step
    real(DP), intent(in ) :: xyz_UN    (0:imax-1, 1:jmax, 1:kmax)
                              ! �±è¥¿é¢���
                              ! Zonal Wind
    real(DP), intent(in ) :: xyz_VN    (0:imax-1, 1:jmax, 1:kmax)
                              ! �������
                              ! Meridional Wind
    real(DP), intent(in ) :: xyr_SigDotN(0:imax-1, 1:jmax, 0:kmax)
                              ! ���´æ���ï¼�SigmaDotï¼�
    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 ) :: xyzf_QMixB(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯�
                              ! Mix ratio of the tracers
    real(DP), intent(out) :: xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯�
                              ! Mix ratio of the tracers

    ! �業��
    ! Work variables
    !
    real(DP) :: f_QMixMax(1:ncmax)
                              ! ����³ªæ··å��æ¯�����大å��
                              ! Maximum of each mix ratio of the tracers
    real(DP) :: f_QMixProcMax(1:ncmax)
                              ! ����³ªæ··å��æ¯������­ã�»ã�¹å����大å��
                              ! Maximum of each mix ratio of the tracers in each process

    integer:: n               ! çµ����¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in dimension of constituents

    real(DP) :: xyz_UTest    (0:imax-1, 1:jmax, 1:kmax)
                              ! �±è¥¿é¢���ï¼����¹ã�����
                              ! Zonal Wind (for test)  
    real(DP) :: xyz_VTest    (0:imax-1, 1:jmax, 1:kmax)
                              ! ����é¢���ï¼����¹ã�����
                              ! Meridional Wind (for test) 
    real(DP) :: xyr_SigDotTest(0:imax-1, 1:jmax, 0:kmax)
                              ! ���´æ���ï¼����¹ã�����;SigmaDot (for test) 
    real(DP) :: xyzf_QMixSave(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)

!!$    real(DP) :: xyrf_QMixA(0:imax-1, 1:jmax, 0:kmax, 1:ncmax)
!!$
!!$    integer :: k


    ! �»ã�����°ã���³ã�¸ã�¥æ���������³ªç§»æ�è¨�ç®�
    ! Semi-Lagrangian method for tracer transport      
!!$!      xyzf_QMixA = xyzf_QMixB !���¹ã����
!!$      xyzf_QMixA = xyzf_QMixB + xyzf_DQMixDtPhy * DelTime
    xyzf_QMixA = xyzf_QMixB + xyzf_DQMixDtPhy * 2.0_DP * DelTime


    ! Save a variable for mass fixer
    xyzf_QMixSave = xyzf_QMixA


    ! Mass fixer
    !   Constituents
    !
!!$!        call MassFixer(                  &
!!$    call MassFixerColumn(                  &
!!$!          & xyr_PressA,                  & ! (in)
!!$      & xyr_PressB,                  & ! (in)
!!$      & xyzf_QMixA,                  & ! (inout)
!!$      & xyr_PressRef = xyr_PressB,   & ! (in) optional
!!$!          & xyzf_QMixRef = ( xyzf_QMixB+xyzf_DQMixDtPhy*DelTime ) & ! (in) optional
!!$!      & xyzf_QMixRef = ( xyzf_QMixB+xyzf_DQMixDtPhy*2.0_DP*DelTime ) & ! (in) optional
!!$      & xyzf_QMixRef = xyzf_QMixSave & ! (in) optional
!!$      & )
    !
!!$      call MassFixer(                   &
      call MassFixerColumn( xyr_PressB, xyzf_QMixA, xyr_PressRef = xyr_PressB, xyzf_QMixRef = xyzf_QMixSave )


    ! Save a variable for mass fixer
    xyzf_QMixSave = xyzf_QMixA


    if (FlagSLTTArcsine) then
      ! ��è²���ä¿�証ã���������� arcsineå¤������£ã����
      ! Arcsine transformation for non-negative filter 

      do n = 1, ncmax
        f_QMixProcMax(n) = maxval( xyzf_QMixA(:,:,:,n) )
      end do
      call MPIWrapperFindMaxVal( ncmax, f_QMixProcMax, f_QMixMax )
      f_QMixMax = f_QMixMax * (1.05_DP) + 1.0e-14_DP
      do n = 1, ncmax
        xyzf_QMixA(:,:,:,n) = 0.5_DP*(asin(2.0_DP*xyzf_QMixA(:,:,:,n)/f_QMixMax(n) - 1.0_DP))
      end do
    end if

    xyzf_QMixA = SLTTHorAdv( xyzf_QMixA, xyz_UN, xyz_VN )           ! æ°´å¹³�»ã������
                                                                    ! Horizontal

!!$    if (FlagSLTTArcsine) then
!!$      ! ��è²���ä¿�証ã���������� arcsineå¤������£ã���¿ï���å¤���ï¼�
!!$      ! Arcsine transformation for non-negative filter
!!$      do n = 1, ncmax
!!$        xyzf_QMixA(:,:,:,n) = &
!!$          & f_QMixMax(n)*(0.5_DP)*(sin(2.0_DP*xyzf_QMixA(:,:,:,n))+1.0_DP) 
!!$      enddo
!!$    endif
!!$
!!$
!!$    call MassFixerLayer(             &
!!$      & xyr_PressA,                  & ! (in)
!!$      & xyzf_QMixA,                  & ! (inout)
!!$!      & xyr_PressRef = xyr_PressA,   & ! (in) optional
!!$      & xyr_PressRef = xyr_PressB,   & ! (in) optional
!!$      & xyzf_QMixRef = xyzf_QMixSave & ! (in) optional
!!$      & )
!!$
!!$
!!$    ! Save a variable for mass fixer
!!$    xyzf_QMixSave = xyzf_QMixA
!!$
!!$
!!$    if (FlagSLTTArcsine) then
!!$      ! ��è²���ä¿�証ã���������� arcsineå¤������£ã����
!!$      ! Arcsine transformation for non-negative filter
!!$
!!$      do n = 1, ncmax
!!$        f_QMixProcMax(n) = maxval( xyzf_QMixA(:,:,:,n) )
!!$      end do
!!$      call MPIWrapperFindMaxVal( &
!!$        & ncmax, f_QMixProcMax,  & ! (in)
!!$        & f_QMixMax              & ! (out)
!!$        & )
!!$      f_QMixMax = f_QMixMax * (1.05_DP) + 1.0e-14_DP
!!$      do n = 1, ncmax
!!$        xyzf_QMixA(:,:,:,n) = &
!!$          & 0.5_DP*(asin(2.0_DP*xyzf_QMixA(:,:,:,n)/f_QMixMax(n) - 1.0_DP))
!!$      end do
!!$    end if


    xyzf_QMixA = SLTTVerAdv( xyr_SigDotN, xyzf_QMixA )              ! ���´ã�»ã������
                                                                    ! Vertical 

    ! Vertical advection by finite difference method
    !
!!$    do n = 1, ncmax
!!$      k = 1
!!$      xyrf_QMixA(:,:,k,n) = 1.0e100_DP
!!$      do k = 1, kmax-1
!!$        xyrf_QMixA(:,:,k,n) = &
!!$          & ( xyzf_QMixA(:,:,k,n) + xyzf_QMixA(:,:,k+1,n) ) / 2.0_DP
!!$      end do
!!$      k = kmax
!!$      xyrf_QMixA(:,:,k,n) = 1.0e100_DP
!!$    end do
!!$    do n = 1, ncmax
!!$      do k = 1, kmax
!!$        xyzf_QMixA(:,:,k,n) = xyzf_QMixA(:,:,k,n)                     &
!!$          & + (                                                       &
!!$          &     - (   xyr_SigDotN(:,:,k-1) * xyrf_QMixA(:,:,k-1,n)    &
!!$          &         - xyr_SigDotN(:,:,k  ) * xyrf_QMixA(:,:,k  ,n) )  &
!!$          &       / z_DelSigma(k)                                     &
!!$          &     + xyzf_QMixA(:,:,k,n)                                 &
!!$          &       * ( xyr_SigDotN(:,:,k-1) - xyr_SigDotN(:,:,k  ) )   &
!!$          &       / z_DelSigma(k)                                     &
!!$          &   ) * 2.0_DP * DelTime
!!$      end do
!!$    end do


    ! 移æ����¹ã��
!    call SLTTTest(xyz_UTest, xyz_VTest, xyr_SigDotTest)
!    xyzf_QMixA = SLTTHorAdv( xyzf_QMixA, xyz_UTest, xyz_VTest )           ! æ°´å¹³�»ã������
!    xyzf_QMixA = SLTTVerAdv( xyr_SigDotTest, xyzf_QMixA )              ! ���´ã�»ã������

    if (FlagSLTTArcsine) then
      ! ��è²���ä¿�証ã���������� arcsineå¤������£ã���¿ï���å¤���ï¼�
      ! Arcsine transformation for non-negative filter
      do n = 1, ncmax
        xyzf_QMixA(:,:,:,n) = f_QMixMax(n)*(0.5_DP)*(sin(2.0_DP*xyzf_QMixA(:,:,:,n))+1.0_DP) 
      enddo
    endif


!!$!      xyzf_QMixA = xyzf_QMixB !���¹ã����
!!$      xyzf_QMixA = xyzf_QMixA + xyzf_DQMixDtPhy * DelTime


    ! Mass fixer
!!$    call MassFixerColumn(               &
!!$      & xyr_PressA,                     & ! (in)
!!$      & xyzf_QMixA,                     & ! (inout)
!!$!      & xyr_PressRef = xyr_PressB,  & ! (in) optional
!!$      & xyr_PressRef = xyr_PressA,      & ! (in) optional
!!$      & xyzf_QMixRef = xyzf_QMixSave    & ! (in) optional
!!$      & )
    call MassFixer( xyr_PressA, xyzf_QMixA, xyr_PressRef = xyr_PressB, xyzf_QMixRef = xyzf_QMixSave )


  end subroutine SLTTMain

Private Instance methods

FlagSLTTArcsine
Variable :
FlagSLTTArcsine :logical, save
: Arcsineå¤�������è²����£ã���¿ã������ Flag for non-negative filter using arcsine trasformation
SLTTHorAdv( xyzf_QMix, xyz_U, xyz_V ) result(xyzf_QMixA)
Function :
xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP)
: 次ã�¹ã����������³ªæ··å��æ¯� Next mix ratio of the tracers

local variables

xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in )
: �¾å�����»ã����³ªæ··å��æ¯� Present mix ratio of the tracers
xyz_U(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
: �±è¥¿é¢��� Zonal Wind
xyz_V(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
: ������� Meridional Wind

�»ã�����°ã���³ã�¸ã�¥æ�������水平移æ����ç®� Calculates tracer transports by Semi-Lagrangian method for horizontal direction

[Source]

  function SLTTHorAdv( xyzf_QMix, xyz_U, xyz_V ) result( xyzf_QMixA )
    ! �»ã�����°ã���³ã�¸ã�¥æ�������水平移æ����ç®�
    ! Calculates tracer transports by Semi-Lagrangian method for horizontal direction

    use timeset    , only : DelTime
                              ! $\Delta t$
    use axesset    , only : x_Lon, y_Lat
                              ! $\lambda, \varphai$ lon and lat
    use sltt_const , only : dtjw, iexmin, iexmax, jexmins, jexmaxs, jexminn, jexmaxn
    use sltt_extarr, only : SLTTExtArrExt, SLTTExtArrExt2
                              ! �����¡å¼µ���¼ã����
                              ! Expansion of arrays
    use sltt_dp    , only : SLTTDPHor
                              ! æ°´å¹³ä¸�æµ��¹æŽ¢ç´�
                              ! Finding departure point in horizontal
    use sltt_lagint, only : SLTTIrrHerIntK13
                              ! 水平�次������
                              ! 2D Interpolation in horizontal 

    ! 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_xya            => wa_xva, xya_wa            => xva_wa, wa_DLon_wa, xya_GradLat_wa => xva_GradLat_wa
#else
    use wa_mpi_module, only: wa_xya            => wa_xva, xya_wa            => xva_wa, wa_DLon_wa, xya_GradLat_wa => xva_GradLat_wa
#endif
#elif AXISYMMETRY
    use wa_zonal_module, only: wa_xya, xya_wa, wa_DLon_wa, xya_GradLat_wa
#elif SJPACK
    use wa_module_sjpack, only: wa_xya, xya_wa, wa_DLon_wa, xya_GradLat_wa 
#elif AXISYMMETRY_SJPACK
    use wa_zonal_module_sjpack, only: wa_xya, xya_wa, wa_DLon_wa, xya_GradLat_wa
#else
    use wa_module, only: wa_xya, xya_wa , wa_DLon_wa, xya_GradLat_wa
#endif


    real(DP), intent(in ) :: xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! �¾å�����»ã����³ªæ··å��æ¯�
                              ! Present mix ratio of the tracers
    real(DP), intent(in ) :: xyz_U    (0:imax-1, 1:jmax, 1:kmax)
                              ! �±è¥¿é¢���
                              ! Zonal Wind    
    real(DP), intent(in ) :: xyz_V    (0:imax-1, 1:jmax, 1:kmax)
                              ! �������
                              ! Meridional Wind    

    real(DP) :: xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! 次ã�¹ã����������³ªæ··å��æ¯�
                              ! Next mix ratio of the tracers
    !
    ! local variables
    !
    real(DP) :: xyzf_ExtQMixS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax)
                              ! �¾å�����»ã����³ªæ··å��æ¯����¡å¼µ����ï¼�������ï¼�
                              ! Extended array (SH) of present mix ratio of the tracers.
    real(DP) :: xyzf_ExtQMixN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax)
                              ! �¾å�����»ã����³ªæ··å��æ¯����¡å¼µ����ï¼�������ï¼�
                              ! Extended array (NH) of present mix ratio of the tracers.
    real(DP) :: xyz_ExtUS    (iexmin:iexmax, jexmins:jexmaxs, 1:kmax)
                              ! �±è¥¿é¢������¡å¼µ����ï¼�������ï¼�
                              ! Extended array (SH) of Zonal Wind        
    real(DP) :: xyz_ExtUN    (iexmin:iexmax, jexminn:jexmaxn, 1:kmax)
                              ! �±è¥¿é¢������¡å¼µ����ï¼�������ï¼�
                              ! Extended array (NH) of Zonal Wind        
    real(DP) :: xyz_ExtVS    (iexmin:iexmax, jexmins:jexmaxs, 1:kmax)
                              ! ����é¢������¡å¼µ����ï¼�������ï¼�
                              ! Extended array (SH) of Meridional Wind
    real(DP) :: xyz_ExtVN    (iexmin:iexmax, jexminn:jexmaxn, 1:kmax)
                              ! ����é¢������¡å¼µ����ï¼�������ï¼�
                              ! Extended array (NH) of Meridional Wind

    integer:: i, ii           ! �±è¥¿�¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in zonal direction
    integer:: j               ! �����¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in meridional direction
    integer:: k               ! ���´æ�¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in vertical direction
    integer:: n               ! çµ����¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in dimension of constituents

    real(DP) :: xyz_DPLonS(0:imax-1, 1:jmax/2, 1:kmax)
                              ! ä¸�æµ��¹ç�åº��������ï¼�
                              ! Lon of the departure point (SH)
    real(DP) :: xyz_DPLonN(0:imax-1, 1:jmax/2, 1:kmax)
                              ! ä¸�æµ��¹ç�åº��������ï¼�
                              ! Lon of the departure point (NH)    
    real(DP) :: xyz_DPLatS(0:imax-1, 1:jmax/2, 1:kmax)
                              ! ä¸�æµ��¹ç·¯åº��������ï¼�
                              ! Lat of the departure point (SH)    
    real(DP) :: xyz_DPLatN(0:imax-1, 1:jmax/2, 1:kmax)
                              ! ä¸�æµ��¹ç·¯åº��������ï¼�
                              ! Lat of the departure point (NH)    

    real(DP) :: xyzf_QMixAS(0:imax-1, 1:jmax/2, 1:kmax, 1:ncmax)
                              ! 次ã�¹ã����������³ªæ··å��æ¯�ï¼�������ï¼�
                              ! Next mix ratio of the tracers (SH)
    real(DP) :: xyzf_QMixAN(0:imax-1, 1:jmax/2, 1:kmax, 1:ncmax)
                              ! 次ã�¹ã����������³ªæ··å��æ¯�ï¼�������ï¼�
                              ! Next mix ratio of the tracers (NH)

!---fx, fy, fxy
    real(DP) :: xyzf_QMix_dlon(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯����åº�¾®��ï¼��°ã������ï¼�
                              ! Zonal derivative of the mix ratio (on grid)
    real(DP) :: xyzf_QMix_dlat(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯���·¯åº�¾®��ï¼��°ã������ï¼�
                              ! Meridional derivative of the mix ratio (on grid)
    real(DP) :: xyzf_QMix_dlonlat(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯���·¯åº��åº�¾®��ï¼��°ã������ï¼�
                              ! Zonal and meridional derivative of the mix ratio (on grid)    
    real(DP) :: xyzf_ExtQMixS_dlon(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯����åº�¾®�����¡å¼µ����ï¼�������ï¼�
                              ! Extended array (SH) of zonal derivative of the mix ratio
    real(DP) :: xyzf_ExtQMixN_dlon(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯����åº�¾®�����¡å¼µ����ï¼�������ï¼�
                              ! Extended array (NH) of zonal derivative of the mix ratio
    real(DP) :: xyzf_ExtQMixS_dlat(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯���·¯åº�¾®�����¡å¼µ����ï¼�������ï¼�
                              ! Extended array (SH) of meridional derivative of the mix ratio
    real(DP) :: xyzf_ExtQMixN_dlat(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax)                              
                              ! ��³ªæ··å��æ¯���·¯åº�¾®�����¡å¼µ����ï¼�������ï¼�
                              ! Extended array (NH) of meridional derivative of the mix ratio
    real(DP) :: xyzf_ExtQMixS_dlonlat(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯���·¯åº��åº�¾®�����¡å¼µ����ï¼�������ï¼�
                              ! Extended array (SH) of zonal and meridional derivative of the mix ratio
    real(DP) :: xyzf_ExtQMixN_dlonlat(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯���·¯åº��åº�¾®�����¡å¼µ����ï¼�������ï¼�
                              ! Extended array (NH) of zonal and meridional derivative of the mix ratio
    real(DP) :: wzf_QMix(1:lmax, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯����åº�¾®��ï¼��¹ã���������
                              ! Zonal derivative of the mix ratio (on grid)    
    real(DP) :: wzf_QMix_dlon(1:lmax, 1:kmax, 1:ncmax)        
                              ! ��³ªæ··å��æ¯����åº�¾®��ï¼��¹ã���������
                              ! Zonal derivative of the mix ratio (on grid)
    real(DP) :: pm            ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                              ! Sign change flag for array extension; -1.0 for sign change over the pole, 1.0 for no sign change

!---fxx, fyy, fxxyy
!    real(DP) :: xyzf_QMix_dlon2(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_QMix_dlat2(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_QMix_dlon2lat2(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixS_dlon2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixN_dlon2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixS_dlat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixN_dlat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixS_dlon2lat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixN_dlon2lat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!----fxxy
!    real(DP) :: xyzf_QMix_dlon2lat(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixS_dlon2lat(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixN_dlon2lat(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!----fxyy
!    real(DP) :: xyzf_QMix_dlonlat2(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixS_dlonlat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!    real(DP) :: xyzf_ExtQMixN_dlonlat2(-2+0:imax-1+3, -jew+1:jmax/2+jew, 1:kmax, 1:ncmax)
!----
!    real(DP) :: wzf_QMix_dlon2(1:lmax, 1:kmax, 1:ncmax)        



    ! ���� ; Executable statement
    !

    ! ������確è�
    ! Initialization check
    !
    if ( .not. sltt_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    ! QMix��¾®��è¨�ç®�ï¼��¹ã����������������
    ! Derivatives of QMix
    do n = 1, ncmax
        wzf_QMix(:,:,n) = wa_xya(xyzf_QMix(:,:,:,n))                     ! �°ã���������¹ã��������
                                                                         ! grid -> spectrum
        xyzf_QMix_dlat(:,:,:,n) = xya_GradLat_wa(wzf_QMix(:,:,n))        ! �¹ã�����������°ã������ç·�º¦å¾���
                                                                         ! spectrum -> grid (dQ/dlat)
        wzf_QMix_dlon(:,:,n) = wa_Dlon_wa(wzf_QMix(:,:,n))               ! �¹ã�����������¹ã���������åº�¾®��
                                                                         ! spectrum -> spectrum (dQ/dlon)        
        xyzf_QMix_dlon(:,:,:,n) = xya_wa(wzf_QMix_dlon(:,:,n))           ! �¹ã���������åº�¾®�����°ã������çµ�åº�¾®��
                                                                         ! spectrum (dQ/dlon) -> grid (dQ/dlon)
        xyzf_QMix_dlonlat(:,:,:,n) = xya_GradLat_wa(wzf_QMix_dlon(:,:,n))! �¹ã���������åº�¾®�����°ã������ç·�º¦çµ�åº�¾®��
                                                                         ! spectrum (dQ/dlon) -> grid (d^2Q/dlon dlat)        

        !---fxx, fyy, fxxy, fxyy, fxxyy ����
        !xyzf_QMix_dlon2(:,:,:,n) = xya_wa(wa_Dlon_wa(wzf_QMix_dlon(:,:,n)))
        !xyzf_QMix_dlat2(:,:,:,n) = xya_GradLat_wa(wa_xya(xyzf_QMix_dlat(:,:,:,n)))
        !xyzf_QMix_dlon2lat(:,:,:,n) = xya_GradLat_wa(wa_xya(xyzf_QMix_dlon2(:,:,:,n)))
        !xyzf_QMix_dlonlat2(:,:,:,n) = xya_GradLat_wa(wa_xya(xyzf_QMix_dlonlat(:,:,:,n)))
        !xyzf_QMix_dlon2lat2(:,:,:,n) = xya_GradLat_wa(wa_xya(xyzf_QMix_dlon2lat(:,:,:,n)))
    enddo


    ! ���������²ã���¡å¼µ
    ! Division and extension of arrays
    !
    ! ���������²ã���¡å¼µ
    ! Division and extension of arrays

    pm = -1.0_DP ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                 ! -1.0 if the sign of value changes over the poles; if not 1.0. 
!!$    call SLTTExtArrExt2(                             &
!!$      & xyzf_QMix_dlon,  pm,                         & ! (in)
!!$      & xyzf_ExtQMixS_dlon, xyzf_ExtQMixN_dlon       & ! (out)
!!$      & )
    call SLTTExtArrExt2( x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix_dlon,  pm, xyzf_ExtQMixS_dlon, xyzf_ExtQMixN_dlon, "Wave1" )

    pm = -1.0_DP ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                 ! -1.0 if the sign of value changes over the poles; if not 1.0.
!!$    call SLTTExtArrExt2(                             &
!!$      & xyzf_QMix_dlat,  pm,                         & ! (in)
!!$      & xyzf_ExtQMixS_dlat, xyzf_ExtQMixN_dlat       & ! (out)
!!$      & )
    call SLTTExtArrExt2( x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix_dlat,  pm, xyzf_ExtQMixS_dlat, xyzf_ExtQMixN_dlat, "Wave1" )

    pm = +1.0_DP ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                 ! -1.0 if the sign of value changes over the poles; if not 1.0.
!!$    call SLTTExtArrExt2(                             &
!!$      & xyzf_QMix_dlonlat, pm,                       & ! (in)
!!$      & xyzf_ExtQMixS_dlonlat, xyzf_ExtQMixN_dlonlat & ! (out)
!!$      & )
    call SLTTExtArrExt2( x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix_dlonlat, pm, xyzf_ExtQMixS_dlonlat, xyzf_ExtQMixN_dlonlat, "Wave1" )

!-----fxx, fyy, fxxy, fxyy, fxxyy �������¡å¼µ
!    pm = +1.0_DP ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                  ! -1.0 if the sign of value changes over the poles; if not 1.0. 
!    call SLTTExtArrExt2(                             &
!      & xyzf_QMix_dlon2,  pm,                        & ! (in)
!      & xyzf_ExtQMixS_dlon2, xyzf_ExtQMixN_dlon2     & ! (out)
!      & )
!    pm = +1.0_DP ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                  ! -1.0 if the sign of value changes over the poles; if not 1.0. 
!    call SLTTExtArrExt2(                             &
!      & xyzf_QMix_dlat2,  pm,                        & ! (in)
!      & xyzf_ExtQMixS_dlat2, xyzf_ExtQMixN_dlat2     & ! (out)
!      & )      
!    pm = -1.0_DP ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                  ! -1.0 if the sign of value changes over the poles; if not 1.0. 
!    call SLTTExtArrExt2(                               &
!      & xyzf_QMix_dlon2lat,  pm,                       & ! (in)
!      & xyzf_ExtQMixS_dlon2lat, xyzf_ExtQMixN_dlon2lat & ! (out)
!      & )      
!    pm = -1.0_DP ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                  ! -1.0 if the sign of value changes over the poles; if not 1.0. 
!    call SLTTExtArrExt2(                               &
!      & xyzf_QMix_dlonlat2,  pm,                       & ! (in)
!      & xyzf_ExtQMixS_dlonlat2, xyzf_ExtQMixN_dlonlat2 & ! (out)
!      & )
!    pm = +1.0_DP ! �����¡å¼µ��������極ã����å¾�����·ã��å¤�����´å���� -1.0��ä¸������������§ã�����´å����1.0��ä¸�������
                  ! -1.0 if the sign of value changes over the poles; if not 1.0. 
!    call SLTTExtArrExt2(                                 &
!      & xyzf_QMix_dlon2lat2,  pm,                        & ! (in)
!      & xyzf_ExtQMixS_dlon2lat2, xyzf_ExtQMixN_dlon2lat2 & ! (out)
!      & )


!!$    call SLTTExtArrExt(                             &
!!$      & x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, & ! (in)
!!$      & xyzf_QMix, xyz_U, xyz_V,                    & ! (in)
!!$      & xyzf_ExtQMixS, xyzf_ExtQMixN,               & ! (out)
!!$      & xyz_ExtUS, xyz_ExtUN,                       & ! (out)
!!$      & xyz_ExtVS, xyz_ExtVN                        & ! (out)
!!$      & )
    call SLTTExtArrExt( y_ExtLatS, y_ExtLatN, x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix, xyz_U, xyz_V, xyzf_ExtQMixS_dlat, xyzf_ExtQMixN_dlat, xyzf_ExtQMixS, xyzf_ExtQMixN, xyz_ExtUS, xyz_ExtUN, xyz_ExtVS, xyz_ExtVN )



    ! ä¸�æµ��¹ã���ç®�
    ! estimation of departure point
    ! ������
    ! south array
    call SLTTDPHor( DelTime, x_LonS, y_LatS, y_SinLatS, y_CosLatS, iexmin, iexmax, jexmins, jexmaxs, x_ExtLonS, y_ExtLatS, xyz_ExtUS, xyz_ExtVS, xyz_DPLonS, xyz_DPLatS )
    ! ������
    ! north array
    call SLTTDPHor( DelTime, x_LonN, y_LatN, y_SinLatN, y_CosLatN, iexmin, iexmax, jexminn, jexmaxn, x_ExtLonN, y_ExtLatN, xyz_ExtUN, xyz_ExtVN, xyz_DPLonN, xyz_DPLatN )



    ! ���
    ! Interpolation
!    do n = 1, ncmax
    call SLTTIrrHerIntK13( iexmin, iexmax, jexmins, jexmaxs, x_ExtLonS, y_ExtLatS, xyz_DPLonS, xyz_DPLatS, xyzf_ExtQMixS(:,:,:,:), xyzf_ExtQMixS_dlon(:,:,:,:), xyzf_ExtQMixS_dlat(:,:,:,:), xyzf_ExtQMixS_dlonlat(:,:,:,:), SLTTIntHor, xyzf_QMixAS(:,:,:,:) )

    call SLTTIrrHerIntK13( iexmin, iexmax, jexminn, jexmaxn, x_ExtLonN, y_ExtLatN, xyz_DPLonN, xyz_DPLatN, xyzf_ExtQMixN(:,:,:,:), xyzf_ExtQMixN_dlon(:,:,:,:), xyzf_ExtQMixN_dlat(:,:,:,:), xyzf_ExtQMixN_dlonlat(:,:,:,:), SLTTIntHor, xyzf_QMixAN(:,:,:,:) )
!    enddo

    ! �������������������
    ! joint of each array
     xyzf_QMixA(:,1:jmax/2,:,:) = xyzf_QMixAS(:,1:jmax/2,:,:)
     xyzf_QMixA(:,jmax/2+1:jmax,:,:) = xyzf_QMixAN(:,1:jmax/2,:,:)


  end function SLTTHorAdv
SLTTIntHor
Variable :
SLTTIntHor :character(TOKEN), save
: æ°´å¹³�¹å��������¹æ�����å®������­ã�¼ã���¼ã�� Keyword for Interpolation Method for Horizontal direction
SLTTIntVer
Variable :
SLTTIntVer :character(TOKEN), save
: ���´æ�¹å��������¹æ�����å®������­ã�¼ã���¼ã�� Keyword for Interpolation Method for Vertical direction
SLTTTest( xyz_UTest, xyz_VTest, xyr_SigDotTest )
Subroutine :
xyz_UTest(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: �±è¥¿é¢��� Zonal Wind
xyz_VTest(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
: ������� Meridional Wind
xyr_SigDotTest(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(out)
: ���´æ���ï¼�SigmaDotï¼�

——-�»ã�����°ã�����¹ã�����������å¸���ä¸�����——— Gives a velocity for Test. Only used for debug.

[Source]

  subroutine SLTTTest( xyz_UTest, xyz_VTest, xyr_SigDotTest )
    !-------�»ã�����°ã�����¹ã�����������å¸���ä¸�����--------
    !Gives a velocity for Test. Only used for debug.
    
    
    use constants0, only : PI
    use axesset   , only : x_Lon, y_Lat, r_Sigma
    use constants, only: RPlanet 
                              ! $ a $ [m]. 
                              ! �������. 
                              ! Radius of planet
    use timeset, only: DelTime, TimeN                 ! �¹ã������ $ t $ ������. Time of step $ t $. 

    real(DP), intent(out) :: xyz_UTest    (0:imax-1, 1:jmax, 1:kmax)
                              ! �±è¥¿é¢���
                              ! Zonal Wind    
    real(DP), intent(out) :: xyz_VTest    (0:imax-1, 1:jmax, 1:kmax)
                              ! �������
                              ! Meridional Wind    
    real(DP), intent(out) :: xyr_SigDotTest(0:imax-1, 1:jmax, 0:kmax)
                              ! ���´æ���ï¼�SigmaDotï¼�

    ! �業��
    ! Work variables
    !
    real(DP) :: u0, t, shape
    real(DP), parameter :: lat0 = 0.5_8*PI, lon0 = 0.0_8*PI
    real(DP), parameter :: tau = 345600.0_DP, p0 = 100000.0_DP
    real(DP), parameter :: omega0 = PI*30000.0_DP/tau

    integer:: i               ! �±è¥¿�¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in zonal direction
    integer:: j               ! �����¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in meridional direction
    integer:: k               ! ���´æ�¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in vertical direction



    ! 水平�������������
    u0 = 2.0_DP*PI*RPlanet/(86400.0_DP*12.0_DP)
    do k = 1, kmax
      do j=1, jmax
        do i=0, imax-1
          xyz_UTest(i,j,k) = u0*(cos(y_Lat(j)) * cos(lat0) + cos(x_Lon(i)) * sin(y_Lat(j)) * sin(lat0))
          xyz_VTest(i,j,k) = -u0*(sin(x_Lon(i))*sin(lat0))
        end do
      end do
    end do
    !���´æ�����å¸���ä¸�����
    t = TimeN
    do k = 0, kmax
      shape = min(1.0_DP, 0.5_DP*(sin((r_Sigma(k)-r_Sigma(kmax))/(1.0_DP - r_Sigma(kmax))*PI)) )
        do j = 1, jmax
          do i = 0, imax-1
            xyr_SigDotTest(i,j,k) = -omega0/p0*cos(2.0_DP*PI/tau*t)*sin(shape*PI*0.5_DP)
          enddo
        enddo
    enddo

    end subroutine SLTTTest
SLTTVerAdv( xyr_SigmaDot, xyzf_QMix ) result(xyzf_QMixA)
Function :
xyzf_QMixA(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP)
: 次ã�¹ã����������³ªæ··å��æ¯� Next mix ratio of the tracers
xyr_SigmaDot(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in)
: ���´æ���ï¼�SigmaDotï¼�
xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) :real(DP), intent(in)
: �¾å�����»ã����³ªæ··å��æ¯� Present mix ratio of the tracers

�»ã�����°ã���³ã�¸ã�¥æ����������´ç§»æµ����ç®� Calculates tracer transports by Semi-Lagrangian method for vertical direction

[Source]

  function SLTTVerAdv( xyr_SigmaDot, xyzf_QMix ) result( xyzf_QMixA )
    ! �»ã�����°ã���³ã�¸ã�¥æ����������´ç§»æµ����ç®�
    ! Calculates tracer transports by Semi-Lagrangian method for vertical direction

    use axesset, only : z_Sigma           ! ���´åº§æ¨�; Sigma coordinate
    use timeset, only : DelTime           ! $\Delta t$
    use sltt_dp, only : SLTTDPVer         ! ���´ä�æµ��¹æŽ¢ç´�; Finding departure point in vertical 
    use sltt_lagint, only : SLTTIrrHerIntQui1DNonUni, SLTTHerIntCub1D                                  

    real(DP), intent(in) :: xyr_SigmaDot(0:imax-1, 1:jmax, 0:kmax)
                              ! ���´æ���ï¼�SigmaDotï¼�
    real(DP), intent(in) :: xyzf_QMix   (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! �¾å�����»ã����³ªæ··å��æ¯�
                              ! Present mix ratio of the tracers
    real(DP)             :: xyzf_QMixA  (0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! 次ã�¹ã����������³ªæ··å��æ¯�
                              ! Next mix ratio of the tracers    

    !
    ! local variables
    !
    real(DP) :: xyz_DPSigma(0:imax-1, 1:jmax, 1:kmax)
                              ! ä¸�æµ��¹é�åº�
                              ! Sigma of the departure point
    integer:: i               ! �±è¥¿�¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in zonal direction
    integer:: j               ! �����¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in meridional direction
    integer:: k, kk           ! ���´æ�¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in vertical direction
    integer:: n               ! çµ����¹å�������� DO ���¼ã�����業å���
                              ! Work variables for DO loop in dimension of constituents
    integer:: xy_kk(0:imax-1, 1:jmax)
                              ! ä¸�æµ��¹ã���ä¸����°ã���������¢ç´¢�����������業å���
                              ! Work variable for finding the grid just above the departure point
                              
    
    real(DP) :: xyzf_QMix_dz(0:imax-1, 1:jmax, 1:kmax, 1:ncmax)
                              ! ��³ªæ··å��æ¯������´å¾®��
                              ! Vertical derivative of the mix ratio    
    real(DP) :: xyzf_ExtQMix(0:imax-1, 1:jmax, 1-2:kmax+2, 1:ncmax)
                              ! ��³ªæ··å��æ¯����¡å¼µ����
                              ! Extended array of the mix ratio    
    real(DP) :: z_ExtSigma(1-2:kmax+2)
                              ! �座æ����¡å¼µ����
                              ! Extended array of the sigma coordinate            
    real(DP) :: xyf_F11(0:imax-1, 1:jmax, 1:ncmax)
                              ! ����������������業��
                              ! work variable for the derivative calculation
    real(DP) :: xyf_F22(0:imax-1, 1:jmax, 1:ncmax)
                              ! ����������������業��
                              ! work variable for the derivative calculation
    real(DP) :: xyf_F12(0:imax-1, 1:jmax, 1:ncmax)
                              ! ����������������業��
                              ! work variable for the derivative calculation
    real(DP) :: xyf_F21(0:imax-1, 1:jmax, 1:ncmax)
                              ! ����������������業��
                              ! work variable for the derivative calculation
    real(DP) :: s1, t1, s2, t2, r1, r2
                              ! ����������������業��
                              ! work variable for the derivative calculation
    





    ! ���� ; Executable statement
    !

    ! ������確è�
    ! Initialization check
    !
    if ( .not. sltt_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if

    ! ä¸�æµ��¹æŽ¢ç´�
    ! estimation of departure point
    !
    call SLTTDPVer( DelTime, xyr_SigmaDot, xyz_DPSigma )

     
     ! �����¡å¼µï¼�z_Sigmaï¼�
     ! Array extension for z_Sigma
     z_ExtSigma(-1) = 2.0_DP - z_Sigma(2)
     z_ExtSigma(0) = 2.0_DP - z_Sigma(1)     
     z_ExtSigma(1:kmax) = z_Sigma(1:kmax)
     z_ExtSigma(kmax+1) = -z_Sigma(kmax)     
     z_ExtSigma(kmax+2) = -z_Sigma(kmax-1)     
          
     ! �����¡å¼µï¼�xyzf_QMixï¼�
     ! Array extension for Q_Mix
     xyzf_ExtQMix(:,:,-1,:) = xyzf_QMix(:,:,2,:)
     xyzf_ExtQMix(:,:,0,:) = xyzf_QMix(:,:,1,:)
     xyzf_ExtQMix(:,:,1:kmax,:) = xyzf_QMix(:,:,1:kmax,:)
     xyzf_ExtQMix(:,:,kmax+1,:) = xyzf_QMix(:,:,kmax,:)
     xyzf_ExtQMix(:,:,kmax+2,:) = xyzf_QMix(:,:,kmax-1,:)

     ! xyzf_QMix_dz������������ 
     ! calculate xyzf_QMix_dz
     do k = 1 , kmax
        s1 = z_ExtSigma(k) - z_ExtSigma(k-1)
        t1 = z_ExtSigma(k+1) - z_ExtSigma(k)
        s2 = z_ExtSigma(k) - z_ExtSigma(k-2)
        t2 = z_ExtSigma(k+2) - z_ExtSigma(k)
        
        if (s1 == t1 .and. s2 == t2 .and. s1 + s1 == s2) then 
          ! �¼å���ç­��������´å��
          ! Uniform depth
          ! 4次精�
          ! 4th order

          xyzf_QMix_dz(:,:,k,:) = ( 8.0_DP*( xyzf_ExtQMix(:,:,k+1,:) - xyzf_ExtQMix(:,:,k-1,:)) - ( xyzf_ExtQMix(:,:,k+2,:) - xyzf_ExtQMix(:,:,k-2,:) ) )/12.0_DP
        else
          ! �¼å���ä¸�å½��������´å��
          ! Non-uniform depth
          xyf_F11 = (s1*s1*xyzf_ExtQMix(:,:,k+1,:) +(t1*t1 - s1*s1)*xyzf_ExtQMix(:,:,k,:) - t1*t1*xyzf_ExtQMix(:,:,k-1,:)) /(s1*t1*(s1+t1))
          xyf_F22 = (s2*s2*xyzf_ExtQMix(:,:,k+2,:) +(t2*t2 - s2*s2)*xyzf_ExtQMix(:,:,k,:) - t2*t2*xyzf_ExtQMix(:,:,k-2,:)) /(s2*t2*(s2+t2))
          xyf_F21 = (s2*s2*xyzf_ExtQMix(:,:,k+1,:) +(t1*t1 - s2*s2)*xyzf_ExtQMix(:,:,k,:) - t1*t1*xyzf_ExtQMix(:,:,k-2,:)) /(s2*t1*(s2+t1))
          xyf_F12 = (s1*s1*xyzf_ExtQMix(:,:,k+2,:) +(t2*t2 - s1*s1)*xyzf_ExtQMix(:,:,k,:) - t2*t2*xyzf_ExtQMix(:,:,k-1,:)) /(s1*t2*(s1+t2))
          
          r1 = t1 - s1 - t2 + s2
          r2 = t1 - s2 - t2 + s1
          !�次精�
          ! 4th order
          xyzf_QMix_dz(:,:,k,:) = ( (xyf_F11*s2*t2 - xyf_F22*s1*t1)*r2 - (xyf_F21*s1*t2 - xyf_F12*s2*t1)*r1 ) / ( (s2*t2-s1*t1)*r2 - (s1*t2-s2*t1)*r1 )
  
          !3次精�
          ! 3rd order
  !        xyzf_QMix_dz(:,:,k,:) = (xyf_F11*s2*t2 - xyf_F22(:,:,:)*s1*t1)/(s2*t2 - s1*t1) 
          
          !2次精�
          ! 2nd order
  !        xyzf_QMix_dz(:,:,k,:) = xyf_F11
        endif





    enddo


    xy_kk = 2
    do k = 1, kmax
    do j = 1, jmax
    do i = 0, imax-1
        if ( xyz_DPSigma(i,j,k) >= z_Sigma(1) ) then     ! DP�� z_Sigma(1) �� �°è¡¨��(sigma = 1.0)�������´å��
                                                         ! if DP is between z_Sigma(1) and the ground (sigma = 1.0)
            xyzf_QMixA(i,j,k,:) = xyzf_QMix(i,j,1,:)     ! Q_1�§ä�å®���������
                                                         ! use Q_1 for interpolated value
                                                                 
        elseif (xyz_DPSigma(i,j,k) <= z_Sigma(kmax)) then! DP�� z_Sigma(kmax) �� 大æ�ä¸�ç«�(sigma = 0.0)����
                                                         ! if DP is between z_Sigma(kmax) and the upper boundary (sigma = 0.0)        
            xyzf_QMixA(i,j,k,:) = xyzf_QMix(i,j,kmax,:)  ! Q_kmax�§ä�å®���������
                                                         ! use Q_kmax for interpolated value            
        else
            do kk = xy_kk(i,j), kmax 
                if ( xyz_DPSigma(i,j,k) > z_Sigma(kk) ) then 
                  select case (SLTTIntVer)
                    case("HQ")    ! å¤����������¼ã��ï¼�次è���; Irregular Hermite Quintic interpolation
                      do n = 1, ncmax 
                          xyzf_QMixA(i,j,k,n) = SLTTIrrHerIntQui1DNonUni(xyzf_ExtQMix(i,j,kk-2,n), xyzf_ExtQMix(i,j,kk-1,n), xyzf_ExtQMix(i,j,kk,n),   xyzf_ExtQMix(i,j,kk+1,n), xyzf_QMix_dz(i,j,kk-1,n), xyzf_QMix_dz(i,j,kk,n), z_ExtSigma(kk-2)-z_ExtSigma(kk-1), z_ExtSigma(kk)-z_ExtSigma(kk-1), z_ExtSigma(kk+1)-z_ExtSigma(kk-1), xyz_DPSigma(i,j,k)-z_ExtSigma(kk-1))
                      enddo

                    case("HC")    ! �������¼ã��ï¼�次è���; Hermitian Cubic interpolation
                      do n = 1, ncmax 
                          xyzf_QMixA(i,j,k,n) = SLTTHerIntCub1D( xyzf_ExtQMix(i,j,kk-1,n), xyzf_ExtQMix(i,j,kk,n), xyzf_QMix_dz(i,j,kk-1,n), xyzf_QMix_dz(i,j,kk,n), z_ExtSigma(kk)-z_ExtSigma(kk-1), xyz_DPSigma(i,j,k)-z_ExtSigma(kk-1))
                      enddo

                    case default
                      write( 6, * ) "ERROR : GIVE CORRECT KEYWORD FOR <SLTTIntVer> IN NAMELIST"
                      stop
                  end select
                    xy_kk(i,j) = kk
                    exit    
                endif
            enddo 
        endif
    enddo
    enddo
    enddo

  end function SLTTVerAdv
module_name
Constant :
module_name = ‘sltt‘ :character(*), parameter
: �¢ã�¸ã�¥ã�¼ã������ç§�. Module name
sltt_inited
Variable :
sltt_inited = .false. :logical, save
: ����設������. Initialization flag
version
Constant :
version = ’$Name: $’ // ’$Id: sltt.F90,v 1.8 2014/06/29 07:21:28 yot Exp $’ :character(*), parameter
: �¢ã�¸ã�¥ã�¼ã�������¼ã�¸ã�§ã�� Module version
x_CosLonN
Variable :
x_CosLonN(:) :real(DP) , save, allocatable
: $cos\lambda_N$
x_CosLonS
Variable :
x_CosLonS(:) :real(DP) , save, allocatable
: $cos\lambda_S$
x_ExtLonN
Variable :
x_ExtLonN(:) :real(DP) , save, allocatable
: $ x_LonN���¡å¼µ������ Extended array of x_LonN.
x_ExtLonS
Variable :
x_ExtLonS(:) :real(DP) , save, allocatable
: $ x_LonS���¡å¼µ������ Extended array of x_LonS.
x_LonN
Variable :
x_LonN(:) :real(DP) , save, allocatable
: $lambda_N$ ������������ longitude in NH.
x_LonS
Variable :
x_LonS(:) :real(DP) , save, allocatable
: $lambda_S$ ������������ longitude in SH.
x_SinLonN
Variable :
x_SinLonN(:) :real(DP) , save, allocatable
: $sin\lambda_N$
x_SinLonS
Variable :
x_SinLonS(:) :real(DP) , save, allocatable
: $sin\lambda_S$
y_CosLatN
Variable :
y_CosLatN(:) :real(DP) , save, allocatable
: $cos\varphai_N$
y_CosLatS
Variable :
y_CosLatS(:) :real(DP) , save, allocatable
: $cos\varphai_S$
y_ExtLatN
Variable :
y_ExtLatN(:) :real(DP) , save, allocatable
: $ x_LatN���¡å¼µ������ Extended array of x_LatN.
y_ExtLatS
Variable :
y_ExtLatS(:) :real(DP) , save, allocatable
: $ x_LatS���¡å¼µ������ Extended array of x_LatS.
y_LatN
Variable :
y_LatN(:) :real(DP) , save, allocatable
: $varphi_N$ ��������·¯åº��� latitude in NH.
y_LatS
Variable :
y_LatS(:) :real(DP) , save, allocatable
: $varphi_S$ ��������·¯åº��� latitude in SH.
y_SinLatN
Variable :
y_SinLatN(:) :real(DP) , save, allocatable
: $sin\varphai_N$
y_SinLatS
Variable :
y_SinLatS(:) :real(DP) , save, allocatable
: $sin\varphai_S$