SCALE-RM
scale_atmos_dyn_common.F90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
10 !-------------------------------------------------------------------------------
11 #include "scalelib.h"
13  !-----------------------------------------------------------------------------
14  !
15  !++ used modules
16  !
17  use scale_precision
18  use scale_io
19  use scale_prof
21  use scale_index
22  use scale_tracer
23 
24 #ifdef DEBUG
25  use scale_debug, only: &
26  check
27  use scale_const, only: &
28  undef => const_undef, &
29  iundef => const_undef2
30 #endif
31  !-----------------------------------------------------------------------------
32  implicit none
33  private
34  !-----------------------------------------------------------------------------
35  !
36  !++ Public procedure
37  !
38  public :: atmos_dyn_wdamp_setup
39  public :: atmos_dyn_fill_halo
40  public :: atmos_dyn_copy_boundary
42  public :: atmos_dyn_divergence
44 
45 
46  !-----------------------------------------------------------------------------
47  !
48  !++ Public parameters & variables
49  !
50  !-----------------------------------------------------------------------------
51  !
52  !++ Private procedure
53  !
54 
55  !-----------------------------------------------------------------------------
56  !
57  !++ Private parameters & variables
58  !
59  !-----------------------------------------------------------------------------
60 
61 contains
62 
63  !-----------------------------------------------------------------------------
65  subroutine atmos_dyn_wdamp_setup( &
66  wdamp_coef, &
67  wdamp_tau, wdamp_height, &
68  FZ )
69  use scale_const, only: &
70  pi => const_pi, &
71  eps => const_eps
72  implicit none
73 
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)
78 
79  real(rp) :: alpha, sw
80 
81  integer :: k
82  !---------------------------------------------------------------------------
83 
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
88  else
89  alpha = 1.0_rp / wdamp_tau
90 
91  do k = ks, ke
92  sw = 0.5_rp + sign( 0.5_rp, fz(k)-wdamp_height )
93 
94  wdamp_coef(k) = alpha * sw &
95  * 0.5_rp * ( 1.0_rp - cos( pi * (fz(k)-wdamp_height) / (fz(ke)-wdamp_height)) )
96  enddo
97  wdamp_coef( 1:ks-1) = wdamp_coef(ks)
98  wdamp_coef(ke+1:ka ) = wdamp_coef(ke)
99 
100  log_newline
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] |'
104  do k = ka, ke+1, -1
105  log_info_cont('(1x,A,I5,F10.2,ES12.4,A)') '| ',k, fz(k), wdamp_coef(k),' |'
106  enddo
107  k = ke
108  log_info_cont('(1x,A,I5,F10.2,ES12.4,A)') '| ',k, fz(k), wdamp_coef(k),' | KE = TOA'
109  do k = ke-1, ks, -1
110  log_info_cont('(1x,A,I5,F10.2,ES12.4,A)') '| ',k, fz(k), wdamp_coef(k),' |'
111  enddo
112  k = ks-1
113  log_info_cont('(1x,A,I5,F10.2,ES12.4,A)') '| ',k, fz(k), wdamp_coef(k),' | KS-1 = surface'
114  do k = ks-2, 1, -1
115  log_info_cont('(1x,A,I5,F10.2,ES12.4,A)') '| ',k, fz(k), wdamp_coef(k),' |'
116  enddo
117  k = 0
118  log_info_cont('(1x,A,I5,F10.2,12x,A)') '| ',k, fz(k), ' |'
119  log_info_cont('(1x,A)') '|=============================|'
120  endif
121 
122  return
123  end subroutine atmos_dyn_wdamp_setup
124 
125  !-----------------------------------------------------------------------------
126 
127  subroutine atmos_dyn_fill_halo( var, &
128  fill_constval, lateral_halo, top_bottom_halo )
129  implicit none
130 
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
135 
136  integer :: i, j, k
137  !----------------------------
138 
139  if (lateral_halo) then
140 !OCL XFILL
141  do j = 1, ja
142  do i = 1, isb-1
143  do k = 1, ka
144  var(k,i,j) = fill_constval
145  enddo
146  enddo
147  do i = ieb+1, ia
148  do k = 1, ka
149  var(k,i,j) = fill_constval
150  enddo
151  enddo
152  enddo
153 !OCL XFILL
154  do j = 1, jsb-1
155  do i = 1, ia
156  do k = 1, ka
157  var(k,i,j) = fill_constval
158  enddo
159  enddo
160  enddo
161 !OCL XFILL
162  do j = jeb+1, ja
163  do i = 1, ia
164  do k = 1, ka
165  var(k,i,j) = fill_constval
166  enddo
167  enddo
168  enddo
169  end if
170 
171  if (top_bottom_halo) then
172  !OCL XFILL
173  do j = js, je
174  do i = is, ie
175  var( 1:ks-1,i,j) = fill_constval
176  var(ke+1:ka ,i,j) = fill_constval
177  enddo
178  enddo
179  end if
180 
181  return
182  end subroutine atmos_dyn_fill_halo
183 
184  subroutine atmos_dyn_copy_boundary( &
185  DENS, MOMZ, MOMX, MOMY, RHOT, PROG, &
186  DENS0, MOMZ0, MOMX0, MOMY0, RHOT0, PROG0, &
187  BND_W, BND_E, BND_S, BND_N, TwoD )
188  implicit none
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)
194  real(rp), intent(inout) :: prog (ka,ia,ja,va)
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)
200  real(rp), intent(in) :: prog0(ka,ia,ja,va)
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
206 
207  integer :: k, i, j, iv
208 
209  if ( bnd_w .and. (.not. twod) ) then
210  !$omp parallel do default(none) private(j,k) OMP_SCHEDULE_ collapse(2) &
211  !$omp private(i,iv) &
212  !$omp shared(JA,IS,KS,KE,DENS,DENS0,MOMZ,MOMZ0,MOMX,MOMX0,MOMY,MOMY0,RHOT,RHOT0,VA,PROG,PROG0)
213 !OCL XFILL
214  do j = 1, ja
215  do i = 1, is-1
216  do k = ks, ke
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)
222  do iv = 1, va
223  prog(k,i,j,iv) = prog0(k,i,j,iv)
224  end do
225  enddo
226  enddo
227  enddo
228  end if
229  if ( bnd_e .and. (.not. twod) ) then
230  !$omp parallel do default(none) private(j,k) OMP_SCHEDULE_ collapse(2) &
231  !$omp private(i,iv) &
232  !$omp shared(JA,IE,IA,KS,KE,DENS,DENS0,MOMZ,MOMZ0,MOMX,MOMX0,MOMY,MOMY0,RHOT,RHOT0,VA,PROG,PROG0)
233 !OCL XFILL
234  do j = 1, ja
235  do i = ie+1, ia
236  do k = ks, ke
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)
242  do iv = 1, va
243  prog(k,i,j,iv) = prog0(k,i,j,iv)
244  end do
245  enddo
246  enddo
247  enddo
248  !$omp parallel do private(j,k) OMP_SCHEDULE_ collapse(2)
249 !OCL XFILL
250  do j = 1, ja
251  do k = ks, ke
252  momx(k,ie,j) = momx0(k,ie,j)
253  enddo
254  enddo
255  end if
256  if ( bnd_s ) then
257  !$omp parallel do default(none) private(j,k) OMP_SCHEDULE_ collapse(2) &
258  !$omp private(i,iv) &
259  !$omp shared(JS,IA,KS,KE,DENS,DENS0,MOMZ,MOMZ0,MOMX,MOMX0,MOMY,MOMY0,RHOT,RHOT0,VA,PROG,PROG0)
260 !OCL XFILL
261  do j = 1, js-1
262  do i = 1, ia
263  do k = ks, ke
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)
269  do iv = 1, va
270  prog(k,i,j,iv) = prog0(k,i,j,iv)
271  end do
272  enddo
273  enddo
274  enddo
275  end if
276  if ( bnd_n ) then
277  !$omp parallel do default(none) private(j,k) OMP_SCHEDULE_ collapse(2) &
278  !$omp private(i,iv) &
279  !$omp shared(JA,JE,IA,KS,KE,DENS,DENS0,MOMZ,MOMZ0,MOMX,MOMX0,MOMY,MOMY0,RHOT,RHOT0,VA,PROG,PROG0)
280 !OCL XFILL
281  do j = je+1, ja
282  do i = 1, ia
283  do k = ks, ke
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)
289  do iv = 1, va
290  prog(k,i,j,iv) = prog0(k,i,j,iv)
291  end do
292  enddo
293  enddo
294  enddo
295  !$omp parallel do private(i,k) OMP_SCHEDULE_ collapse(2)
296 !OCL XFILL
297  do i = 1, ia
298  do k = ks, ke
299  momy(k,i,je) = momy0(k,i,je)
300  enddo
301  enddo
302  end if
303 
304  return
305  end subroutine atmos_dyn_copy_boundary
306 
307  !-----------------------------------------------------------------------------
308  subroutine atmos_dyn_copy_boundary_tracer( &
309  QTRC, QTRC0, &
310  BND_W, BND_E, BND_S, BND_N, TwoD )
311  implicit none
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
319 
320  integer :: k, i, j
321 
322  if ( bnd_w .and. (.not. twod) ) then
323  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
324  !$omp shared(JA,IS,KS,KE,QTRC,QTRC0)
325 !OCL XFILL
326  do j = 1, ja
327  do i = 1, is-1
328  do k = ks, ke
329  qtrc(k,i,j) = qtrc0(k,i,j)
330  enddo
331  enddo
332  enddo
333  end if
334  if ( bnd_e .and. (.not. twod) ) then
335  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
336  !$omp shared(JA,IE,IA,KS,KE,QTRC,QTRC0)
337 !OCL XFILL
338  do j = 1, ja
339  do i = ie+1, ia
340  do k = ks, ke
341  qtrc(k,i,j) = qtrc0(k,i,j)
342  enddo
343  enddo
344  enddo
345  end if
346  if ( bnd_s ) then
347  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
348  !$omp shared(JS,IA,KS,KE,QTRC,QTRC0)
349 !OCL XFILL
350  do j = 1, js-1
351  do i = 1, ia
352  do k = ks, ke
353  qtrc(k,i,j) = qtrc0(k,i,j)
354  enddo
355  enddo
356  enddo
357  end if
358  if ( bnd_n ) then
359  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
360  !$omp shared(JA,JE,IA,KS,KE,QTRC,QTRC0)
361 !OCL XFILL
362  do j = je+1, ja
363  do i = 1, ia
364  do k = ks, ke
365  qtrc(k,i,j) = qtrc0(k,i,j)
366  enddo
367  enddo
368  enddo
369  end if
370 
371  return
372  end subroutine atmos_dyn_copy_boundary_tracer
373 
374  !-----------------------------------------------------------------------------
375 
376  subroutine atmos_dyn_divergence( &
377  DDIV, &
378  MOMZ, MOMX, MOMY, &
379  GSQRT, J13G, J23G, J33G, MAPF, &
380  TwoD, &
381  RCDZ, RCDX, RCDY, RFDZ, FDZ )
382  implicit none
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)
398 
399  integer :: k, i, j
400 
401  call prof_rapstart("DYN_divercence", 2)
402 
403  ! 3D divergence
404 
405  if ( twod ) then
406  !$omp parallel do private(j,k) OMP_SCHEDULE_
407  do j = js, je+1
408  do k = ks-1, ke+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) ) &
412  + mapf(is,j,2,i_xy) &
413  * ( momy(k,is,j ) * gsqrt(k,is,j ,i_xvz) / mapf(is,j ,1,i_xv) &
414  - momy(k,is,j-1) * gsqrt(k,is,j-1,i_xvz) / mapf(is,j-1,1,i_xv) ) * rcdy(j)
415  enddo
416  enddo
417 #ifdef DEBUG
418  k = iundef; i = iundef; j = iundef
419 #endif
420  !$omp parallel do private(j) OMP_SCHEDULE_
421  do j = js, je+1
422  ddiv(ks,is,j) = j33g * ( momz(ks,is,j) ) * rcdz(ks) &
423  + ( ( momy(ks+1,is,j) + momy(ks+1,is,j-1) ) * j23g(ks+1,is,j,i_xyw) &
424  - ( momy(ks ,is,j) + momy(ks ,is,j-1) ) * j23g(ks ,is,j,i_xyw) ) * rfdz(ks) &
425  + mapf(is,j,2,i_xy) &
426  * ( momy(ks,is,j ) * gsqrt(ks,is,j ,i_xvz) / mapf(is,j ,1,i_xv) &
427  - momy(ks,is,j-1) * gsqrt(ks,is,j-1,i_xvz) / mapf(is,j-1,1,i_xv) ) * rcdy(j)
428  ddiv(ke,is,j) = j33g * ( - momz(ke-1,is,j ) ) * rcdz(ke) &
429  + ( ( momy(ke ,is,j) + momy(ke ,is,j-1) ) * j23g(ke ,is,j,i_xyw) &
430  - ( momy(ke-1,is,j) + momy(ke-1,is,j-1) ) * j23g(ke-1,is,j,i_xyw) ) * rfdz(ke) &
431  + mapf(is,j,2,i_xy) &
432  * ( momy(ke,is,j ) * gsqrt(ke,is,j ,i_xvz) / mapf(is,j ,1,i_xv) &
433  - momy(ke,is,j-1) * gsqrt(ke,is,j-1,i_xvz) / mapf(is,j-1,1,i_xv) ) * rcdy(j)
434  enddo
435 #ifdef DEBUG
436  k = iundef; i = iundef; j = iundef
437 #endif
438  else
439  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
440  do j = js, je+1
441  do i = is, ie+1
442  do k = ks-1, ke+1
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) )
453  enddo
454  enddo
455  enddo
456 #ifdef DEBUG
457  k = iundef; i = iundef; j = iundef
458 #endif
459  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
460  do j = js, je+1
461  do i = is, ie+1
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) )
482  enddo
483  enddo
484 #ifdef DEBUG
485  k = iundef; i = iundef; j = iundef
486 #endif
487  end if
488  call prof_rapend ("DYN_divercence", 2)
489 
490  return
491  end subroutine atmos_dyn_divergence
492 
493 
494  !------------------------------------------------------------------------
495  ! prepare thermodynamical data
496  ! specific heat
497  ! pressure data ( linearization )
498  !
499  ! pres = P0 * ( R * rhot / P0 )**(CP/CV)
500  ! d pres / d rhot ~ CP*R/CV * ( R * rhot / P0 )**(R/CV)
501  ! = CP*R/CV * ( pres / P0 )**(R/CP)
502  ! = CP*R/CV * temp / pott
503  ! = CP/CV * pres / rhot
504  ! pres ~ P0 * ( R * rhot0 / P0 ) ** (CP/CV) + CV*R/CP * ( pres / P0 )**(R/CP) * rhot'
505  !------------------------------------------------------------------------
506 
508  DPRES, RT2P, REF_rhot, & ! (out)
509  rhot, qtrc, ref_pres, aq_r, aq_cv, aq_cp, aq_mass ) ! (in)
510 
511  use scale_const, only: &
512  p0 => const_pre00, &
513  rdry => const_rdry, &
514  cvdry => const_cvdry, &
515  cpdry => const_cpdry
516  implicit none
517 
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)
522  real(rp), intent(in) :: qtrc(ka,ia,ja,qa)
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)
528 
529  integer :: i, j, k
530  integer :: iq
531  real(rp) :: qdry ! dry air
532  real(rp) :: rtot ! total R
533  real(rp) :: cvtot ! total CV
534  real(rp) :: cptot ! total CP
535  real(rp) :: pres ! pressure
536 
537 #ifdef DRY
538  real(rp) :: cpovcv
539 #endif
540 
541  !--------------------------------------
542 
543 #ifdef DRY
544  cpovcv = cpdry / cvdry
545 #endif
546 
547 !OCL XFILL
548  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
549  !$omp shared(JA,IA,KS,KE) &
550  !$omp shared(P0,Rdry,RHOT,AQ_R,AQ_CV,AQ_CP,QTRC,AQ_MASS,REF_rhot,REF_pres,CPdry,CVdry,QA,RT2P,DPRES) &
551 #ifdef DRY
552  !$omp shared(CPovCV) &
553 #endif
554  !$omp private(i,j,k,iq) &
555  !$omp private(PRES,Rtot,CVtot,CPtot,QDRY)
556  do j = 1, ja
557  do i = 1, ia
558  do k = ks, ke
559 #ifdef DRY
560  pres = p0 * ( rdry * rhot(k,i,j) / p0 )**cpovcv
561  rt2p(k,i,j) = cpovcv * pres / rhot(k,i,j)
562 #else
563  rtot = 0.0_rp
564  cvtot = 0.0_rp
565  cptot = 0.0_rp
566  qdry = 1.0_rp
567  do iq = 1, qa
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)
572  enddo
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)
578 #endif
579  dpres(k,i,j) = pres - ref_pres(k,i,j)
580  ref_rhot(k,i,j) = rhot(k,i,j)
581  end do
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) )
584  end do
585  end do
586 
587  return
588  end subroutine atmos_dyn_prep_pres_linearization
589 
590 end module scale_atmos_dyn_common
scale_atmos_grid_cartesc_index::isb
integer, public isb
Definition: scale_atmos_grid_cartesC_index.F90:63
scale_atmos_grid_cartesc_index::i_uy
integer, public i_uy
Definition: scale_atmos_grid_cartesC_index.F90:99
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_atmos_grid_cartesc_index::i_xv
integer, public i_xv
Definition: scale_atmos_grid_cartesC_index.F90:100
scale_tracer::qa
integer, public qa
Definition: scale_tracer.F90:34
scale_index
module Index
Definition: scale_index.F90:11
scale_const::const_undef2
integer, parameter, public const_undef2
undefined value (INT2)
Definition: scale_const.F90:38
scale_atmos_dyn_common
module Atmosphere / Dynamics common
Definition: scale_atmos_dyn_common.F90:12
scale_precision
module PRECISION
Definition: scale_precision.F90:14
scale_atmos_grid_cartesc_index::ka
integer, public ka
Definition: scale_atmos_grid_cartesC_index.F90:47
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:33
scale_prof::prof_rapstart
subroutine, public prof_rapstart(rapname_base, level, disable_barrier)
Start raptime.
Definition: scale_prof.F90:159
scale_atmos_grid_cartesc_index::jeb
integer, public jeb
Definition: scale_atmos_grid_cartesC_index.F90:66
scale_index::va
integer, public va
Definition: scale_index.F90:35
scale_const::const_pi
real(rp), public const_pi
pi
Definition: scale_const.F90:31
scale_atmos_grid_cartesc_index::i_xy
integer, public i_xy
Definition: scale_atmos_grid_cartesC_index.F90:98
scale_precision::rp
integer, parameter, public rp
Definition: scale_precision.F90:41
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_atmos_grid_cartesc_index::i_uyz
integer, public i_uyz
Definition: scale_atmos_grid_cartesC_index.F90:94
scale_io
module STDIO
Definition: scale_io.F90:10
scale_tracer::k
real(rp), public k
Definition: scale_tracer.F90:44
scale_atmos_grid_cartesc_index
module atmosphere / grid / cartesC index
Definition: scale_atmos_grid_cartesC_index.F90:12
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_debug::check
subroutine, public check(current_line, v)
Undefined value checker.
Definition: scale_debug.F90:56
scale_atmos_dyn_common::atmos_dyn_wdamp_setup
subroutine, public atmos_dyn_wdamp_setup(wdamp_coef, wdamp_tau, wdamp_height, FZ)
Setup.
Definition: scale_atmos_dyn_common.F90:69
scale_atmos_grid_cartesc_index::i_xvz
integer, public i_xvz
Definition: scale_atmos_grid_cartesC_index.F90:95
scale_const::const_cvdry
real(rp), public const_cvdry
specific heat (dry air,constant volume) [J/kg/K]
Definition: scale_const.F90:57
scale_atmos_dyn_common::atmos_dyn_fill_halo
subroutine, public atmos_dyn_fill_halo(var, fill_constval, lateral_halo, top_bottom_halo)
Definition: scale_atmos_dyn_common.F90:129
scale_const::const_cpdry
real(rp), public const_cpdry
specific heat (dry air,constant pressure) [J/kg/K]
Definition: scale_const.F90:56
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_prof
module profiler
Definition: scale_prof.F90:11
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_tracer
module TRACER
Definition: scale_tracer.F90:12
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_debug
module DEBUG
Definition: scale_debug.F90:11
scale_atmos_dyn_common::atmos_dyn_divergence
subroutine, public atmos_dyn_divergence(DDIV, MOMZ, MOMX, MOMY, GSQRT, J13G, J23G, J33G, MAPF, TwoD, RCDZ, RCDX, RCDY, RFDZ, FDZ)
Definition: scale_atmos_dyn_common.F90:382
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_atmos_grid_cartesc_index::ieb
integer, public ieb
Definition: scale_atmos_grid_cartesC_index.F90:64
scale_const::const_rdry
real(rp), public const_rdry
specific gas constant (dry air) [J/kg/K]
Definition: scale_const.F90:55
scale_prof::prof_rapend
subroutine, public prof_rapend(rapname_base, level, disable_barrier)
Save raptime.
Definition: scale_prof.F90:217
scale_atmos_dyn_common::atmos_dyn_prep_pres_linearization
subroutine, public atmos_dyn_prep_pres_linearization(DPRES, RT2P, REF_rhot, RHOT, QTRC, REF_pres, AQ_R, AQ_CV, AQ_CP, AQ_MASS)
Definition: scale_atmos_dyn_common.F90:510
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:41
scale_atmos_grid_cartesc_index::jsb
integer, public jsb
Definition: scale_atmos_grid_cartesC_index.F90:65
scale_atmos_dyn_common::atmos_dyn_copy_boundary
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)
Definition: scale_atmos_dyn_common.F90:188
scale_const::const_pre00
real(rp), public const_pre00
pressure reference [Pa]
Definition: scale_const.F90:88
scale_atmos_grid_cartesc_index::je
integer, public je
end point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:56
scale_atmos_grid_cartesc_index::i_xyw
integer, public i_xyw
Definition: scale_atmos_grid_cartesC_index.F90:91