Go to the documentation of this file.
67 wdamp_tau, wdamp_height, &
74 real(
rp),
intent(inout) :: wdamp_coef(
ka)
75 real(
rp),
intent(in) :: wdamp_tau
76 real(
rp),
intent(in) :: wdamp_height
77 real(
rp),
intent(in) :: fz(0:
ka)
84 if ( wdamp_height < 0.0_rp )
then
85 wdamp_coef(:) = 0.0_rp
86 elseif( fz(
ke)-wdamp_height < eps )
then
87 wdamp_coef(:) = 0.0_rp
89 alpha = 1.0_rp / wdamp_tau
92 sw = 0.5_rp + sign( 0.5_rp, fz(
k)-wdamp_height )
94 wdamp_coef(
k) = alpha * sw &
95 * 0.5_rp * ( 1.0_rp - cos( pi * (fz(
k)-wdamp_height) / (fz(
ke)-wdamp_height)) )
97 wdamp_coef( 1:
ks-1) = wdamp_coef(
ks)
98 wdamp_coef(
ke+1:
ka ) = wdamp_coef(
ke)
101 log_info(
"ATMOS_DYN_wdamp_setup",*)
'Setup Rayleigh damping coefficient'
102 log_info_cont(
'(1x,A)')
'|=== Rayleigh Damping Coef ===|'
103 log_info_cont(
'(1x,A)')
'| k zh[m] coef[/s] |'
105 log_info_cont(
'(1x,A,I5,F10.2,ES12.4,A)')
'| ',
k, fz(
k), wdamp_coef(
k),
' |'
108 log_info_cont(
'(1x,A,I5,F10.2,ES12.4,A)')
'| ',
k, fz(
k), wdamp_coef(
k),
' | KE = TOA'
110 log_info_cont(
'(1x,A,I5,F10.2,ES12.4,A)')
'| ',
k, fz(
k), wdamp_coef(
k),
' |'
113 log_info_cont(
'(1x,A,I5,F10.2,ES12.4,A)')
'| ',
k, fz(
k), wdamp_coef(
k),
' | KS-1 = surface'
115 log_info_cont(
'(1x,A,I5,F10.2,ES12.4,A)')
'| ',
k, fz(
k), wdamp_coef(
k),
' |'
118 log_info_cont(
'(1x,A,I5,F10.2,12x,A)')
'| ',
k, fz(
k),
' |'
119 log_info_cont(
'(1x,A)')
'|=============================|'
128 fill_constval, lateral_halo, top_bottom_halo )
131 real(
rp),
intent(inout) :: var(
ka,
ia,
ja)
132 real(
rp),
intent(in) :: fill_constval
133 logical,
intent(in) :: lateral_halo
134 logical,
intent(in) :: top_bottom_halo
139 if (lateral_halo)
then
144 var(
k,i,j) = fill_constval
149 var(
k,i,j) = fill_constval
157 var(
k,i,j) = fill_constval
165 var(
k,i,j) = fill_constval
171 if (top_bottom_halo)
then
175 var( 1:
ks-1,i,j) = fill_constval
176 var(
ke+1:
ka ,i,j) = fill_constval
185 DENS, MOMZ, MOMX, MOMY, RHOT, PROG, &
186 DENS0, MOMZ0, MOMX0, MOMY0, RHOT0, PROG0, &
187 BND_W, BND_E, BND_S, BND_N, TwoD )
189 real(
rp),
intent(inout) :: dens (
ka,
ia,
ja)
190 real(
rp),
intent(inout) :: momz (
ka,
ia,
ja)
191 real(
rp),
intent(inout) :: momx (
ka,
ia,
ja)
192 real(
rp),
intent(inout) :: momy (
ka,
ia,
ja)
193 real(
rp),
intent(inout) :: rhot (
ka,
ia,
ja)
195 real(
rp),
intent(in) :: dens0(
ka,
ia,
ja)
196 real(
rp),
intent(in) :: momz0(
ka,
ia,
ja)
197 real(
rp),
intent(in) :: momx0(
ka,
ia,
ja)
198 real(
rp),
intent(in) :: momy0(
ka,
ia,
ja)
199 real(
rp),
intent(in) :: rhot0(
ka,
ia,
ja)
201 logical,
intent(in) :: bnd_w
202 logical,
intent(in) :: bnd_e
203 logical,
intent(in) :: bnd_s
204 logical,
intent(in) :: bnd_n
205 logical,
intent(in) :: twod
207 integer ::
k, i, j, iv
209 if ( bnd_w .and. (.not. twod) )
then
217 dens(
k,i,j) = dens0(
k,i,j)
218 momz(
k,i,j) = momz0(
k,i,j)
219 momx(
k,i,j) = momx0(
k,i,j)
220 momy(
k,i,j) = momy0(
k,i,j)
221 rhot(
k,i,j) = rhot0(
k,i,j)
223 prog(
k,i,j,iv) = prog0(
k,i,j,iv)
229 if ( bnd_e .and. (.not. twod) )
then
237 dens(
k,i,j) = dens0(
k,i,j)
238 momz(
k,i,j) = momz0(
k,i,j)
239 momx(
k,i,j) = momx0(
k,i,j)
240 momy(
k,i,j) = momy0(
k,i,j)
241 rhot(
k,i,j) = rhot0(
k,i,j)
243 prog(
k,i,j,iv) = prog0(
k,i,j,iv)
252 momx(
k,
ie,j) = momx0(
k,
ie,j)
264 dens(
k,i,j) = dens0(
k,i,j)
265 momz(
k,i,j) = momz0(
k,i,j)
266 momx(
k,i,j) = momx0(
k,i,j)
267 momy(
k,i,j) = momy0(
k,i,j)
268 rhot(
k,i,j) = rhot0(
k,i,j)
270 prog(
k,i,j,iv) = prog0(
k,i,j,iv)
284 dens(
k,i,j) = dens0(
k,i,j)
285 momz(
k,i,j) = momz0(
k,i,j)
286 momx(
k,i,j) = momx0(
k,i,j)
287 momy(
k,i,j) = momy0(
k,i,j)
288 rhot(
k,i,j) = rhot0(
k,i,j)
290 prog(
k,i,j,iv) = prog0(
k,i,j,iv)
299 momy(
k,i,
je) = momy0(
k,i,
je)
310 BND_W, BND_E, BND_S, BND_N, TwoD )
312 real(
rp),
intent(inout) :: qtrc (
ka,
ia,
ja)
313 real(
rp),
intent(in) :: qtrc0(
ka,
ia,
ja)
314 logical,
intent(in) :: bnd_w
315 logical,
intent(in) :: bnd_e
316 logical,
intent(in) :: bnd_s
317 logical,
intent(in) :: bnd_n
318 logical,
intent(in) :: twod
322 if ( bnd_w .and. (.not. twod) )
then
329 qtrc(
k,i,j) = qtrc0(
k,i,j)
334 if ( bnd_e .and. (.not. twod) )
then
341 qtrc(
k,i,j) = qtrc0(
k,i,j)
353 qtrc(
k,i,j) = qtrc0(
k,i,j)
365 qtrc(
k,i,j) = qtrc0(
k,i,j)
379 GSQRT, J13G, J23G, J33G, MAPF, &
381 RCDZ, RCDX, RCDY, RFDZ, FDZ )
383 real(
rp),
intent(out) :: ddiv(
ka,
ia,
ja)
384 real(
rp),
intent(in) :: momz(
ka,
ia,
ja)
385 real(
rp),
intent(in) :: momx(
ka,
ia,
ja)
386 real(
rp),
intent(in) :: momy(
ka,
ia,
ja)
387 real(
rp),
intent(in) :: gsqrt(
ka,
ia,
ja,7)
388 real(
rp),
intent(in) :: j13g(
ka,
ia,
ja,7)
389 real(
rp),
intent(in) :: j23g(
ka,
ia,
ja,7)
390 real(
rp),
intent(in) :: j33g
391 real(
rp),
intent(in) :: mapf(
ia,
ja,2,7)
392 logical,
intent(in) :: twod
393 real(
rp),
intent(in) :: rcdz(
ka)
394 real(
rp),
intent(in) :: rcdx(
ia)
395 real(
rp),
intent(in) :: rcdy(
ja)
396 real(
rp),
intent(in) :: rfdz(
ka-1)
397 real(
rp),
intent(in) :: fdz(
ka-1)
409 ddiv(
k,
is,j) = j33g * ( momz(
k,
is,j) - momz(
k-1,
is,j) ) * rcdz(
k) &
410 + ( ( momy(
k+1,
is,j) + momy(
k+1,
is,j-1) ) * j23g(
k+1,
is,j,
i_xyw) &
411 - ( momy(
k-1,
is,j) + momy(
k-1,
is,j-1) ) * j23g(
k-1,
is,j,
i_xyw) ) / ( fdz(
k)+fdz(
k-1) ) &
418 k = iundef; i = iundef; j = iundef
422 ddiv(
ks,
is,j) = j33g * ( momz(
ks,
is,j) ) * rcdz(
ks) &
428 ddiv(
ke,
is,j) = j33g * ( - momz(
ke-1,
is,j ) ) * rcdz(
ke) &
436 k = iundef; i = iundef; j = iundef
443 ddiv(
k,i,j) = j33g * ( momz(
k,i,j) - momz(
k-1,i ,j ) ) * rcdz(
k) &
444 + ( ( momx(
k+1,i,j) + momx(
k+1,i-1,j ) ) * j13g(
k+1,i,j,
i_xyw) &
445 - ( momx(
k-1,i,j) + momx(
k-1,i-1,j ) ) * j13g(
k-1,i,j,
i_xyw) &
446 + ( momy(
k+1,i,j) + momy(
k+1,i ,j-1) ) * j23g(
k+1,i,j,
i_xyw) &
447 - ( momy(
k-1,i,j) + momy(
k-1,i ,j-1) ) * j23g(
k-1,i,j,
i_xyw) ) / ( fdz(
k)+fdz(
k-1) ) &
448 + mapf(i,j,1,
i_xy) * mapf(i,j,2,
i_xy) &
449 * ( ( momx(
k,i ,j ) * gsqrt(
k,i ,j ,
i_uyz) / mapf(i ,j ,2,
i_uy) &
450 - momx(
k,i-1,j ) * gsqrt(
k,i-1,j ,
i_uyz) / mapf(i-1,j ,2,
i_uy) ) * rcdx(i) &
451 + ( momy(
k,i ,j ) * gsqrt(
k,i ,j ,
i_xvz) / mapf(i ,j ,1,
i_xv) &
452 - momy(
k,i, j-1) * gsqrt(
k,i ,j-1,
i_xvz) / mapf(i ,j-1,1,
i_xv) ) * rcdy(j) )
457 k = iundef; i = iundef; j = iundef
462 ddiv(
ks,i,j) = j33g * ( momz(
ks,i,j) ) * rcdz(
ks) &
463 + ( ( momx(
ks+1,i,j) + momx(
ks+1,i-1,j ) ) * j13g(
ks+1,i,j,
i_xyw) &
464 - ( momx(
ks-1,i,j) + momx(
ks ,i-1,j ) ) * j13g(
ks ,i,j,
i_xyw) &
465 + ( momy(
ks+1,i,j) + momy(
ks+1,i ,j-1) ) * j23g(
ks+1,i,j,
i_xyw) &
466 - ( momy(
ks ,i,j) + momy(
ks ,i ,j-1) ) * j23g(
ks ,i,j,
i_xyw) ) * rfdz(
ks) &
467 + mapf(i,j,1,
i_xy) * mapf(i,j,2,
i_xy) &
468 * ( ( momx(
ks,i ,j ) * gsqrt(
ks,i ,j ,
i_uyz) / mapf(i ,j ,2,
i_uy) &
469 - momx(
ks,i-1,j ) * gsqrt(
ks,i-1,j ,
i_uyz) / mapf(i-1,j ,2,
i_uy) ) * rcdx(i) &
470 + ( momy(
ks,i ,j ) * gsqrt(
ks,i ,j ,
i_xvz) / mapf(i ,j ,1,
i_xv) &
471 - momy(
ks,i, j-1) * gsqrt(
ks,i ,j-1,
i_xvz) / mapf(i ,j-1,1,
i_xv) ) * rcdy(j) )
472 ddiv(
ke,i,j) = j33g * ( - momz(
ke-1,i ,j ) ) * rcdz(
ke) &
473 + ( ( momx(
ke ,i,j) + momx(
ke ,i-1,j ) ) * j13g(
ke ,i,j,
i_xyw) &
474 - ( momx(
ke-1,i,j) + momx(
ke-1,i-1,j ) ) * j13g(
ke-1,i,j,
i_xyw) &
475 + ( momy(
ke ,i,j) + momy(
ke ,i ,j-1) ) * j23g(
ke ,i,j,
i_xyw) &
476 - ( momy(
ke-1,i,j) + momy(
ke-1,i ,j-1) ) * j23g(
ke-1,i,j,
i_xyw) ) * rfdz(
ke) &
477 + mapf(i,j,1,
i_xy) * mapf(i,j,2,
i_xy) &
478 * ( ( momx(
ke,i ,j ) * gsqrt(
ke,i ,j ,
i_uyz) / mapf(i ,j ,2,
i_uy) &
479 - momx(
ke,i-1,j ) * gsqrt(
ke,i-1,j ,
i_uyz) / mapf(i-1,j ,2,
i_uy) ) * rcdx(i) &
480 + ( momy(
ke,i ,j ) * gsqrt(
ke,i ,j ,
i_xvz) / mapf(i ,j ,1,
i_xv) &
481 - momy(
ke,i, j-1) * gsqrt(
ke,i ,j-1,
i_xvz) / mapf(i ,j-1,1,
i_xv) ) * rcdy(j) )
485 k = iundef; i = iundef; j = iundef
508 DPRES, RT2P, REF_rhot, & ! (out)
509 rhot, qtrc, ref_pres, aq_r, aq_cv, aq_cp, aq_mass )
518 real(
rp),
intent(out) :: dpres(
ka,
ia,
ja)
519 real(
rp),
intent(out) :: rt2p(
ka,
ia,
ja)
520 real(
rp),
intent(out) :: ref_rhot(
ka,
ia,
ja)
521 real(
rp),
intent(in) :: rhot(
ka,
ia,
ja)
523 real(
rp),
intent(in) :: ref_pres(
ka,
ia,
ja)
524 real(
rp),
intent(in) :: aq_r(
qa)
525 real(
rp),
intent(in) :: aq_cv(
qa)
526 real(
rp),
intent(in) :: aq_cp(
qa)
527 real(
rp),
intent(in) :: aq_mass(
qa)
544 cpovcv = cpdry / cvdry
560 pres = p0 * ( rdry * rhot(
k,i,j) / p0 )**cpovcv
561 rt2p(
k,i,j) = cpovcv * pres / rhot(
k,i,j)
568 rtot = rtot + aq_r(iq) * qtrc(
k,i,j,iq)
569 cvtot = cvtot + aq_cv(iq) * qtrc(
k,i,j,iq)
570 cptot = cptot + aq_cp(iq) * qtrc(
k,i,j,iq)
571 qdry = qdry - qtrc(
k,i,j,iq) * aq_mass(iq)
573 rtot = rtot + rdry * qdry
574 cvtot = cvtot + cvdry * qdry
575 cptot = cptot + cpdry * qdry
576 pres = p0 * ( rtot * rhot(
k,i,j) / p0 )**( cptot / cvtot )
577 rt2p(
k,i,j) = cptot / cvtot * pres / rhot(
k,i,j)
579 dpres(
k,i,j) = pres - ref_pres(
k,i,j)
580 ref_rhot(
k,i,j) = rhot(
k,i,j)
582 dpres(
ks-1,i,j) = dpres(
ks+1,i,j) + ( ref_pres(
ks+1,i,j) - ref_pres(
ks-1,i,j) )
583 dpres(
ke+1,i,j) = dpres(
ke-1,i,j) + ( ref_pres(
ke-1,i,j) - ref_pres(
ke+1,i,j) )
integer, public ke
end point of inner domain: z, local
integer, parameter, public const_undef2
undefined value (INT2)
module Atmosphere / Dynamics common
real(rp), public const_eps
small number
subroutine, public prof_rapstart(rapname_base, level, disable_barrier)
Start raptime.
real(rp), public const_pi
pi
integer, parameter, public rp
integer, public ie
end point of inner domain: x, local
module atmosphere / grid / cartesC index
subroutine, public check(current_line, v)
Undefined value checker.
subroutine, public atmos_dyn_wdamp_setup(wdamp_coef, wdamp_tau, wdamp_height, FZ)
Setup.
real(rp), public const_cvdry
specific heat (dry air,constant volume) [J/kg/K]
subroutine, public atmos_dyn_fill_halo(var, fill_constval, lateral_halo, top_bottom_halo)
real(rp), public const_cpdry
specific heat (dry air,constant pressure) [J/kg/K]
subroutine, public atmos_dyn_copy_boundary_tracer(QTRC, QTRC0, BND_W, BND_E, BND_S, BND_N, TwoD)
integer, public is
start point of inner domain: x, local
integer, public ks
start point of inner domain: z, local
subroutine, public atmos_dyn_divergence(DDIV, MOMZ, MOMX, MOMY, GSQRT, J13G, J23G, J33G, MAPF, TwoD, RCDZ, RCDX, RCDY, RFDZ, FDZ)
integer, public js
start point of inner domain: y, local
real(rp), public const_rdry
specific gas constant (dry air) [J/kg/K]
subroutine, public prof_rapend(rapname_base, level, disable_barrier)
Save raptime.
subroutine, public atmos_dyn_prep_pres_linearization(DPRES, RT2P, REF_rhot, RHOT, QTRC, REF_pres, AQ_R, AQ_CV, AQ_CP, AQ_MASS)
real(rp), public const_undef
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, TwoD)
real(rp), public const_pre00
pressure reference [Pa]
integer, public je
end point of inner domain: y, local