SCALE-RM
Functions/Subroutines
scale_atmos_dyn_tinteg_large_rk3 Module Reference

module Atmosphere / Dyn Tinteg More...

Functions/Subroutines

subroutine, public atmos_dyn_tinteg_large_rk3_setup (tinteg_type)
 Setup. More...
 
subroutine, public atmos_dyn_tinteg_large_rk3 (DENS, MOMZ, MOMX, MOMY, RHOT, QTRC, PROG, DENS_av, MOMZ_av, MOMX_av, MOMY_av, RHOT_av, QTRC_av, mflx_hi, tflx_hi, num_diff, num_diff_q, DENS_tp, MOMZ_tp, MOMX_tp, MOMY_tp, RHOT_tp, RHOQ_tp, CORIOLI, CDZ, CDX, CDY, FDZ, FDX, FDY, RCDZ, RCDX, RCDY, RFDZ, RFDX, RFDY, PHI, GSQRT, J13G, J23G, J33G, MAPF, AQ_R, AQ_CV, AQ_CP, AQ_MASS, REF_dens, REF_pott, REF_qv, REF_pres, BND_W, BND_E, BND_S, BND_N, ND_COEF, ND_COEF_Q, ND_ORDER, ND_SFC_FACT, ND_USE_RS, BND_QA, BND_SMOOTHER_FACT, DAMP_DENS, DAMP_VELZ, DAMP_VELX, DAMP_VELY, DAMP_POTT, DAMP_QTRC, DAMP_alpha_DENS, DAMP_alpha_VELZ, DAMP_alpha_VELX, DAMP_alpha_VELY, DAMP_alpha_POTT, DAMP_alpha_QTRC, wdamp_coef, divdmp_coef, FLAG_TRACER_SPLIT_TEND, FLAG_FCT_MOMENTUM, FLAG_FCT_T, FLAG_FCT_TRACER, FLAG_FCT_ALONG_STREAM, USE_AVERAGE, I_QV, DTL, DTS)
 RK3. More...
 

Detailed Description

module Atmosphere / Dyn Tinteg

Description
Temporal integration in Dynamical core for Atmospheric process three step Runge-Kutta scheme
Author
Team SCALE

Function/Subroutine Documentation

◆ atmos_dyn_tinteg_large_rk3_setup()

subroutine, public scale_atmos_dyn_tinteg_large_rk3::atmos_dyn_tinteg_large_rk3_setup ( character(len=*)  tinteg_type)

Setup.

Definition at line 61 of file scale_atmos_dyn_tinteg_large_rk3.F90.

References scale_prc::prc_abort().

Referenced by scale_atmos_dyn_tinteg_large::atmos_dyn_tinteg_large_setup().

61  use scale_prc, only: &
62  prc_abort
63  implicit none
64 
65  character(len=*) :: tinteg_type
66 
67  integer :: iv
68  !---------------------------------------------------------------------------
69 
70  if ( tinteg_type /= 'RK3' ) then
71  log_error("ATMOS_DYN_Tinteg_large_rk3_setup",*) 'TINTEG_LARGE_TYPE is not RK3. Check!'
72  call prc_abort
73  end if
74 
75  return
module PROCESS
Definition: scale_prc.F90:11
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:338
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_tinteg_large_rk3()

subroutine, public scale_atmos_dyn_tinteg_large_rk3::atmos_dyn_tinteg_large_rk3 ( real(rp), dimension(ka,ia,ja), intent(inout)  DENS,
real(rp), dimension(ka,ia,ja), intent(inout)  MOMZ,
real(rp), dimension(ka,ia,ja), intent(inout)  MOMX,
real(rp), dimension(ka,ia,ja), intent(inout)  MOMY,
real(rp), dimension(ka,ia,ja), intent(inout)  RHOT,
real(rp), dimension(ka,ia,ja,qa), intent(inout)  QTRC,
real(rp), dimension(ka,ia,ja,va), intent(inout)  PROG,
real(rp), dimension(ka,ia,ja), intent(inout)  DENS_av,
real(rp), dimension(ka,ia,ja), intent(inout)  MOMZ_av,
real(rp), dimension(ka,ia,ja), intent(inout)  MOMX_av,
real(rp), dimension(ka,ia,ja), intent(inout)  MOMY_av,
real(rp), dimension(ka,ia,ja), intent(inout)  RHOT_av,
real(rp), dimension(ka,ia,ja,qa), intent(inout)  QTRC_av,
real(rp), dimension(ka,ia,ja,3), intent(out)  mflx_hi,
real(rp), dimension(ka,ia,ja,3), intent(out)  tflx_hi,
real(rp), dimension(ka,ia,ja,5,3), intent(out)  num_diff,
real(rp), dimension(ka,ia,ja,3), intent(out)  num_diff_q,
real(rp), dimension(ka,ia,ja), intent(in)  DENS_tp,
real(rp), dimension(ka,ia,ja), intent(in)  MOMZ_tp,
real(rp), dimension(ka,ia,ja), intent(in)  MOMX_tp,
real(rp), dimension(ka,ia,ja), intent(in)  MOMY_tp,
real(rp), dimension(ka,ia,ja), intent(in)  RHOT_tp,
real(rp), dimension(ka,ia,ja,qa), intent(in)  RHOQ_tp,
real(rp), dimension(ia,ja), intent(in)  CORIOLI,
real(rp), dimension (ka), intent(in)  CDZ,
real(rp), dimension (ia), intent(in)  CDX,
real(rp), dimension (ja), intent(in)  CDY,
real(rp), dimension (ka-1), intent(in)  FDZ,
real(rp), dimension (ia-1), intent(in)  FDX,
real(rp), dimension (ja-1), intent(in)  FDY,
real(rp), dimension(ka), intent(in)  RCDZ,
real(rp), dimension(ia), intent(in)  RCDX,
real(rp), dimension(ja), intent(in)  RCDY,
real(rp), dimension(ka-1), intent(in)  RFDZ,
real(rp), dimension(ia-1), intent(in)  RFDX,
real(rp), dimension(ja-1), intent(in)  RFDY,
real(rp), dimension (ka,ia,ja), intent(in)  PHI,
real(rp), dimension(ka,ia,ja,7), intent(in)  GSQRT,
real(rp), dimension (ka,ia,ja,7), intent(in)  J13G,
real(rp), dimension (ka,ia,ja,7), intent(in)  J23G,
real(rp), intent(in)  J33G,
real(rp), dimension (ia,ja,2,4), intent(in)  MAPF,
real(rp), dimension (qa), intent(in)  AQ_R,
real(rp), dimension (qa), intent(in)  AQ_CV,
real(rp), dimension (qa), intent(in)  AQ_CP,
real(rp), dimension(qa), intent(in)  AQ_MASS,
real(rp), dimension(ka,ia,ja), intent(in)  REF_dens,
real(rp), dimension(ka,ia,ja), intent(in)  REF_pott,
real(rp), dimension (ka,ia,ja), intent(in)  REF_qv,
real(rp), dimension(ka,ia,ja), intent(in)  REF_pres,
logical, intent(in)  BND_W,
logical, intent(in)  BND_E,
logical, intent(in)  BND_S,
logical, intent(in)  BND_N,
real(rp), intent(in)  ND_COEF,
real(rp), intent(in)  ND_COEF_Q,
integer, intent(in)  ND_ORDER,
real(rp), intent(in)  ND_SFC_FACT,
logical, intent(in)  ND_USE_RS,
integer, intent(in)  BND_QA,
real(rp), intent(in)  BND_SMOOTHER_FACT,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_DENS,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_VELZ,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_VELX,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_VELY,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_POTT,
real(rp), dimension(ka,ia,ja,bnd_qa), intent(in)  DAMP_QTRC,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_alpha_DENS,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_alpha_VELZ,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_alpha_VELX,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_alpha_VELY,
real(rp), dimension(ka,ia,ja), intent(in)  DAMP_alpha_POTT,
real(rp), dimension(ka,ia,ja,bnd_qa), intent(in)  DAMP_alpha_QTRC,
real(rp), dimension(ka), intent(in)  wdamp_coef,
real(rp), intent(in)  divdmp_coef,
logical, intent(in)  FLAG_TRACER_SPLIT_TEND,
logical, intent(in)  FLAG_FCT_MOMENTUM,
logical, intent(in)  FLAG_FCT_T,
logical, intent(in)  FLAG_FCT_TRACER,
logical, intent(in)  FLAG_FCT_ALONG_STREAM,
logical, intent(in)  USE_AVERAGE,
integer, intent(in)  I_QV,
real(dp), intent(in)  DTL,
real(dp), intent(in)  DTS 
)

RK3.

Parameters
[in]phigeopotential
[in]gsqrtvertical metrics {G}^1/2
[in]j13g(1,3) element of Jacobian matrix
[in]j23g(2,3) element of Jacobian matrix
[in]j33g(3,3) element of Jacobian matrix
[in]mapfmap factor
[in]ref_presreference pressure

Definition at line 109 of file scale_atmos_dyn_tinteg_large_rk3.F90.

References scale_atmos_dyn_common::atmos_dyn_fct(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_xyz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_xyz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_xyz_ud1(), scale_atmos_dyn_common::atmos_dyn_numfilter_coef(), scale_atmos_dyn_common::atmos_dyn_numfilter_coef_q(), scale_atmos_dyn_tstep_large::atmos_dyn_tstep_large, scale_const::const_cvdry, scale_const::const_rdry, scale_const::const_rvap, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), and scale_index::va.

Referenced by scale_atmos_dyn_tinteg_large::atmos_dyn_tinteg_large_setup().

109  use scale_const, only: &
110  rdry => const_rdry, &
111  rvap => const_rvap, &
112  cvdry => const_cvdry
113  use scale_comm_cartesc, only: &
114  comm_vars8, &
115  comm_wait
116  use scale_atmos_dyn_common, only: &
120  use scale_atmos_dyn_fvm_flux_ud1, only: &
124  use scale_atmos_dyn_fvm_flux, only: &
128  use scale_atmos_dyn_tstep_large, only: &
130 #ifdef HIST_TEND
131  use scale_file_history, only: &
132  file_history_in
133 #endif
134  implicit none
135 
136  real(RP), intent(inout) :: dens(ka,ia,ja)
137  real(RP), intent(inout) :: momz(ka,ia,ja)
138  real(RP), intent(inout) :: momx(ka,ia,ja)
139  real(RP), intent(inout) :: momy(ka,ia,ja)
140  real(RP), intent(inout) :: rhot(ka,ia,ja)
141  real(RP), intent(inout) :: qtrc(ka,ia,ja,qa)
142 
143  real(RP), intent(inout) :: prog(ka,ia,ja,va)
144 
145  real(RP), intent(inout) :: dens_av(ka,ia,ja)
146  real(RP), intent(inout) :: momz_av(ka,ia,ja)
147  real(RP), intent(inout) :: momx_av(ka,ia,ja)
148  real(RP), intent(inout) :: momy_av(ka,ia,ja)
149  real(RP), intent(inout) :: rhot_av(ka,ia,ja)
150  real(RP), intent(inout) :: qtrc_av(ka,ia,ja,qa)
151 
152  real(RP), intent(out) :: mflx_hi(ka,ia,ja,3)
153  real(RP), intent(out) :: tflx_hi(ka,ia,ja,3)
154 
155  real(RP), intent(out) :: num_diff(ka,ia,ja,5,3)
156  real(RP), intent(out) :: num_diff_q(ka,ia,ja,3)
157 
158  real(RP), intent(in) :: dens_tp(ka,ia,ja)
159  real(RP), intent(in) :: momz_tp(ka,ia,ja)
160  real(RP), intent(in) :: momx_tp(ka,ia,ja)
161  real(RP), intent(in) :: momy_tp(ka,ia,ja)
162  real(RP), intent(in) :: rhot_tp(ka,ia,ja)
163  real(RP), intent(in) :: rhoq_tp(ka,ia,ja,qa)
164 
165  real(RP), intent(in) :: corioli(ia,ja)
166 
167  real(RP), intent(in) :: cdz (ka)
168  real(RP), intent(in) :: cdx (ia)
169  real(RP), intent(in) :: cdy (ja)
170  real(RP), intent(in) :: fdz (ka-1)
171  real(RP), intent(in) :: fdx (ia-1)
172  real(RP), intent(in) :: fdy (ja-1)
173  real(RP), intent(in) :: rcdz(ka)
174  real(RP), intent(in) :: rcdx(ia)
175  real(RP), intent(in) :: rcdy(ja)
176  real(RP), intent(in) :: rfdz(ka-1)
177  real(RP), intent(in) :: rfdx(ia-1)
178  real(RP), intent(in) :: rfdy(ja-1)
179 
180  real(RP), intent(in) :: phi (ka,ia,ja)
181  real(RP), intent(in) :: gsqrt(ka,ia,ja,7)
182  real(RP), intent(in) :: j13g (ka,ia,ja,7)
183  real(RP), intent(in) :: j23g (ka,ia,ja,7)
184  real(RP), intent(in) :: j33g
185  real(RP), intent(in) :: mapf (ia,ja,2,4)
186 
187  real(RP), intent(in) :: aq_r (qa)
188  real(RP), intent(in) :: aq_cv (qa)
189  real(RP), intent(in) :: aq_cp (qa)
190  real(RP), intent(in) :: aq_mass(qa)
191 
192  real(RP), intent(in) :: ref_dens(ka,ia,ja)
193  real(RP), intent(in) :: ref_pott(ka,ia,ja)
194  real(RP), intent(in) :: ref_qv (ka,ia,ja)
195  real(RP), intent(in) :: ref_pres(ka,ia,ja)
196 
197  logical, intent(in) :: bnd_w
198  logical, intent(in) :: bnd_e
199  logical, intent(in) :: bnd_s
200  logical, intent(in) :: bnd_n
201 
202  real(RP), intent(in) :: nd_coef
203  real(RP), intent(in) :: nd_coef_q
204  integer, intent(in) :: nd_order
205  real(RP), intent(in) :: nd_sfc_fact
206  logical, intent(in) :: nd_use_rs
207 
208  integer, intent(in) :: bnd_qa
209  real(RP), intent(in) :: bnd_smoother_fact
210 
211  real(RP), intent(in) :: damp_dens(ka,ia,ja)
212  real(RP), intent(in) :: damp_velz(ka,ia,ja)
213  real(RP), intent(in) :: damp_velx(ka,ia,ja)
214  real(RP), intent(in) :: damp_vely(ka,ia,ja)
215  real(RP), intent(in) :: damp_pott(ka,ia,ja)
216  real(RP), intent(in) :: damp_qtrc(ka,ia,ja,bnd_qa)
217 
218  real(RP), intent(in) :: damp_alpha_dens(ka,ia,ja)
219  real(RP), intent(in) :: damp_alpha_velz(ka,ia,ja)
220  real(RP), intent(in) :: damp_alpha_velx(ka,ia,ja)
221  real(RP), intent(in) :: damp_alpha_vely(ka,ia,ja)
222  real(RP), intent(in) :: damp_alpha_pott(ka,ia,ja)
223  real(RP), intent(in) :: damp_alpha_qtrc(ka,ia,ja,bnd_qa)
224 
225  real(RP), intent(in) :: wdamp_coef(ka)
226  real(RP), intent(in) :: divdmp_coef
227 
228  logical, intent(in) :: flag_tracer_split_tend
229  logical, intent(in) :: flag_fct_momentum
230  logical, intent(in) :: flag_fct_t
231  logical, intent(in) :: flag_fct_tracer
232  logical, intent(in) :: flag_fct_along_stream
233 
234  logical, intent(in) :: use_average
235 
236  integer, intent(in) :: i_qv
237 
238  real(DP), intent(in) :: dtl
239  real(DP), intent(in) :: dts
240 
241  ! for time integartion
242  real(RP) :: dens0 (ka,ia,ja)
243  real(RP) :: momz0 (ka,ia,ja)
244  real(RP) :: momx0 (ka,ia,ja)
245  real(RP) :: momy0 (ka,ia,ja)
246  real(RP) :: rhot0 (ka,ia,ja)
247  real(RP) :: qtrc0 (ka,ia,ja,qa)
248  real(RP) :: prog0 (ka,ia,ja,va)
249 
250  real(DP) :: dtrk
251 
252  logical :: last
253 
254  integer :: n, iq
255 
256  !---------------------------------------------------------------------------
257 
258  call prof_rapstart("DYN_Large_RK3_Prep", 3)
259 
260  !##### SAVE #####
261 !OCL XFILL
262  dens0(:,:,:) = dens(:,:,:)
263 !OCL XFILL
264  momz0(:,:,:) = momz(:,:,:)
265 !OCL XFILL
266  momx0(:,:,:) = momx(:,:,:)
267 !OCL XFILL
268  momy0(:,:,:) = momy(:,:,:)
269 !OCL XFILL
270  rhot0(:,:,:) = rhot(:,:,:)
271 !OCL XFILL
272  qtrc0(:,:,:,:) = qtrc(:,:,:,:)
273 
274  if ( va > 0 ) then
275 !OCL XFILL
276  prog0(:,:,:,:) = prog(:,:,:,:)
277  end if
278 
279  call prof_rapend ("DYN_Large_RK3_Prep", 3)
280 
281  !------------------------------------------------------------------------
282  ! Start RK
283  !------------------------------------------------------------------------
284 
285  call prof_rapstart("DYN_Large_RK3", 3)
286 
287  do n = 1, 3
288 
289  dtrk = dtl / real( 4-n, kind=dp )
290  last = n == 3
291 
292  if ( n > 1 ) then
293  dens = dens0
294  momz = momz0
295  momx = momx0
296  momy = momy0
297  rhot = rhot0
298  if ( va > 0 ) prog = prog0
299  end if
300 
301  call atmos_dyn_tstep_large( &
302  dens, momz, momx, momy, rhot, qtrc, prog, & ! (inout)
303  dens_av, momz_av, momx_av, momy_av, rhot_av, qtrc_av, & ! (inout)
304  mflx_hi, tflx_hi, & ! (out)
305  num_diff, num_diff_q, & ! (out)
306  qtrc0, & ! (in)
307  dens_tp, momz_tp, momx_tp, momy_tp, rhot_tp, rhoq_tp, & ! (in)
308  corioli, & ! (in)
309  cdz, cdx, cdy, fdz, fdx, fdy, & ! (in)
310  rcdz, rcdx, rcdy, rfdz, rfdx, rfdy, & ! (in)
311  phi, gsqrt, & ! (in)
312  j13g, j23g, j33g, mapf, & ! (in)
313  aq_r, aq_cv, aq_cp, aq_mass, & ! (in)
314  ref_dens, ref_pott, ref_qv, ref_pres, & ! (in)
315  bnd_w, bnd_e, bnd_s, bnd_n, & ! (in)
316  nd_coef, nd_coef_q, nd_order, nd_sfc_fact, nd_use_rs, & ! (in)
317  bnd_qa, bnd_smoother_fact, & ! (in)
318  damp_dens, damp_velz, damp_velx, & ! (in)
319  damp_vely, damp_pott, damp_qtrc, & ! (in)
320  damp_alpha_dens, damp_alpha_velz, damp_alpha_velx, & ! (in)
321  damp_alpha_vely, damp_alpha_pott, damp_alpha_qtrc, & ! (in)
322  wdamp_coef, & ! (in)
323  divdmp_coef, & ! (in)
324  flag_tracer_split_tend, & ! (in)
325  flag_fct_momentum, flag_fct_t, flag_fct_tracer, & ! (in)
326  flag_fct_along_stream, & ! (in)
327  use_average .AND. last, & ! (in)
328  i_qv, & ! (in)
329  dtrk, dts, last ) ! (in)
330 
331  end do
332 
333  call prof_rapend ("DYN_Large_RK3", 3)
334 
335  return
real(rp), public const_cvdry
specific heat (dry air,constant volume) [J/kg/K]
Definition: scale_const.F90:57
subroutine, public atmos_dyn_fvm_fluxy_xyz_ud1(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation Y-flux at XYZ
integer, public va
Definition: scale_index.F90:35
subroutine, public atmos_dyn_fvm_fluxz_xyz_ud1(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation z-flux at XYZ
subroutine, public atmos_dyn_numfilter_coef(num_diff, DENS, MOMZ, MOMX, MOMY, RHOT, CDZ, CDX, CDY, FDZ, FDX, FDY, DT, REF_dens, REF_pott, ND_COEF, ND_ORDER, ND_SFC_FACT, ND_USE_RS)
Calc coefficient of numerical filter.
integer, public ia
of whole cells: x, local, with HALO
procedure(flux_phi), pointer, public atmos_dyn_fvm_fluxx_xyz
integer, public ja
of whole cells: y, local, with HALO
module Atmosphere / Dynamical scheme
real(rp), public const_rdry
specific gas constant (dry air) [J/kg/K]
Definition: scale_const.F90:55
module COMMUNICATION
subroutine, public atmos_dyn_numfilter_coef_q(num_diff_q, DENS, QTRC, is_qv, CDZ, CDX, CDY, dt, REF_qv, iq, ND_COEF, ND_ORDER, ND_SFC_FACT, ND_USE_RS)
Calc coefficient of numerical filter.
module Atmosphere / Dynamics common
real(rp), parameter, public const_rvap
specific gas constant (water vapor) [J/kg/K]
Definition: scale_const.F90:63
module CONSTANT
Definition: scale_const.F90:11
module scale_atmos_dyn_fvm_flux
integer, public ka
of whole cells: z, local, with HALO
subroutine, public atmos_dyn_fct(qflx_anti, phi_in, DENS0, DENS, qflx_hi, qflx_lo, mflx_hi, rdz, rdx, rdy, GSQRT, MAPF, dt, flag_vect)
Flux Correction Transport Limiter.
procedure(flux_phi), pointer, public atmos_dyn_fvm_fluxz_xyz
integer, parameter, public dp
procedure(flux_phi), pointer, public atmos_dyn_fvm_fluxy_xyz
procedure(large), pointer, public atmos_dyn_tstep_large
module file_history
subroutine, public atmos_dyn_fvm_fluxx_xyz_ud1(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation X-flux at XYZ
module scale_atmos_dyn_fvm_flux_ud1
Here is the call graph for this function:
Here is the caller graph for this function: