SCALE-RM
Functions/Subroutines
scale_atmos_dyn_tinteg_short_rk4 Module Reference

module Atmosphere / Dyn Tinteg More...

Functions/Subroutines

subroutine, public atmos_dyn_tinteg_short_rk4_setup (tinteg_type)
 Setup. More...
 
subroutine, public atmos_dyn_tinteg_short_rk4 (DENS, MOMZ, MOMX, MOMY, RHOT, PROG, mflx_hi, tflx_hi, DENS_t, MOMZ_t, MOMX_t, MOMY_t, RHOT_t, Rtot, CVtot, CORIOLI, num_diff, divdmp_coef, DDIV, FLAG_FCT_MOMENTUM, FLAG_FCT_T, FLAG_FCT_ALONG_STREAM, CDZ, FDZ, FDX, FDY, RCDZ, RCDX, RCDY, RFDZ, RFDX, RFDY, PHI, GSQRT, J13G, J23G, J33G, MAPF, REF_pres, REF_dens, BND_W, BND_E, BND_S, BND_N, dt)
 RK3. More...
 

Detailed Description

module Atmosphere / Dyn Tinteg

Description
Temporal integration in Dynamical core for Atmospheric process four step Runge-Kutta scheme
Author
Team SCALE
History
  • 2016-04-30 (S.Nishizawa) [new]

Function/Subroutine Documentation

◆ atmos_dyn_tinteg_short_rk4_setup()

subroutine, public scale_atmos_dyn_tinteg_short_rk4::atmos_dyn_tinteg_short_rk4_setup ( character(len=*)  tinteg_type)

Setup.

Definition at line 103 of file scale_atmos_dyn_tinteg_short_rk4.F90.

References scale_comm::comm_vars8_init(), scale_const::const_undef, scale_grid_index::ia, scale_grid_index::ja, scale_grid_index::ka, scale_process::prc_mpistop(), and scale_index::va.

Referenced by scale_atmos_dyn_tinteg_short::atmos_dyn_tinteg_short_setup().

103  use scale_process, only: &
105  use scale_const, only: &
106  undef => const_undef
107  use scale_comm, only: &
109  implicit none
110 
111  character(len=*) :: tinteg_type
112 
113  integer :: iv
114  !---------------------------------------------------------------------------
115 
116  if ( tinteg_type /= 'RK4' ) then
117  write(*,*) 'xxx TINTEG_TYPE is not RK4. Check!'
118  call prc_mpistop
119  end if
120 
121  allocate( dens_rk1(ka,ia,ja) )
122  allocate( momz_rk1(ka,ia,ja) )
123  allocate( momx_rk1(ka,ia,ja) )
124  allocate( momy_rk1(ka,ia,ja) )
125  allocate( rhot_rk1(ka,ia,ja) )
126 
127  allocate( dens_rk2(ka,ia,ja) )
128  allocate( momz_rk2(ka,ia,ja) )
129  allocate( momx_rk2(ka,ia,ja) )
130  allocate( momy_rk2(ka,ia,ja) )
131  allocate( rhot_rk2(ka,ia,ja) )
132 
133  allocate( dens_rk3(ka,ia,ja) )
134  allocate( momz_rk3(ka,ia,ja) )
135  allocate( momx_rk3(ka,ia,ja) )
136  allocate( momy_rk3(ka,ia,ja) )
137  allocate( rhot_rk3(ka,ia,ja) )
138 
139  allocate( prog_rk1(ka,ia,ja,max(va,1)) )
140  allocate( prog_rk2(ka,ia,ja,max(va,1)) )
141  allocate( prog_rk3(ka,ia,ja,max(va,1)) )
142  allocate( i_comm_prog_rk1(max(va,1)) )
143  allocate( i_comm_prog_rk2(max(va,1)) )
144  allocate( i_comm_prog_rk3(max(va,1)) )
145 
146  call comm_vars8_init( dens_rk1, i_comm_dens_rk1 )
147  call comm_vars8_init( momz_rk1, i_comm_momz_rk1 )
148  call comm_vars8_init( momx_rk1, i_comm_momx_rk1 )
149  call comm_vars8_init( momy_rk1, i_comm_momy_rk1 )
150  call comm_vars8_init( rhot_rk1, i_comm_rhot_rk1 )
151  do iv = 1, va
152  i_comm_prog_rk1(iv) = 5 + iv
153  call comm_vars8_init( prog_rk1(:,:,:,iv), i_comm_prog_rk1(iv) )
154  enddo
155 
156  call comm_vars8_init( dens_rk2, i_comm_dens_rk2 )
157  call comm_vars8_init( momz_rk2, i_comm_momz_rk2 )
158  call comm_vars8_init( momx_rk2, i_comm_momx_rk2 )
159  call comm_vars8_init( momy_rk2, i_comm_momy_rk2 )
160  call comm_vars8_init( rhot_rk2, i_comm_rhot_rk2 )
161  do iv = 1, va
162  i_comm_prog_rk2(iv) = 5 + iv
163  call comm_vars8_init( prog_rk2(:,:,:,iv), i_comm_prog_rk2(iv) )
164  enddo
165 
166  call comm_vars8_init( dens_rk3, i_comm_dens_rk3 )
167  call comm_vars8_init( momz_rk3, i_comm_momz_rk3 )
168  call comm_vars8_init( momx_rk3, i_comm_momx_rk3 )
169  call comm_vars8_init( momy_rk3, i_comm_momy_rk3 )
170  call comm_vars8_init( rhot_rk3, i_comm_rhot_rk3 )
171  do iv = 1, va
172  i_comm_prog_rk3(iv) = 5 + iv
173  call comm_vars8_init( prog_rk3(:,:,:,iv), i_comm_prog_rk3(iv) )
174  enddo
175 
176  dens_rk1(:,:,:) = undef
177  momz_rk1(:,:,:) = undef
178  momx_rk1(:,:,:) = undef
179  momy_rk1(:,:,:) = undef
180  rhot_rk1(:,:,:) = undef
181  if ( va > 0 ) prog_rk1(:,:,:,:) = undef
182 
183  dens_rk2(:,:,:) = undef
184  momz_rk2(:,:,:) = undef
185  momx_rk2(:,:,:) = undef
186  momy_rk2(:,:,:) = undef
187  rhot_rk2(:,:,:) = undef
188  if ( va > 0 ) prog_rk2(:,:,:,:) = undef
189 
190  dens_rk3(:,:,:) = undef
191  momz_rk3(:,:,:) = undef
192  momx_rk3(:,:,:) = undef
193  momy_rk3(:,:,:) = undef
194  rhot_rk3(:,:,:) = undef
195  if ( va > 0 ) prog_rk3(:,:,:,:) = undef
196 
197  return
subroutine, public prc_mpistop
Abort MPI.
integer, public va
Definition: scale_index.F90:38
real(rp), public const_undef
Definition: scale_const.F90:43
integer, public ia
of x whole cells (local, with HALO)
subroutine, public comm_vars8_init(var, vid)
Register variables.
Definition: scale_comm.F90:328
integer, public ka
of z whole cells (local, with HALO)
module COMMUNICATION
Definition: scale_comm.F90:23
module PROCESS
module CONSTANT
Definition: scale_const.F90:14
integer, public ja
of y whole cells (local, with HALO)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_tinteg_short_rk4()

subroutine, public scale_atmos_dyn_tinteg_short_rk4::atmos_dyn_tinteg_short_rk4 ( 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,va), intent(inout)  PROG,
real(rp), dimension(ka,ia,ja,3), intent(inout)  mflx_hi,
real(rp), dimension(ka,ia,ja,3), intent(inout)  tflx_hi,
real(rp), dimension(ka,ia,ja), intent(in)  DENS_t,
real(rp), dimension(ka,ia,ja), intent(in)  MOMZ_t,
real(rp), dimension(ka,ia,ja), intent(in)  MOMX_t,
real(rp), dimension(ka,ia,ja), intent(in)  MOMY_t,
real(rp), dimension(ka,ia,ja), intent(in)  RHOT_t,
real(rp), dimension(ka,ia,ja), intent(in)  Rtot,
real(rp), dimension(ka,ia,ja), intent(in)  CVtot,
real(rp), dimension(ia,ja), intent(in)  CORIOLI,
real(rp), dimension(ka,ia,ja,5,3), intent(in)  num_diff,
real(rp), intent(in)  divdmp_coef,
real(rp), dimension(ka,ia,ja), intent(in)  DDIV,
logical, intent(in)  FLAG_FCT_MOMENTUM,
logical, intent(in)  FLAG_FCT_T,
logical, intent(in)  FLAG_FCT_ALONG_STREAM,
real(rp), dimension (ka), intent(in)  CDZ,
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(ka,ia,ja), intent(in)  REF_pres,
real(rp), dimension(ka,ia,ja), intent(in)  REF_dens,
logical, intent(in)  BND_W,
logical, intent(in)  BND_E,
logical, intent(in)  BND_S,
logical, intent(in)  BND_N,
real(rp), intent(in)  dt 
)

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 216 of file scale_atmos_dyn_tinteg_short_rk4.F90.

References scale_atmos_dyn_common::atmos_dyn_copy_boundary(), scale_atmos_dyn_tstep_short::atmos_dyn_tstep_short, scale_grid_index::ie, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ke, scale_grid_index::ks, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), and scale_index::va.

Referenced by scale_atmos_dyn_tinteg_short::atmos_dyn_tinteg_short_setup().

216  use scale_comm, only: &
217  comm_vars8, &
218  comm_wait
219  use scale_atmos_dyn_tstep_short, only: &
220  atmos_dyn_tstep => atmos_dyn_tstep_short
221  use scale_atmos_dyn_common, only: &
223  implicit none
224 
225  real(RP), intent(inout) :: dens(ka,ia,ja)
226  real(RP), intent(inout) :: momz(ka,ia,ja)
227  real(RP), intent(inout) :: momx(ka,ia,ja)
228  real(RP), intent(inout) :: momy(ka,ia,ja)
229  real(RP), intent(inout) :: rhot(ka,ia,ja)
230  real(RP), intent(inout) :: prog(ka,ia,ja,va)
231 
232  real(RP), intent(inout) :: mflx_hi(ka,ia,ja,3)
233  real(RP), intent(inout) :: tflx_hi(ka,ia,ja,3)
234 
235  real(RP), intent(in) :: dens_t(ka,ia,ja)
236  real(RP), intent(in) :: momz_t(ka,ia,ja)
237  real(RP), intent(in) :: momx_t(ka,ia,ja)
238  real(RP), intent(in) :: momy_t(ka,ia,ja)
239  real(RP), intent(in) :: rhot_t(ka,ia,ja)
240 
241  real(RP), intent(in) :: rtot(ka,ia,ja)
242  real(RP), intent(in) :: cvtot(ka,ia,ja)
243  real(RP), intent(in) :: corioli(ia,ja)
244  real(RP), intent(in) :: num_diff(ka,ia,ja,5,3)
245  real(RP), intent(in) :: divdmp_coef
246  real(RP), intent(in) :: ddiv(ka,ia,ja)
247 
248  logical, intent(in) :: flag_fct_momentum
249  logical, intent(in) :: flag_fct_t
250  logical, intent(in) :: flag_fct_along_stream
251 
252  real(RP), intent(in) :: cdz (ka)
253  real(RP), intent(in) :: fdz (ka-1)
254  real(RP), intent(in) :: fdx (ia-1)
255  real(RP), intent(in) :: fdy (ja-1)
256  real(RP), intent(in) :: rcdz(ka)
257  real(RP), intent(in) :: rcdx(ia)
258  real(RP), intent(in) :: rcdy(ja)
259  real(RP), intent(in) :: rfdz(ka-1)
260  real(RP), intent(in) :: rfdx(ia-1)
261  real(RP), intent(in) :: rfdy(ja-1)
262 
263  real(RP), intent(in) :: phi (ka,ia,ja)
264  real(RP), intent(in) :: gsqrt(ka,ia,ja,7)
265  real(RP), intent(in) :: j13g (ka,ia,ja,7)
266  real(RP), intent(in) :: j23g (ka,ia,ja,7)
267  real(RP), intent(in) :: j33g
268  real(RP), intent(in) :: mapf (ia,ja,2,4)
269 
270  real(RP), intent(in) :: ref_pres(ka,ia,ja)
271  real(RP), intent(in) :: ref_dens(ka,ia,ja)
272 
273  logical, intent(in) :: bnd_w
274  logical, intent(in) :: bnd_e
275  logical, intent(in) :: bnd_s
276  logical, intent(in) :: bnd_n
277 
278  real(RP), intent(in) :: dt
279 
280  real(RP) :: dens0(ka,ia,ja)
281  real(RP) :: momz0(ka,ia,ja)
282  real(RP) :: momx0(ka,ia,ja)
283  real(RP) :: momy0(ka,ia,ja)
284  real(RP) :: rhot0(ka,ia,ja)
285  real(RP) :: prog0(ka,ia,ja,va)
286 
287  real(RP) :: mflx_hi_rk(ka,ia,ja,3,3)
288  real(RP) :: tflx_hi_rk(ka,ia,ja,3,3)
289 
290  real(RP) :: dtrk
291 
292  integer :: i, j, k, iv, n
293  !---------------------------------------------------------------------------
294 
295  call prof_rapstart("DYN_RK4_Prep",3)
296 
297 #ifdef DEBUG
298  dens_rk1(:,:,:) = undef
299  momz_rk1(:,:,:) = undef
300  momx_rk1(:,:,:) = undef
301  momy_rk1(:,:,:) = undef
302  rhot_rk1(:,:,:) = undef
303  if ( va > 0 ) prog_rk1(:,:,:,:) = undef
304 
305  dens_rk2(:,:,:) = undef
306  momz_rk2(:,:,:) = undef
307  momx_rk2(:,:,:) = undef
308  momy_rk2(:,:,:) = undef
309  rhot_rk2(:,:,:) = undef
310  if ( va > 0 ) prog_rk2(:,:,:,:) = undef
311 
312  dens_rk3(:,:,:) = undef
313  momz_rk3(:,:,:) = undef
314  momx_rk3(:,:,:) = undef
315  momy_rk3(:,:,:) = undef
316  rhot_rk3(:,:,:) = undef
317  if ( va > 0 ) prog_rk3(:,:,:,:) = undef
318 
319  mflx_hi(:,:,:,:) = undef
320  tflx_hi(:,:,:,:) = undef
321 
322  mflx_hi_rk(:,:,:,:,:) = undef
323  tflx_hi_rk(:,:,:,:,:) = undef
324 #endif
325 
326 !OCL XFILL
327  dens0 = dens
328 !OCL XFILL
329  momz0 = momz
330 !OCL XFILL
331  momx0 = momx
332 !OCL XFILL
333  momy0 = momy
334 !OCL XFILL
335  rhot0 = rhot
336 !OCL XFILL
337  if ( va > 0 ) prog0 = prog
338 
339  if ( bnd_w ) then
340  do j = js, je
341  do k = ks, ke
342  mflx_hi_rk(k,is-1,j,2,:) = mflx_hi(k,is-1,j,2)
343  end do
344  end do
345  end if
346  if ( bnd_e ) then
347  do j = js, je
348  do k = ks, ke
349  mflx_hi_rk(k,ie,j,2,:) = mflx_hi(k,ie,j,2)
350  end do
351  end do
352  end if
353  if ( bnd_s ) then
354  do i = is, ie
355  do k = ks, ke
356  mflx_hi_rk(k,i,js-1,3,:) = mflx_hi(k,i,js-1,3)
357  end do
358  end do
359  end if
360  if ( bnd_n ) then
361  do i = is, ie
362  do k = ks, ke
363  mflx_hi_rk(k,i,je,3,:) = mflx_hi(k,i,je,3)
364  end do
365  end do
366  end if
367 
368  call prof_rapend ("DYN_RK4_Prep",3)
369 
370  !------------------------------------------------------------------------
371  ! Start RK
372  !------------------------------------------------------------------------
373 
374  !##### RK1 : PROG0,PROG->PROG_RK1 #####
375 
376  call prof_rapstart("DYN_RK4",3)
377 
378  dtrk = dt / 2.0_rp
379 
380  call atmos_dyn_tstep( dens_rk1, momz_rk1, momx_rk1, momy_rk1, rhot_rk1, & ! [OUT]
381  prog_rk1, & ! [OUT]
382  mflx_hi_rk(:,:,:,:,1), tflx_hi_rk(:,:,:,:,1), & ! [INOUT]
383 
384  dens0, momz0, momx0, momy0, rhot0, & ! [IN]
385  dens, momz, momx, momy, rhot, & ! [IN]
386  dens_t, momz_t, momx_t, momy_t, rhot_t, & ! [IN]
387  prog0, prog, & ! [IN]
388  rtot, cvtot, corioli, & ! [IN]
389  num_diff, divdmp_coef, ddiv, & ! [IN]
390  flag_fct_momentum, flag_fct_t, & ! [IN]
391  flag_fct_along_stream, & ! [IN]
392  cdz, fdz, fdx, fdy, & ! [IN]
393  rcdz, rcdx, rcdy, rfdz, rfdx, rfdy, & ! [IN]
394  phi, gsqrt, j13g, j23g, j33g, mapf, & ! [IN]
395  ref_pres, ref_dens, & ! [IN]
396  bnd_w, bnd_e, bnd_s, bnd_n, & ! [IN]
397  dtrk, dt ) ! [IN]
398 
399  call prof_rapend ("DYN_RK4",3)
400  call prof_rapstart("DYN_RK4_BND",3)
401 
402  call atmos_dyn_copy_boundary( dens_rk1, momz_rk1, momx_rk1, momy_rk1, rhot_rk1, & ! [INOUT]
403  prog_rk1, & ! [INOUT]
404  dens0, momz0, momx0, momy0, rhot0, & ! [IN]
405  prog0, & ! [IN]
406  bnd_w, bnd_e, bnd_s, bnd_n ) ! [IN]
407 
408  call prof_rapend ("DYN_RK4_BND",3)
409 
410  call comm_vars8( dens_rk1(:,:,:), i_comm_dens_rk1 )
411  call comm_vars8( momz_rk1(:,:,:), i_comm_momz_rk1 )
412  call comm_vars8( momx_rk1(:,:,:), i_comm_momx_rk1 )
413  call comm_vars8( momy_rk1(:,:,:), i_comm_momy_rk1 )
414  call comm_vars8( rhot_rk1(:,:,:), i_comm_rhot_rk1 )
415  do iv = 1, va
416  call comm_vars8( prog_rk1(:,:,:,iv), i_comm_prog_rk1(iv) )
417  enddo
418 
419  call comm_wait ( dens_rk1(:,:,:), i_comm_dens_rk1, .false. )
420  call comm_wait ( momz_rk1(:,:,:), i_comm_momz_rk1, .false. )
421  call comm_wait ( momx_rk1(:,:,:), i_comm_momx_rk1, .false. )
422  call comm_wait ( momy_rk1(:,:,:), i_comm_momy_rk1, .false. )
423  call comm_wait ( rhot_rk1(:,:,:), i_comm_rhot_rk1, .false. )
424  do iv = 1, va
425  call comm_wait ( prog_rk1(:,:,:,iv), i_comm_prog_rk1(iv), .false. )
426  enddo
427 
428  !##### RK2 : PROG0,PROG_RK1->PROG_RK2 #####
429 
430  call prof_rapstart("DYN_RK4",3)
431 
432  dtrk = dt / 2.0_rp
433 
434  call atmos_dyn_tstep( dens_rk2, momz_rk2, momx_rk2, momy_rk2, rhot_rk2, & ! [OUT]
435  prog_rk2, & ! [OUT]
436  mflx_hi_rk(:,:,:,:,2), tflx_hi_rk(:,:,:,:,2), & ! [INOUT]
437  dens0, momz0, momx0, momy0, rhot0, & ! [IN]
438  dens_rk1, momz_rk1, momx_rk1, momy_rk1, rhot_rk1, & ! [IN]
439  dens_t, momz_t, momx_t, momy_t, rhot_t, & ! [IN]
440  prog0, prog_rk1, & ! [IN]
441  rtot, cvtot, corioli, & ! [IN]
442  num_diff, divdmp_coef, ddiv, & ! [IN]
443  flag_fct_momentum, flag_fct_t, & ! [IN]
444  flag_fct_along_stream, & ! [IN]
445  cdz, fdz, fdx, fdy, & ! [IN]
446  rcdz, rcdx, rcdy, rfdz, rfdx, rfdy, & ! [IN]
447  phi, gsqrt, j13g, j23g, j33g, mapf, & ! [IN]
448  ref_pres, ref_dens, & ! [IN]
449  bnd_w, bnd_e, bnd_s, bnd_n, & ! [IN]
450  dtrk, dt ) ! [IN]
451 
452  call prof_rapend ("DYN_RK4",3)
453  call prof_rapstart("DYN_RK4_BND",3)
454 
455  call atmos_dyn_copy_boundary( dens_rk2, momz_rk2, momx_rk2, momy_rk2, rhot_rk2, & ! [INOUT]
456  prog_rk2, & ! [INOUT]
457  dens0, momz0, momx0, momy0, rhot0, & ! [IN]
458  prog0, & ! [IN]
459  bnd_w, bnd_e, bnd_s, bnd_n ) ! [IN]
460 
461  call prof_rapend ("DYN_RK4_BND",3)
462 
463  call comm_vars8( dens_rk2(:,:,:), i_comm_dens_rk2 )
464  call comm_vars8( momz_rk2(:,:,:), i_comm_momz_rk2 )
465  call comm_vars8( momx_rk2(:,:,:), i_comm_momx_rk2 )
466  call comm_vars8( momy_rk2(:,:,:), i_comm_momy_rk2 )
467  call comm_vars8( rhot_rk2(:,:,:), i_comm_rhot_rk2 )
468  do iv = 1, va
469  call comm_vars8( prog_rk2(:,:,:,iv), i_comm_prog_rk2(iv) )
470  enddo
471 
472  call comm_wait ( dens_rk2(:,:,:), i_comm_dens_rk2, .false. )
473  call comm_wait ( momz_rk2(:,:,:), i_comm_momz_rk2, .false. )
474  call comm_wait ( momx_rk2(:,:,:), i_comm_momx_rk2, .false. )
475  call comm_wait ( momy_rk2(:,:,:), i_comm_momy_rk2, .false. )
476  call comm_wait ( rhot_rk2(:,:,:), i_comm_rhot_rk2, .false. )
477  do iv = 1, va
478  call comm_wait ( prog_rk2(:,:,:,iv), i_comm_prog_rk2(iv), .false. )
479  enddo
480 
481  !##### RK3 : PROG0,PROG_RK2->PROG_RK3 #####
482 
483  call prof_rapstart("DYN_RK4",3)
484 
485  dtrk = dt
486 
487  call atmos_dyn_tstep( dens_rk3, momz_rk3, momx_rk3, momy_rk3, rhot_rk3, & ! [OUT]
488  prog_rk3, & ! [OUT]
489  mflx_hi_rk(:,:,:,:,3), tflx_hi_rk(:,:,:,:,3), & ! [INOUT]
490  dens0, momz0, momx0, momy0, rhot0, & ! [IN]
491  dens_rk2, momz_rk2, momx_rk2, momy_rk2, rhot_rk2, & ! [IN]
492  dens_t, momz_t, momx_t, momy_t, rhot_t, & ! [IN]
493  prog0, prog_rk2, & ! [IN]
494  rtot, cvtot, corioli, & ! [IN]
495  num_diff, divdmp_coef, ddiv, & ! [IN]
496  flag_fct_momentum, flag_fct_t, & ! [IN]
497  flag_fct_along_stream, & ! [IN]
498  cdz, fdz, fdx, fdy, & ! [IN]
499  rcdz, rcdx, rcdy, rfdz, rfdx, rfdy, & ! [IN]
500  phi, gsqrt, j13g, j23g, j33g, mapf, & ! [IN]
501  ref_pres, ref_dens, & ! [IN]
502  bnd_w, bnd_e, bnd_s, bnd_n, & ! [IN]
503  dtrk, dt ) ! [IN]
504 
505  call prof_rapend ("DYN_RK4",3)
506  call prof_rapstart("DYN_RK4_BND",3)
507 
508  call atmos_dyn_copy_boundary( dens_rk3, momz_rk3, momx_rk3, momy_rk3, rhot_rk3, & ! [INOUT]
509  prog_rk3, & ! [INOUT]
510  dens0, momz0, momx0, momy0, rhot0, & ! [IN]
511  prog0, & ! [IN]
512  bnd_w, bnd_e, bnd_s, bnd_n ) ! [IN]
513 
514  call prof_rapend ("DYN_RK4_BND",3)
515 
516  call comm_vars8( dens_rk3(:,:,:), i_comm_dens_rk3 )
517  call comm_vars8( momz_rk3(:,:,:), i_comm_momz_rk3 )
518  call comm_vars8( momx_rk3(:,:,:), i_comm_momx_rk3 )
519  call comm_vars8( momy_rk3(:,:,:), i_comm_momy_rk3 )
520  call comm_vars8( rhot_rk3(:,:,:), i_comm_rhot_rk3 )
521  do iv = 1, va
522  call comm_vars8( prog_rk3(:,:,:,iv), i_comm_prog_rk3(iv) )
523  enddo
524 
525  call comm_wait ( dens_rk3(:,:,:), i_comm_dens_rk3, .false. )
526  call comm_wait ( momz_rk3(:,:,:), i_comm_momz_rk3, .false. )
527  call comm_wait ( momx_rk3(:,:,:), i_comm_momx_rk3, .false. )
528  call comm_wait ( momy_rk3(:,:,:), i_comm_momy_rk3, .false. )
529  call comm_wait ( rhot_rk3(:,:,:), i_comm_rhot_rk3, .false. )
530  do iv = 1, va
531  call comm_wait ( prog_rk3(:,:,:,iv), i_comm_prog_rk3(iv), .false. )
532  enddo
533 
534  !##### RK4 : PROG0,PROG_RK3->PROG #####
535 
536  call prof_rapstart("DYN_RK4",3)
537 
538  dtrk = dt
539 
540  call atmos_dyn_tstep( dens, momz, momx, momy, rhot, & ! [OUT]
541  prog, & ! [OUT]
542  mflx_hi, tflx_hi, & ! [INOUT]
543  dens0, momz0, momx0, momy0, rhot0, & ! [IN]
544  dens_rk3, momz_rk3, momx_rk3, momy_rk3, rhot_rk3, & ! [IN]
545  dens_t, momz_t, momx_t, momy_t, rhot_t, & ! [IN]
546  prog0, prog_rk3, & ! [IN]
547  rtot, cvtot, corioli, & ! [IN]
548  num_diff, divdmp_coef, ddiv, & ! [IN]
549  flag_fct_momentum, flag_fct_t, & ! [IN]
550  flag_fct_along_stream, & ! [IN]
551  cdz, fdz, fdx, fdy, & ! [IN]
552  rcdz, rcdx, rcdy, rfdz, rfdx, rfdy, & ! [IN]
553  phi, gsqrt, j13g, j23g, j33g, mapf, & ! [IN]
554  ref_pres, ref_dens, & ! [IN]
555  bnd_w, bnd_e, bnd_s, bnd_n, & ! [IN]
556  dtrk, dt ) ! [IN]
557 
558  do j = js, je
559  do i = is, ie
560  do k = ks, ke
561  dens(k,i,j) = ( dens_rk1(k,i,j) * 2.0_rp &
562  + dens_rk2(k,i,j) * 4.0_rp &
563  + dens_rk3(k,i,j) * 2.0_rp &
564  + dens(k,i,j) &
565  - dens0(k,i,j) * 3.0_rp ) / 6.0_rp
566  enddo
567  enddo
568  enddo
569 
570  do j = js, je
571  do i = is, ie
572  do k = ks, ke-1
573  momz(k,i,j) = ( momz_rk1(k,i,j) * 2.0_rp &
574  + momz_rk2(k,i,j) * 4.0_rp &
575  + momz_rk3(k,i,j) * 2.0_rp &
576  + momz(k,i,j) &
577  - momz0(k,i,j) * 3.0_rp ) / 6.0_rp
578  enddo
579  enddo
580  enddo
581 
582  do j = js, je
583  do i = is, ie
584  do k = ks, ke
585  momx(k,i,j) = ( momx_rk1(k,i,j) * 2.0_rp &
586  + momx_rk2(k,i,j) * 4.0_rp &
587  + momx_rk3(k,i,j) * 2.0_rp &
588  + momx(k,i,j) &
589  - momx0(k,i,j) * 3.0_rp ) / 6.0_rp
590  enddo
591  enddo
592  enddo
593 
594  do j = js, je
595  do i = is, ie
596  do k = ks, ke
597  momy(k,i,j) = ( momy_rk1(k,i,j) * 2.0_rp &
598  + momy_rk2(k,i,j) * 4.0_rp &
599  + momy_rk3(k,i,j) * 2.0_rp &
600  + momy(k,i,j) &
601  - momy0(k,i,j) * 3.0_rp ) / 6.0_rp
602  enddo
603  enddo
604  enddo
605 
606  do j = js, je
607  do i = is, ie
608  do k = ks, ke
609  rhot(k,i,j) = ( rhot_rk1(k,i,j) * 2.0_rp &
610  + rhot_rk2(k,i,j) * 4.0_rp &
611  + rhot_rk3(k,i,j) * 2.0_rp &
612  + rhot(k,i,j) &
613  - rhot0(k,i,j) * 3.0_rp ) / 6.0_rp
614  enddo
615  enddo
616  enddo
617 
618  do iv = 1, va
619  do j = js, je
620  do i = is, ie
621  do k = ks, ke
622  prog(k,i,j,iv) = ( prog_rk1(k,i,j,iv) * 2.0_rp &
623  + prog_rk2(k,i,j,iv) * 4.0_rp &
624  + prog_rk3(k,i,j,iv) * 2.0_rp &
625  + prog(k,i,j,iv) &
626  - prog0(k,i,j,iv) * 3.0_rp ) / 6.0_rp
627  enddo
628  enddo
629  enddo
630  enddo
631 
632 
633  do n = 1, 3
634  do j = js, je
635  do i = is, ie
636  do k = ks, ke
637  mflx_hi(k,i,j,n) = ( mflx_hi_rk(k,i,j,n,1) &
638  + mflx_hi_rk(k,i,j,n,2) * 2.0_rp &
639  + mflx_hi_rk(k,i,j,n,3) * 2.0_rp &
640  + mflx_hi(k,i,j,n ) ) / 6.0_rp
641  enddo
642  enddo
643  enddo
644  enddo
645 
646  do n = 1, 3
647  do j = js, je
648  do i = is, ie
649  do k = ks, ke
650  tflx_hi(k,i,j,n) = ( tflx_hi_rk(k,i,j,n,1) &
651  + tflx_hi_rk(k,i,j,n,2) * 2.0_rp &
652  + tflx_hi_rk(k,i,j,n,3) * 2.0_rp &
653  + tflx_hi(k,i,j,n ) ) / 6.0_rp
654  enddo
655  enddo
656  enddo
657  enddo
658 
659  call prof_rapend ("DYN_RK4",3)
660 
661  return
integer, public is
start point of inner domain: x, local
integer, public je
end point of inner domain: y, local
integer, public va
Definition: scale_index.F90:38
module Atmosphere / Dynamical scheme
integer, public ke
end point of inner domain: z, local
integer, public ia
of x whole cells (local, with HALO)
integer, public ka
of z whole cells (local, with HALO)
procedure(short), pointer, public atmos_dyn_tstep_short
module COMMUNICATION
Definition: scale_comm.F90:23
integer, public js
start point of inner domain: y, local
module Atmosphere / Dynamics common
integer, public ks
start point of inner domain: z, local
subroutine, public prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:132
subroutine, public atmos_dyn_copy_boundary(DENS, MOMZ, MOMX, MOMY, RHOT, PROG, DENS0, MOMZ0, MOMX0, MOMY0, RHOT0, PROG0, BND_W, BND_E, BND_S, BND_N)
integer, public ie
end point of inner domain: x, local
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:178
integer, public ja
of y whole cells (local, with HALO)
Here is the call graph for this function:
Here is the caller graph for this function: