15 #include "inc_openmp.h" 27 #if defined DEBUG || defined QUICKDEBUG 56 real(RP),
private,
allocatable :: DENS_RK1(:,:,:)
57 real(RP),
private,
allocatable :: MOMZ_RK1(:,:,:)
58 real(RP),
private,
allocatable :: MOMX_RK1(:,:,:)
59 real(RP),
private,
allocatable :: MOMY_RK1(:,:,:)
60 real(RP),
private,
allocatable :: RHOT_RK1(:,:,:)
61 real(RP),
private,
allocatable :: PROG_RK1(:,:,:,:)
62 real(RP),
private,
allocatable :: DENS_RK2(:,:,:)
63 real(RP),
private,
allocatable :: MOMZ_RK2(:,:,:)
64 real(RP),
private,
allocatable :: MOMX_RK2(:,:,:)
65 real(RP),
private,
allocatable :: MOMY_RK2(:,:,:)
66 real(RP),
private,
allocatable :: RHOT_RK2(:,:,:)
67 real(RP),
private,
allocatable :: PROG_RK2(:,:,:,:)
70 integer :: I_COMM_DENS_RK1 = 1
71 integer :: I_COMM_MOMZ_RK1 = 2
72 integer :: I_COMM_MOMX_RK1 = 3
73 integer :: I_COMM_MOMY_RK1 = 4
74 integer :: I_COMM_RHOT_RK1 = 5
75 integer,
allocatable :: I_COMM_PROG_RK1(:)
77 integer :: I_COMM_DENS_RK2 = 1
78 integer :: I_COMM_MOMZ_RK2 = 2
79 integer :: I_COMM_MOMX_RK2 = 3
80 integer :: I_COMM_MOMY_RK2 = 4
81 integer :: I_COMM_RHOT_RK2 = 5
82 integer,
allocatable :: I_COMM_PROG_RK2(:)
84 logical :: FLAG_WS2002
102 character(len=*) :: tinteg_type
107 select case( tinteg_type )
117 flag_ws2002 = .false.
118 fact_dt1 = 1.0_rp / 3.0_rp
119 fact_dt2 = 2.0_rp / 3.0_rp
121 if(
io_l )
write(
io_fid_log,*)
"*** RK3: Wichere and Skamarock (2002) is used" 128 fact_dt1 = 1.0_rp / 3.0_rp
129 fact_dt2 = 1.0_rp / 2.0_rp
131 write(*,*)
'xxx TINTEG_TYPE is not RK3. Check!' 135 allocate( dens_rk1(
ka,
ia,
ja) )
136 allocate( momz_rk1(
ka,
ia,
ja) )
137 allocate( momx_rk1(
ka,
ia,
ja) )
138 allocate( momy_rk1(
ka,
ia,
ja) )
139 allocate( rhot_rk1(
ka,
ia,
ja) )
141 allocate( dens_rk2(
ka,
ia,
ja) )
142 allocate( momz_rk2(
ka,
ia,
ja) )
143 allocate( momx_rk2(
ka,
ia,
ja) )
144 allocate( momy_rk2(
ka,
ia,
ja) )
145 allocate( rhot_rk2(
ka,
ia,
ja) )
147 allocate( prog_rk1(
ka,
ia,
ja,max(
va,1)) )
148 allocate( prog_rk2(
ka,
ia,
ja,max(
va,1)) )
149 allocate( i_comm_prog_rk1(max(
va,1)) )
150 allocate( i_comm_prog_rk2(max(
va,1)) )
158 i_comm_prog_rk1(iv) = 5 + iv
159 call comm_vars8_init(
'PROG_RK1', prog_rk1(:,:,:,iv), i_comm_prog_rk1(iv) )
168 i_comm_prog_rk2(iv) = 5 + iv
169 call comm_vars8_init(
'PROG_RK2', prog_rk2(:,:,:,iv), i_comm_prog_rk2(iv) )
172 dens_rk1(:,:,:) = undef
173 momz_rk1(:,:,:) = undef
174 momx_rk1(:,:,:) = undef
175 momy_rk1(:,:,:) = undef
176 rhot_rk1(:,:,:) = undef
177 if (
va > 0 ) prog_rk1(:,:,:,:) = undef
179 dens_rk2(:,:,:) = undef
180 momz_rk2(:,:,:) = undef
181 momx_rk2(:,:,:) = undef
182 momy_rk2(:,:,:) = undef
183 rhot_rk2(:,:,:) = undef
184 if (
va > 0 ) prog_rk2(:,:,:,:) = undef
192 DENS, MOMZ, MOMX, MOMY, RHOT, PROG, &
194 DENS_t, MOMZ_t, MOMX_t, MOMY_t, RHOT_t, &
195 Rtot, CVtot, CORIOLI, &
196 num_diff, wdamp_coef, divdmp_coef, DDIV, &
197 FLAG_FCT_MOMENTUM, FLAG_FCT_T, &
198 FLAG_FCT_ALONG_STREAM, &
199 CDZ, FDZ, FDX, FDY, &
200 RCDZ, RCDX, RCDY, RFDZ, RFDX, RFDY, &
201 PHI, GSQRT, J13G, J23G, J33G, MAPF, &
202 REF_pres, REF_dens, &
203 BND_W, BND_E, BND_S, BND_N, &
214 real(RP),
intent(inout) :: dens(
ka,
ia,
ja)
215 real(RP),
intent(inout) :: momz(
ka,
ia,
ja)
216 real(RP),
intent(inout) :: momx(
ka,
ia,
ja)
217 real(RP),
intent(inout) :: momy(
ka,
ia,
ja)
218 real(RP),
intent(inout) :: rhot(
ka,
ia,
ja)
219 real(RP),
intent(inout) :: prog(
ka,
ia,
ja,
va)
221 real(RP),
intent(inout) :: mflx_hi(
ka,
ia,
ja,3)
222 real(RP),
intent(inout) :: tflx_hi(
ka,
ia,
ja,3)
224 real(RP),
intent(in) :: dens_t(
ka,
ia,
ja)
225 real(RP),
intent(in) :: momz_t(
ka,
ia,
ja)
226 real(RP),
intent(in) :: momx_t(
ka,
ia,
ja)
227 real(RP),
intent(in) :: momy_t(
ka,
ia,
ja)
228 real(RP),
intent(in) :: rhot_t(
ka,
ia,
ja)
230 real(RP),
intent(in) :: rtot(
ka,
ia,
ja)
231 real(RP),
intent(in) :: cvtot(
ka,
ia,
ja)
232 real(RP),
intent(in) :: corioli(
ia,
ja)
233 real(RP),
intent(in) :: num_diff(
ka,
ia,
ja,5,3)
234 real(RP),
intent(in) :: wdamp_coef(
ka)
235 real(RP),
intent(in) :: divdmp_coef
236 real(RP),
intent(in) :: ddiv(
ka,
ia,
ja)
238 logical,
intent(in) :: flag_fct_momentum
239 logical,
intent(in) :: flag_fct_t
240 logical,
intent(in) :: flag_fct_along_stream
242 real(RP),
intent(in) :: cdz (
ka)
243 real(RP),
intent(in) :: fdz (
ka-1)
244 real(RP),
intent(in) :: fdx (
ia-1)
245 real(RP),
intent(in) :: fdy (
ja-1)
246 real(RP),
intent(in) :: rcdz(
ka)
247 real(RP),
intent(in) :: rcdx(
ia)
248 real(RP),
intent(in) :: rcdy(
ja)
249 real(RP),
intent(in) :: rfdz(
ka-1)
250 real(RP),
intent(in) :: rfdx(
ia-1)
251 real(RP),
intent(in) :: rfdy(
ja-1)
253 real(RP),
intent(in) :: phi (
ka,
ia,
ja)
254 real(RP),
intent(in) :: gsqrt(
ka,
ia,
ja,7)
255 real(RP),
intent(in) :: j13g (
ka,
ia,
ja,7)
256 real(RP),
intent(in) :: j23g (
ka,
ia,
ja,7)
257 real(RP),
intent(in) :: j33g
258 real(RP),
intent(in) :: mapf (
ia,
ja,2,4)
260 real(RP),
intent(in) :: ref_pres(
ka,
ia,
ja)
261 real(RP),
intent(in) :: ref_dens(
ka,
ia,
ja)
263 logical,
intent(in) :: bnd_w
264 logical,
intent(in) :: bnd_e
265 logical,
intent(in) :: bnd_s
266 logical,
intent(in) :: bnd_n
268 real(RP),
intent(in) :: dt
270 real(RP) :: dens0(
ka,
ia,
ja)
271 real(RP) :: momz0(
ka,
ia,
ja)
272 real(RP) :: momx0(
ka,
ia,
ja)
273 real(RP) :: momy0(
ka,
ia,
ja)
274 real(RP) :: rhot0(
ka,
ia,
ja)
277 real(RP) :: mflx_hi_rk(
ka,
ia,
ja,3,2)
278 real(RP) :: tflx_hi_rk(
ka,
ia,
ja,3,2)
282 integer :: i, j, k, iv, n
288 dens_rk1(:,:,:) = undef
289 momz_rk1(:,:,:) = undef
290 momx_rk1(:,:,:) = undef
291 momy_rk1(:,:,:) = undef
292 rhot_rk1(:,:,:) = undef
293 if (
va > 0 ) prog_rk1(:,:,:,:) = undef
295 dens_rk2(:,:,:) = undef
296 momz_rk2(:,:,:) = undef
297 momx_rk2(:,:,:) = undef
298 momy_rk2(:,:,:) = undef
299 rhot_rk2(:,:,:) = undef
300 if (
va > 0 ) prog_rk2(:,:,:,:) = undef
302 mflx_hi_rk(:,:,:,:,:) = undef
303 tflx_hi_rk(:,:,:,:,:) = undef
307 mflx_hi( 1:
ks-1,:,:,:) = undef
308 mflx_hi(
ke+1:
ka ,:,:,:) = undef
322 if (
va > 0 ) prog0 = prog
327 mflx_hi_rk(k,
is-1,j,2,:) = mflx_hi(k,
is-1,j,2)
334 mflx_hi_rk(k,
ie,j,2,:) = mflx_hi(k,
ie,j,2)
341 mflx_hi_rk(k,i,
js-1,3,:) = mflx_hi(k,i,
js-1,3)
348 mflx_hi_rk(k,i,
je,3,:) = mflx_hi(k,i,
je,3)
365 call atmos_dyn_tstep( dens_rk1, momz_rk1, momx_rk1, momy_rk1, rhot_rk1, &
367 mflx_hi_rk(:,:,:,:,1), tflx_hi_rk(:,:,:,:,1), &
368 dens0, momz0, momx0, momy0, rhot0, &
369 dens, momz, momx, momy, rhot, &
370 dens_t, momz_t, momx_t, momy_t, rhot_t, &
372 rtot, cvtot, corioli, &
373 num_diff, wdamp_coef, divdmp_coef, ddiv, &
374 flag_fct_momentum, flag_fct_t, &
375 flag_fct_along_stream, &
376 cdz, fdz, fdx, fdy, &
377 rcdz, rcdx, rcdy, rfdz, rfdx, rfdy, &
378 phi, gsqrt, j13g, j23g, j33g, mapf, &
379 ref_pres, ref_dens, &
380 bnd_w, bnd_e, bnd_s, bnd_n, &
388 dens0, momz0, momx0, momy0, rhot0, &
390 bnd_w, bnd_e, bnd_s, bnd_n )
394 call comm_vars8( dens_rk1(:,:,:), i_comm_dens_rk1 )
395 call comm_vars8( momz_rk1(:,:,:), i_comm_momz_rk1 )
396 call comm_vars8( momx_rk1(:,:,:), i_comm_momx_rk1 )
397 call comm_vars8( momy_rk1(:,:,:), i_comm_momy_rk1 )
398 call comm_vars8( rhot_rk1(:,:,:), i_comm_rhot_rk1 )
400 call comm_vars8( prog_rk1(:,:,:,iv), i_comm_prog_rk1(iv) )
403 call comm_wait ( dens_rk1(:,:,:), i_comm_dens_rk1, .false. )
404 call comm_wait ( momz_rk1(:,:,:), i_comm_momz_rk1, .false. )
405 call comm_wait ( momx_rk1(:,:,:), i_comm_momx_rk1, .false. )
406 call comm_wait ( momy_rk1(:,:,:), i_comm_momy_rk1, .false. )
407 call comm_wait ( rhot_rk1(:,:,:), i_comm_rhot_rk1, .false. )
409 call comm_wait ( prog_rk1(:,:,:,iv), i_comm_prog_rk1(iv), .false. )
418 call atmos_dyn_tstep( dens_rk2, momz_rk2, momx_rk2, momy_rk2, rhot_rk2, &
420 mflx_hi_rk(:,:,:,:,2), tflx_hi_rk(:,:,:,:,2), &
421 dens0, momz0, momx0, momy0, rhot0, &
422 dens_rk1, momz_rk1, momx_rk1, momy_rk1, rhot_rk1, &
423 dens_t, momz_t, momx_t, momy_t, rhot_t, &
425 rtot, cvtot, corioli, &
426 num_diff, wdamp_coef, divdmp_coef, ddiv, &
427 flag_fct_momentum, flag_fct_t, &
428 flag_fct_along_stream, &
429 cdz, fdz, fdx, fdy, &
430 rcdz, rcdx, rcdy, rfdz, rfdx, rfdy, &
431 phi, gsqrt, j13g, j23g, j33g, mapf, &
432 ref_pres, ref_dens, &
433 bnd_w, bnd_e, bnd_s, bnd_n, &
441 dens0, momz0, momx0, momy0, rhot0, &
443 bnd_w, bnd_e, bnd_s, bnd_n )
447 call comm_vars8( dens_rk2(:,:,:), i_comm_dens_rk2 )
448 call comm_vars8( momz_rk2(:,:,:), i_comm_momz_rk2 )
449 call comm_vars8( momx_rk2(:,:,:), i_comm_momx_rk2 )
450 call comm_vars8( momy_rk2(:,:,:), i_comm_momy_rk2 )
451 call comm_vars8( rhot_rk2(:,:,:), i_comm_rhot_rk2 )
453 call comm_vars8( prog_rk2(:,:,:,iv), i_comm_prog_rk2(iv) )
456 call comm_wait ( dens_rk2(:,:,:), i_comm_dens_rk2, .false. )
457 call comm_wait ( momz_rk2(:,:,:), i_comm_momz_rk2, .false. )
458 call comm_wait ( momx_rk2(:,:,:), i_comm_momx_rk2, .false. )
459 call comm_wait ( momy_rk2(:,:,:), i_comm_momy_rk2, .false. )
460 call comm_wait ( rhot_rk2(:,:,:), i_comm_rhot_rk2, .false. )
462 call comm_wait ( prog_rk2(:,:,:,iv), i_comm_prog_rk2(iv), .false. )
471 call atmos_dyn_tstep( dens, momz, momx, momy, rhot, &
474 dens0, momz0, momx0, momy0, rhot0, &
475 dens_rk2, momz_rk2, momx_rk2, momy_rk2, rhot_rk2, &
476 dens_t, momz_t, momx_t, momy_t, rhot_t, &
478 rtot, cvtot, corioli, &
479 num_diff, wdamp_coef, divdmp_coef, ddiv, &
480 flag_fct_momentum, flag_fct_t, &
481 flag_fct_along_stream, &
482 cdz, fdz, fdx, fdy, &
483 rcdz, rcdx, rcdy, rfdz, rfdx, rfdy, &
484 phi, gsqrt, j13g, j23g, j33g, mapf, &
485 ref_pres, ref_dens, &
486 bnd_w, bnd_e, bnd_s, bnd_n, &
489 if ( .NOT. flag_ws2002 )
then 493 dens(k,i,j) = ( dens_rk1(k,i,j) * 3.0_rp &
494 + dens(k,i,j) * 3.0_rp &
495 - dens0(k,i,j) * 2.0_rp ) / 4.0_rp
503 momz(k,i,j) = ( momz_rk1(k,i,j) * 3.0_rp &
504 + momz(k,i,j) * 3.0_rp &
505 - momz0(k,i,j) * 2.0_rp ) / 4.0_rp
513 momx(k,i,j) = ( momx_rk1(k,i,j) * 3.0_rp &
514 + momx(k,i,j) * 3.0_rp &
515 - momx0(k,i,j) * 2.0_rp ) / 4.0_rp
523 momy(k,i,j) = ( momy_rk1(k,i,j) * 3.0_rp &
524 + momy(k,i,j) * 3.0_rp &
525 - momy0(k,i,j) * 2.0_rp ) / 4.0_rp
533 rhot(k,i,j) = ( rhot_rk1(k,i,j) * 3.0_rp &
534 + rhot(k,i,j) * 3.0_rp &
535 - rhot0(k,i,j) * 2.0_rp ) / 4.0_rp
544 prog(k,i,j,iv) = ( prog_rk1(k,i,j,iv) * 3.0_rp &
545 + prog(k,i,j,iv) * 3.0_rp &
546 - prog0(k,i,j,iv) * 2.0_rp ) / 4.0_rp
556 mflx_hi(k,i,j,n) = ( mflx_hi_rk(k,i,j,n,1) &
557 + mflx_hi(k,i,j,n ) * 3.0_rp ) / 4.0_rp
567 tflx_hi(k,i,j,n) = ( tflx_hi_rk(k,i,j,n,1) &
568 + tflx_hi(k,i,j,n ) * 3.0_rp ) / 4.0_rp
integer, public is
start point of inner domain: x, local
integer, public je
end point of inner domain: y, local
subroutine, public prc_mpistop
Abort MPI.
subroutine, public comm_vars8_init(varname, var, vid)
Register variables.
module Atmosphere / Dynamical scheme
logical, public io_l
output log or not? (this process)
integer, public ke
end point of inner domain: z, local
subroutine, public check(current_line, v)
Undefined value checker.
real(rp), public const_undef
integer, public ia
of whole cells: x, local, with HALO
integer, public ka
of whole cells: z, local, with HALO
procedure(short), pointer, public atmos_dyn_tstep_short
subroutine, public atmos_dyn_tinteg_short_rk3(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, wdamp_coef, 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.
integer, public js
start point of inner domain: y, local
integer, parameter, public const_undef2
undefined value (INT2)
module Atmosphere / Dynamics common
module Atmosphere / Dyn Tinteg
integer, public ks
start point of inner domain: z, local
subroutine, public prof_rapstart(rapname_base, level)
Start raptime.
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)
subroutine, public atmos_dyn_tinteg_short_rk3_setup(tinteg_type)
Setup.
integer, public ie
end point of inner domain: x, local
integer, public io_fid_log
Log file ID.
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
integer, public ja
of whole cells: y, local, with HALO