SCALE-RM
Functions/Subroutines
scale_atmos_dyn_tinteg_tracer_rk3 Module Reference

module Atmosphere / Dyn Tinteg More...

Functions/Subroutines

subroutine, public atmos_dyn_tinteg_tracer_rk3_setup (tinteg_type)
 Setup. More...
 
subroutine, public atmos_dyn_tinteg_tracer_rk3 (QTRC, qflx, QTRC0, RHOQ_t, DENS0, DENS, mflx_hi, num_diff, GSQRT, MAPF, CDZ, RCDZ, RCDX, RCDY, BND_W, BND_E, BND_S, BND_N, TwoD, dtl, FLAG_FCT_TRACER, FLAG_FCT_ALONG_STREAM)
 RK3. More...
 

Detailed Description

module Atmosphere / Dyn Tinteg

Description
Temporal integration for tracer advection for Atmospheric process three stage Runge-Kutta scheme
Author
Team SCALE

Function/Subroutine Documentation

◆ atmos_dyn_tinteg_tracer_rk3_setup()

subroutine, public scale_atmos_dyn_tinteg_tracer_rk3::atmos_dyn_tinteg_tracer_rk3_setup ( character(len=*)  tinteg_type)

Setup.

Definition at line 65 of file scale_atmos_dyn_tinteg_tracer_rk3.F90.

65  use scale_const, only: &
66  undef => const_undef
67  use scale_prc, only: &
68  prc_abort
69  use scale_comm_cartesc, only: &
71  implicit none
72 
73  character(len=*) :: tinteg_type
74 
75  integer :: iv
76  !---------------------------------------------------------------------------
77 
78  if ( tinteg_type /= 'RK3WS2002' ) then
79  log_error("ATMOS_DYN_Tinteg_tracer_rk3_setup",*) 'TINTEG_TRACER_TYPE is not RK3WS2002. Check!'
80  call prc_abort
81  end if
82 
83  allocate( qtrc_rk1(ka,ia,ja) )
84  allocate( qtrc_rk2(ka,ia,ja) )
85  qtrc_rk1(:,:,:) = undef
86  qtrc_rk2(:,:,:) = undef
87 
88  call comm_vars8_init( 'QTRC_RK1', qtrc_rk1, i_comm_rk1 )
89  call comm_vars8_init( 'QTRC_RK2', qtrc_rk2, i_comm_rk2 )
90 
91  return

References scale_comm_cartesc::comm_vars8_init(), scale_const::const_undef, scale_atmos_grid_cartesc_index::ia, scale_atmos_grid_cartesc_index::ja, scale_atmos_grid_cartesc_index::ka, and scale_prc::prc_abort().

Referenced by scale_atmos_dyn_tinteg_tracer::atmos_dyn_tinteg_tracer_setup().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_tinteg_tracer_rk3()

subroutine, public scale_atmos_dyn_tinteg_tracer_rk3::atmos_dyn_tinteg_tracer_rk3 ( real(rp), dimension (ka,ia,ja), intent(inout)  QTRC,
real(rp), dimension (ka,ia,ja,3), intent(out)  qflx,
real(rp), dimension (ka,ia,ja), intent(in)  QTRC0,
real(rp), dimension (ka,ia,ja), intent(in)  RHOQ_t,
real(rp), dimension (ka,ia,ja), intent(in)  DENS0,
real(rp), dimension (ka,ia,ja), intent(in)  DENS,
real(rp), dimension (ka,ia,ja,3), intent(in)  mflx_hi,
real(rp), dimension(ka,ia,ja,3), intent(in)  num_diff,
real(rp), dimension (ka,ia,ja,7), intent(in)  GSQRT,
real(rp), dimension (ia,ja), intent(in)  MAPF,
real(rp), dimension(ka), intent(in)  CDZ,
real(rp), dimension(ka), intent(in)  RCDZ,
real(rp), dimension(ia), intent(in)  RCDX,
real(rp), dimension(ja), intent(in)  RCDY,
logical, intent(in)  BND_W,
logical, intent(in)  BND_E,
logical, intent(in)  BND_S,
logical, intent(in)  BND_N,
logical, intent(in)  TwoD,
real(rp), intent(in)  dtl,
logical, intent(in)  FLAG_FCT_TRACER,
logical, intent(in)  FLAG_FCT_ALONG_STREAM 
)

RK3.

Definition at line 109 of file scale_atmos_dyn_tinteg_tracer_rk3.F90.

109  use scale_comm_cartesc, only: &
110  comm_vars8, &
111  comm_wait
112  use scale_atmos_dyn_tstep_tracer, only: &
114  use scale_atmos_dyn_common, only: &
116  implicit none
117  real(RP), intent(inout) :: QTRC (KA,IA,JA)
118  real(RP), intent(out) :: qflx (KA,IA,JA,3)
119  real(RP), intent(in) :: QTRC0 (KA,IA,JA)
120  real(RP), intent(in) :: RHOQ_t (KA,IA,JA)
121  real(RP), intent(in) :: DENS0 (KA,IA,JA)
122  real(RP), intent(in) :: DENS (KA,IA,JA)
123  real(RP), intent(in) :: mflx_hi (KA,IA,JA,3)
124  real(RP), intent(in) :: num_diff(KA,IA,JA,3)
125  real(RP), intent(in) :: GSQRT (KA,IA,JA,7)
126  real(RP), intent(in) :: MAPF (IA,JA)
127  real(RP), intent(in) :: CDZ(KA)
128  real(RP), intent(in) :: RCDZ(KA)
129  real(RP), intent(in) :: RCDX(IA)
130  real(RP), intent(in) :: RCDY(JA)
131  logical, intent(in) :: BND_W
132  logical, intent(in) :: BND_E
133  logical, intent(in) :: BND_S
134  logical, intent(in) :: BND_N
135  logical, intent(in) :: TwoD
136  real(RP), intent(in) :: dtl
137  logical, intent(in) :: FLAG_FCT_TRACER
138  logical, intent(in) :: FLAG_FCT_ALONG_STREAM
139 
140  real(RP) :: DENS_RK(KA,IA,JA)
141  real(RP) :: dtrk
142  integer :: k, i, j
143 
144  !------------------------------------------------------------------------
145  ! Start RK
146  !------------------------------------------------------------------------
147 
148  do j = js-1, je+1
149  do i = max(is-1,1), min(ie+1,ia)
150  do k = ks, ke
151  dens_rk(k,i,j) = dens0(k,i,j) &
152  + ( dens(k,i,j) - dens0(k,i,j) ) / 3.0_rp
153  end do
154  end do
155  end do
156 
157  dtrk = dtl / 3.0_rp
158  call atmos_dyn_tstep_tracer( &
159  qtrc_rk1, & ! (out)
160  qflx, & ! (out)
161  qtrc, qtrc0, rhoq_t, &! (in)
162  dens0, dens_rk, & ! (in)
163  mflx_hi, num_diff, & ! (in)
164  gsqrt, mapf, & ! (in)
165  cdz, rcdz, rcdx, rcdy, & ! (in)
166  twod, dtrk, & ! (in)
167  .false., flag_fct_along_stream ) ! (in)
168 
169  call atmos_dyn_copy_boundary_tracer( qtrc_rk1, & ! [INOUT]
170  qtrc0, & ! [IN]
171  bnd_w, bnd_e, bnd_s, bnd_n, & ! [IN]
172  twod ) ! [IN]
173 
174  call comm_vars8( qtrc_rk1(:,:,:), i_comm_rk1 )
175  call comm_wait ( qtrc_rk1(:,:,:), i_comm_rk1, .false. )
176 
177 
178  do j = js-1, je+1
179  do i = max(is-1,1), min(ie+1,ia)
180  do k = ks, ke
181  dens_rk(k,i,j) = dens0(k,i,j) &
182  + ( dens(k,i,j) - dens0(k,i,j) ) * 0.5_rp
183  end do
184  end do
185  end do
186 
187  dtrk = dtl / 2.0_rp
188  call atmos_dyn_tstep_tracer( &
189  qtrc_rk2, & ! (out)
190  qflx, & ! (out)
191  qtrc_rk1, qtrc0, rhoq_t, &! (in)
192  dens0, dens_rk, & ! (in)
193  mflx_hi, num_diff, & ! (in)
194  gsqrt, mapf, & ! (in)
195  cdz, rcdz, rcdx, rcdy, & ! (in)
196  twod, dtrk, & ! (in)
197  .false., flag_fct_along_stream ) ! (in)
198 
199  call atmos_dyn_copy_boundary_tracer( qtrc_rk2, & ! [INOUT]
200  qtrc0, & ! [IN]
201  bnd_w, bnd_e, bnd_s, bnd_n, & ! [IN]
202  twod ) ! [IN]
203 
204  call comm_vars8( qtrc_rk2(:,:,:), i_comm_rk2 )
205  call comm_wait ( qtrc_rk2(:,:,:), i_comm_rk2, .false. )
206 
207 
208  dtrk = dtl
209  call atmos_dyn_tstep_tracer( &
210  qtrc, & ! (out)
211  qflx, & ! (out)
212  qtrc_rk2, qtrc0, rhoq_t, &! (in)
213  dens0, dens, & ! (in)
214  mflx_hi, num_diff, & ! (in)
215  gsqrt, mapf, & ! (in)
216  cdz, rcdz, rcdx, rcdy, & ! (in)
217  twod, & ! (in)
218  dtrk, & ! (in)
219  flag_fct_tracer, flag_fct_along_stream ) ! (in)
220 
221  return

References scale_atmos_dyn_common::atmos_dyn_copy_boundary_tracer(), scale_atmos_dyn_tstep_tracer::atmos_dyn_tstep_tracer, scale_atmos_grid_cartesc_index::ia, scale_atmos_grid_cartesc_index::ie, scale_atmos_grid_cartesc_index::is, scale_atmos_grid_cartesc_index::je, scale_atmos_grid_cartesc_index::js, scale_tracer::k, scale_atmos_grid_cartesc_index::ke, and scale_atmos_grid_cartesc_index::ks.

Referenced by scale_atmos_dyn_tinteg_tracer::atmos_dyn_tinteg_tracer_setup().

Here is the call graph for this function:
Here is the caller graph for this function:
scale_atmos_grid_cartesc_index::ke
integer, public ke
end point of inner domain: z, local
Definition: scale_atmos_grid_cartesC_index.F90:52
scale_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:342
scale_atmos_dyn_common
module Atmosphere / Dynamics common
Definition: scale_atmos_dyn_common.F90:12
scale_atmos_grid_cartesc_index::ka
integer, public ka
Definition: scale_atmos_grid_cartesC_index.F90:47
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_atmos_grid_cartesc_index::ie
integer, public ie
end point of inner domain: x, local
Definition: scale_atmos_grid_cartesC_index.F90:54
scale_const
module CONSTANT
Definition: scale_const.F90:11
scale_atmos_grid_cartesc_index::ia
integer, public ia
Definition: scale_atmos_grid_cartesC_index.F90:48
scale_comm_cartesc::comm_vars8_init
subroutine, public comm_vars8_init(varname, var, vid)
Register variables.
Definition: scale_comm_cartesC.F90:294
scale_atmos_dyn_common::atmos_dyn_copy_boundary_tracer
subroutine, public atmos_dyn_copy_boundary_tracer(QTRC, QTRC0, BND_W, BND_E, BND_S, BND_N, TwoD)
Definition: scale_atmos_dyn_common.F90:311
scale_atmos_dyn_tstep_tracer
module Atmosphere / Dynamical scheme
Definition: scale_atmos_dyn_tstep_tracer.F90:12
scale_atmos_grid_cartesc_index::is
integer, public is
start point of inner domain: x, local
Definition: scale_atmos_grid_cartesC_index.F90:53
scale_atmos_grid_cartesc_index::ja
integer, public ja
Definition: scale_atmos_grid_cartesC_index.F90:49
scale_atmos_grid_cartesc_index::ks
integer, public ks
start point of inner domain: z, local
Definition: scale_atmos_grid_cartesC_index.F90:51
scale_comm_cartesc
module COMMUNICATION
Definition: scale_comm_cartesC.F90:11
scale_atmos_dyn_tstep_tracer::atmos_dyn_tstep_tracer
procedure(step), pointer, public atmos_dyn_tstep_tracer
Definition: scale_atmos_dyn_tstep_tracer.F90:72
scale_atmos_grid_cartesc_index::js
integer, public js
start point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:55
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:41
scale_atmos_grid_cartesc_index::je
integer, public je
end point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:56