SCALE-RM
scale_atmos_hydrostatic.F90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
10 !-------------------------------------------------------------------------------
11 #include "scalelib.h"
12 module scale_atmos_hydrostatic
13  !-----------------------------------------------------------------------------
14  !
15  !++ used modules
16  !
17  use scale_precision
18  use scale_io
19  use scale_prof
20  use scale_const, only: &
21  grav => const_grav, &
22  rdry => const_rdry, &
23  rvap => const_rvap, &
24  cvdry => const_cvdry, &
25  cvvap => const_cvvap, &
26  cpdry => const_cpdry, &
27  cpvap => const_cpvap, &
28  laps => const_laps, &
29  lapsdry => const_lapsdry, &
30  p00 => const_pre00, &
31  epstvap => const_epstvap
32  !-----------------------------------------------------------------------------
33  implicit none
34  private
35  !-----------------------------------------------------------------------------
36  !
37  !++ Public procedure
38  !
39  public :: atmos_hydrostatic_setup
40  public :: atmos_hydrostatic_buildrho
41  public :: atmos_hydrostatic_buildrho_real
42  public :: atmos_hydrostatic_buildrho_atmos
43  public :: atmos_hydrostatic_buildrho_atmos_rev
44  public :: atmos_hydrostatic_buildrho_bytemp
45  public :: atmos_hydrostatic_buildrho_bytemp_atmos
46 
47  public :: atmos_hydrostatic_barometric_law_mslp
48  public :: atmos_hydrostatic_barometric_law_pres
49 
50  interface atmos_hydrostatic_buildrho
51 #ifdef _OPENACC
52  module procedure atmos_hydrostatic_buildrho_1d_cpu
53 #endif
54  module procedure atmos_hydrostatic_buildrho_1d
55  module procedure atmos_hydrostatic_buildrho_3d
56  end interface atmos_hydrostatic_buildrho
57 
58  interface atmos_hydrostatic_buildrho_real
59  module procedure atmos_hydrostatic_buildrho_real_3d
60  end interface atmos_hydrostatic_buildrho_real
61 
62  interface atmos_hydrostatic_buildrho_atmos
63  module procedure atmos_hydrostatic_buildrho_atmos_0d
64  module procedure atmos_hydrostatic_buildrho_atmos_1d
65  module procedure atmos_hydrostatic_buildrho_atmos_3d
66  end interface atmos_hydrostatic_buildrho_atmos
67 
68  interface atmos_hydrostatic_buildrho_atmos_rev
69  module procedure atmos_hydrostatic_buildrho_atmos_rev_1d
70  module procedure atmos_hydrostatic_buildrho_atmos_rev_2d
71  module procedure atmos_hydrostatic_buildrho_atmos_rev_3d
72  end interface atmos_hydrostatic_buildrho_atmos_rev
73 
74  interface atmos_hydrostatic_buildrho_bytemp
75  module procedure atmos_hydrostatic_buildrho_bytemp_1d
76  module procedure atmos_hydrostatic_buildrho_bytemp_3d
77  end interface atmos_hydrostatic_buildrho_bytemp
78 
79  interface atmos_hydrostatic_buildrho_bytemp_atmos
80  module procedure atmos_hydrostatic_buildrho_bytemp_atmos_1d
81  module procedure atmos_hydrostatic_buildrho_bytemp_atmos_3d
82  end interface atmos_hydrostatic_buildrho_bytemp_atmos
83 
84  interface atmos_hydrostatic_barometric_law_mslp
85  module procedure atmos_hydrostatic_barometric_law_mslp_0d
86  module procedure atmos_hydrostatic_barometric_law_mslp_2d
87  end interface atmos_hydrostatic_barometric_law_mslp
88 
89  interface atmos_hydrostatic_barometric_law_pres
90  module procedure atmos_hydrostatic_barometric_law_pres_0d
91  module procedure atmos_hydrostatic_barometric_law_pres_2d
92  end interface atmos_hydrostatic_barometric_law_pres
93 
94  !-----------------------------------------------------------------------------
95  !
96  !++ Public parameters & variables
97  !
98  !-----------------------------------------------------------------------------
99  !
100  !++ Private procedure
101  !
102  private :: atmos_hydrostatic_buildrho_atmos_2d
103 
104  !-----------------------------------------------------------------------------
105  !
106  !++ Private parameters & variables
107  !
108  integer, private, parameter :: itelim = 100
109  real(RP), private :: criteria
110  logical, private :: HYDROSTATIC_uselapserate = .false.
111  integer, private :: HYDROSTATIC_buildrho_real_kref = 1
112  integer, private :: HYDROSTATIC_barometric_law_mslp_kref = 1
113  !$acc declare create(criteria, HYDROSTATIC_uselapserate, HYDROSTATIC_barometric_law_mslp_kref)
114 
115  !-----------------------------------------------------------------------------
116 contains
117  !-----------------------------------------------------------------------------
119  subroutine atmos_hydrostatic_setup
120  use scale_prc, only: &
121  prc_abort
122  use scale_const, only: &
123  const_eps
124  implicit none
125 
126  namelist / param_atmos_hydrostatic / &
127  hydrostatic_uselapserate, &
128  hydrostatic_buildrho_real_kref, &
129  hydrostatic_barometric_law_mslp_kref
130 
131  integer :: ierr
132  !---------------------------------------------------------------------------
133 
134  log_newline
135  log_info("ATMOS_HYDROSTATIC_setup",*) 'Setup'
136 
137  !--- read namelist
138  rewind(io_fid_conf)
139  read(io_fid_conf,nml=param_atmos_hydrostatic,iostat=ierr)
140  if( ierr < 0 ) then !--- missing
141  log_info("ATMOS_HYDROSTATIC_setup",*) 'Not found namelist. Default used.'
142  elseif( ierr > 0 ) then !--- fatal error
143  log_error("ATMOS_HYDROSTATIC_setup",*) 'Not appropriate names in namelist PARAM_ATMOS_HYDROSTATIC. Check!'
144  call prc_abort
145  endif
146  log_nml(param_atmos_hydrostatic)
147 
148  criteria = sqrt( const_eps )
149 
150  log_newline
151  log_info("ATMOS_HYDROSTATIC_setup",*) 'Use lapse rate for estimation of surface temperature? : ', hydrostatic_uselapserate
152  log_info("ATMOS_HYDROSTATIC_setup",*) 'Buildrho conversion criteria : ', criteria
153 
154  !$acc update device(criteria, HYDROSTATIC_uselapserate, HYDROSTATIC_barometric_law_mslp_kref)
155 
156  return
157  end subroutine atmos_hydrostatic_setup
158 
159  !-----------------------------------------------------------------------------
161 #ifdef _OPENACC
162 !OCL SERIAL
163  subroutine atmos_hydrostatic_buildrho_1d_cpu( &
164  KA, KS, KE, &
165  pott, qv, qc, &
166  pres_sfc, pott_sfc, qv_sfc, qc_sfc, &
167  cz, fz, &
168  dens, temp, pres, &
169  temp_sfc, &
170  converged )
171  implicit none
172 
173  integer, intent(in) :: ka, ks, ke
174  real(rp), intent(in) :: pott(ka)
175  real(rp), intent(in) :: qv (ka)
176  real(rp), intent(in) :: qc (ka)
177  real(rp), intent(in) :: pres_sfc
178  real(rp), intent(in) :: pott_sfc
179  real(rp), intent(in) :: qv_sfc
180  real(rp), intent(in) :: qc_sfc
181  real(rp), intent(in) :: cz (ka)
182  real(rp), intent(in) :: fz (0:ka)
183  real(rp), intent(out) :: dens(ka)
184  real(rp), intent(out) :: temp(ka)
185  real(rp), intent(out) :: pres(ka)
186  real(rp), intent(out) :: temp_sfc
187  logical, intent(out) :: converged
188 
189  real(rp) :: dz(ka)
190  real(rp) :: work1(ka)
191  real(rp) :: work2(ka)
192 
193  call atmos_hydrostatic_buildrho_1d( &
194  ka, ks, ke, &
195  pott, qv, qc, &
196  pres_sfc, pott_sfc, qv_sfc, qc_sfc, &
197  cz, fz, &
198  dz, work1, work2, &
199  dens, temp, pres, &
200  temp_sfc, &
201  converged )
202 
203  return
204  end subroutine atmos_hydrostatic_buildrho_1d_cpu
205 #endif
206 !OCL SERIAL
207  subroutine atmos_hydrostatic_buildrho_1d( &
208  KA, KS, KE, &
209  pott, qv, qc, &
210  pres_sfc, pott_sfc, qv_sfc, qc_sfc, &
211  cz, fz, &
212 #ifdef _OPENACC
213  dz, work1, work2, &
214 #endif
215  dens, temp, pres, &
216  temp_sfc, &
217  converged )
218  !$acc routine seq
219  use scale_prc, only: &
220  prc_abort
221  use scale_atmos_hydrometeor, only: &
222  cv_vapor, &
223  cv_water, &
224  cp_vapor, &
225  cp_water
226  implicit none
227 
228  integer, intent(in) :: KA, KS, KE
229  real(RP), intent(in) :: pott(KA)
230  real(RP), intent(in) :: qv (KA)
231  real(RP), intent(in) :: qc (KA)
232  real(RP), intent(in) :: pres_sfc
233  real(RP), intent(in) :: pott_sfc
234  real(RP), intent(in) :: qv_sfc
235  real(RP), intent(in) :: qc_sfc
236  real(RP), intent(in) :: cz (KA)
237  real(RP), intent(in) :: fz (0:KA)
238  real(RP), intent(out) :: dens(KA)
239  real(RP), intent(out) :: temp(KA)
240  real(RP), intent(out) :: pres(KA)
241  real(RP), intent(out) :: temp_sfc
242  logical, intent(out) :: converged
243 
244 #ifdef _OPENACC
245  real(RP), intent(out) :: dz(KA)
246  real(RP), intent(out) :: work1(KA)
247  real(RP), intent(out) :: work2(KA)
248 #else
249  real(RP) :: dz(KA)
250 #endif
251 
252  real(RP) :: dens_sfc
253  real(RP) :: Rtot_sfc
254  real(RP) :: CVtot_sfc
255  real(RP) :: CPtot_sfc
256  real(RP) :: CPovCV_sfc
257  real(RP) :: Rtot
258  real(RP) :: CVtot
259  real(RP) :: CPtot
260  real(RP) :: CPovCV
261 
262  integer :: k
263  !---------------------------------------------------------------------------
264 
265  !--- from surface to lowermost atmosphere
266 
267  rtot_sfc = rdry * ( 1.0_rp - qv_sfc - qc_sfc ) &
268  + rvap * qv_sfc
269  cvtot_sfc = cvdry * ( 1.0_rp - qv_sfc - qc_sfc ) &
270  + cv_vapor * qv_sfc &
271  + cv_water * qc_sfc
272  cptot_sfc = cpdry * ( 1.0_rp - qv_sfc - qc_sfc ) &
273  + cp_vapor * qv_sfc &
274  + cp_water * qc_sfc
275  cpovcv_sfc = cptot_sfc / cvtot_sfc
276 
277  rtot = rdry * ( 1.0_rp - qv(ks) - qc(ks) ) &
278  + rvap * qv(ks)
279  cvtot = cvdry * ( 1.0_rp - qv(ks) - qc(ks) ) &
280  + cv_vapor * qv(ks) &
281  + cv_water * qc(ks)
282  cptot = cpdry * ( 1.0_rp - qv(ks) - qc(ks) ) &
283  + cp_vapor * qv(ks) &
284  + cp_water * qc(ks)
285  cpovcv = cptot / cvtot
286 
287  ! density at surface
288  dens_sfc = p00 / rtot_sfc / pott_sfc * ( pres_sfc/p00 )**(cvtot_sfc/cptot_sfc)
289  temp_sfc = pres_sfc / ( dens_sfc * rtot_sfc )
290 
291  dz(ks-1) = cz(ks) - fz(ks-1)
292  do k = ks, ke-1
293  dz(k) = cz(k+1) - cz(k)
294  end do
295 
296  ! make density at lowermost cell center
297  if ( hydrostatic_uselapserate ) then
298 
299  temp(ks) = pott_sfc - lapsdry * dz(ks-1) ! use dry lapse rate
300  pres(ks) = p00 * ( temp(ks)/pott(ks) )**(cptot/rtot)
301  dens(ks) = p00 / rtot / pott(ks) * ( pres(ks)/p00 )**(cvtot/cptot)
302 
303  converged = .true.
304 
305  else ! use itelation
306 
307  call atmos_hydrostatic_buildrho_atmos_0d( pott(ks), qv(ks), qc(ks), & ! [IN]
308  dens_sfc, pott_sfc, qv_sfc, qc_sfc, & ! [IN]
309  dz(ks-1), ks-1, & ! [IN]
310  dens(ks), temp(ks), pres(ks), & ! [OUT]
311  converged ) ! [OUT]
312 
313  endif
314 
315  if ( converged ) then
316  !--- from lowermost atmosphere to top of atmosphere
317  call atmos_hydrostatic_buildrho_atmos_1d( ka, ks, ke, &
318  pott(:), qv(:), qc(:), & ! [IN]
319  dz(:), & ! [IN]
320 #ifdef _OPENACC
321  work1(:), work2(:), & ! [WORK]
322 #endif
323  dens(:), & ! [INOUT]
324  temp(:), pres(:), & ! [OUT]
325  converged ) ! [OUT]
326  if ( converged ) then
327  ! fill dummy
328  dens( 1:ks-1) = 0.0_rp
329  dens(ke+1:ka ) = 0.0_rp
330  end if
331  end if
332 
333  return
334  end subroutine atmos_hydrostatic_buildrho_1d
335 
336  !-----------------------------------------------------------------------------
338  subroutine atmos_hydrostatic_buildrho_3d( &
339  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
340  pott, qv, qc, &
341  pres_sfc, pott_sfc, &
342  qv_sfc, qc_sfc, &
343  cz, fz, area, &
344  dens, temp, pres, &
345  temp_sfc )
346  use scale_prc, only: &
347  prc_abort
348  use scale_statistics, only: &
349  statistics_horizontal_mean
350  implicit none
351 
352  integer, intent(in) :: KA, KS, KE
353  integer, intent(in) :: IA, IS, IE
354  integer, intent(in) :: JA, JS, JE
355 
356  real(RP), intent(in) :: pott(KA,IA,JA)
357  real(RP), intent(in) :: qv (KA,IA,JA)
358  real(RP), intent(in) :: qc (KA,IA,JA)
359  real(RP), intent(in) :: pres_sfc(IA,JA)
360  real(RP), intent(in) :: pott_sfc(IA,JA)
361  real(RP), intent(in) :: qv_sfc (IA,JA)
362  real(RP), intent(in) :: qc_sfc (IA,JA)
363  real(RP), intent(in) :: cz( KA,IA,JA)
364  real(RP), intent(in) :: fz(0:KA,IA,JA)
365  real(RP), intent(in) :: area(IA,JA)
366 
367  real(RP), intent(out) :: dens(KA,IA,JA)
368  real(RP), intent(out) :: temp(KA,IA,JA)
369  real(RP), intent(out) :: pres(KA,IA,JA)
370  real(RP), intent(out) :: temp_sfc(IA,JA)
371 
372  real(RP) :: dz(KA,IA,JA), dz_top(IA,JA)
373 
374  ! TOA
375  real(RP) :: pott_toa(IA,JA)
376  real(RP) :: qv_toa (IA,JA)
377  real(RP) :: qc_toa (IA,JA)
378  real(RP) :: dens_toa(IA,JA)
379  real(RP) :: temp_toa(IA,JA)
380  real(RP) :: pres_toa(IA,JA)
381  ! k = KE
382  real(RP) :: pott_ke(IA,JA)
383  real(RP) :: qv_ke (IA,JA)
384  real(RP) :: qc_ke (IA,JA)
385  real(RP) :: dens_ke(IA,JA)
386  real(RP) :: temp_ke(IA,JA)
387  real(RP) :: pres_ke(IA,JA)
388 
389  real(RP) :: dens_mean
390 
391 #ifdef _OPENACC
392  real(RP) :: work1(KA,IA,JA)
393  real(RP) :: work2(KA,IA,JA)
394  real(RP) :: work3(KA,IA,JA)
395 #endif
396 
397  logical :: converged, flag
398  integer :: k, i, j
399  !---------------------------------------------------------------------------
400 
401  !$acc data copyin(pott, qv, qc, pres_sfc, pott_sfc, qv_sfc, qc_sfc, cz, fz, area) &
402  !$acc copyout(dens, temp, pres, temp_sfc) &
403  !$acc create(dz, dz_top, pott_toa, qv_toa, qc_toa, dens_toa, temp_toa, pres_toa, pott_ke, qv_ke, qc_ke, dens_ke, temp_ke, pres_ke)
404 
405 
406  converged = .true.
407 
408  !--- from surface to lowermost atmosphere
409  !$omp parallel do OMP_SCHEDULE_ collapse(2) private(flag) reduction(.and.:converged)
410  !$acc kernels
411  !$acc loop independent collapse(2) &
412  !$acc private(flag) reduction(.and.: converged)
413  do j = js, je
414  do i = is, ie
415  call atmos_hydrostatic_buildrho_1d( ka, ks, ke, &
416  pott(:,i,j), qv(:,i,j), qc(:,i,j), & ! [IN]
417  pres_sfc(i,j), pott_sfc(i,j), qv_sfc(i,j), qc_sfc(i,j), & ! [IN]
418  cz(:,i,j), fz(:,i,j), & ! [IN]
419 #ifdef _OPENACC
420  work1(:,i,j), work2(:,i,j), work3(:,i,j), & ! [WORK]
421 #endif
422  dens(:,i,j), temp(:,i,j), pres(:,i,j), & ! [OUT]
423  temp_sfc(i,j), & ! [OUT]
424  flag ) ! [OUT]
425  converged = converged .and. flag
426  end do
427  end do
428  !$acc end kernels
429 
430  if ( .not. converged ) then
431  log_error("ATMOS_HYDROSTATIC_buildrho_3D",*) "not converged"
432  call prc_abort
433  end if
434 
435  !$omp parallel do OMP_SCHEDULE_ collapse(2)
436  !$acc kernels
437  do j = js, je
438  do i = is, ie
439  do k = ks, ke-1
440  dz(k,i,j) = cz(k+1,i,j) - cz(k,i,j)
441  end do
442  dz_top(i,j) = fz(ke,i,j) - cz(ke,i,j) ! distance from cell center to TOA
443 
444  ! value at TOA
445  pott_toa(i,j) = pott(ke,i,j)
446  qv_toa(i,j) = qv(ke,i,j)
447  qc_toa(i,j) = qc(ke,i,j)
448  ! value at k=KE
449  dens_ke(i,j) = dens(ke,i,j)
450  pott_ke(i,j) = pott(ke,i,j)
451  qv_ke(i,j) = qv(ke,i,j)
452  qc_ke(i,j) = qc(ke,i,j)
453  end do
454  end do
455  !$acc end kernels
456 
457  call atmos_hydrostatic_buildrho_atmos_2d( ia, is, ie, ja, js, je, &
458  pott_toa(:,:), qv_toa(:,:), qc_toa(:,:), & ! [IN]
459  dens_ke(:,:), pott_ke(:,:), qv_ke(:,:), qc_ke(:,:), & ! [IN]
460  dz_top(:,:), ke+1, & ! [IN]
461  dens_toa(:,:), temp_toa(:,:), pres_toa(:,:) ) ! [OUT]
462 
463  call statistics_horizontal_mean( ia, is, ie, ja, js, je, &
464  dens_toa(:,:), area(:,:), & ! [IN]
465  dens_mean ) ! [OUT]
466 
467  !$omp parallel do OMP_SCHEDULE_ collapse(2)
468  !$acc kernels
469  do j = js, je
470  do i = is, ie
471  dens_toa(i,j) = dens_mean
472  enddo
473  enddo
474  !$acc end kernels
475 
476  call atmos_hydrostatic_buildrho_atmos_rev_2d( ia, is, ie, ja, js, je, &
477  pott_ke(:,:), qv_ke(:,:), qc_ke(:,:), & ! [IN]
478  dens_toa(:,:), pott_toa(:,:), qv_toa(:,:), qc_toa(:,:), & ! [IN]
479  dz_top(:,:), ke+1, & ! [IN]
480  dens_ke(:,:), temp_ke(:,:), pres_ke(:,:) ) ! [OUT]
481 
482  !$omp parallel do OMP_SCHEDULE_
483  !$acc kernels
484  do j = js, je
485  do i = is, ie
486  dens(ke,i,j) = dens_ke(i,j)
487  end do
488  end do
489  !$acc end kernels
490 
491  !--- from top of atmosphere to lowermost atmosphere
492  call atmos_hydrostatic_buildrho_atmos_rev_3d( ka, ks, ke, ia, is, ie, ja, js, je, &
493  pott(:,:,:), qv(:,:,:), qc(:,:,:), & ! [IN]
494  dz(:,:,:), & ! [IN]
495  dens(:,:,:), & ! [INOUT]
496  temp(:,:,:), pres(:,:,:) ) ! [OUT]
497 
498  !$acc end data
499 
500  return
501  end subroutine atmos_hydrostatic_buildrho_3d
502 
503  !-----------------------------------------------------------------------------
505  subroutine atmos_hydrostatic_buildrho_real_3d( &
506  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
507  pott, qv, qc, &
508  cz, &
509  pres, &
510  dens, temp )
511  use scale_atmos_hydrometeor, only: &
512  cv_vapor, &
513  cv_water, &
514  cp_vapor, &
515  cp_water
516  use scale_const, only: &
517  undef => const_undef
518  implicit none
519  integer, intent(in) :: KA, KS, KE
520  integer, intent(in) :: IA, IS, IE
521  integer, intent(in) :: JA, JS, JE
522 
523  real(RP), intent(in) :: pott(KA,IA,JA)
524  real(RP), intent(in) :: qv (KA,IA,JA)
525  real(RP), intent(in) :: qc (KA,IA,JA)
526  real(RP), intent(in) :: cz (KA,IA,JA)
527 
528  real(RP), intent(inout) :: pres(KA,IA,JA)
529 
530  real(RP), intent(out) :: dens(KA,IA,JA)
531  real(RP), intent(out) :: temp(KA,IA,JA)
532 
533  real(RP) :: dz(KA,IA,JA)
534 
535  real(RP) :: pott_toa(IA,JA)
536  real(RP) :: qv_toa (IA,JA)
537  real(RP) :: qc_toa (IA,JA)
538 
539  real(RP) :: Rtot
540  real(RP) :: CVtot
541  real(RP) :: CPtot
542 
543  integer :: kref(IA,JA)
544  integer :: k, i, j
545  !---------------------------------------------------------------------------
546 
547  !$acc data copyin(pott, qv, qc, cz) copy(pres), copyout(dens, temp) create(dz, pott_toa, qv_toa, qc_toa, kref)
548 
549  !--- from surface to lowermost atmosphere
550 
551  !$omp parallel do OMP_SCHEDULE_ collapse(2)
552  !$acc kernels
553  do j = js, je
554  do i = is, ie
555  do k = ks, ke-1
556  dz(k,i,j) = cz(k+1,i,j) - cz(k,i,j) ! distance from cell center to cell center
557  enddo
558  enddo
559  enddo
560  !$acc end kernels
561 
562  !$omp parallel do OMP_SCHEDULE_ collapse(2)
563  !$acc kernels
564  do j = js, je
565  do i = is, ie
566  pott_toa(i,j) = pott(ke,i,j)
567  qv_toa(i,j) = qv(ke,i,j)
568  qc_toa(i,j) = qc(ke,i,j)
569  enddo
570  enddo
571  !$acc end kernels
572 
573  !$omp parallel do OMP_SCHEDULE_ collapse(2)
574  !$acc kernels
575  do j = js, je
576  do i = is, ie
577  kref(i,j) = hydrostatic_buildrho_real_kref + ks - 1
578  !$acc loop seq
579  do k = kref(i,j), ke
580  if ( pres(k,i,j) .ne. undef ) then
581  kref(i,j) = k
582  exit
583  end if
584  end do
585  end do
586  end do
587  !$acc end kernels
588 
589  ! calc density at reference level
590  !$omp parallel do OMP_SCHEDULE_ collapse(2) &
591  !$omp private(Rtot,CVtot,CPtot,k)
592  !$acc kernels
593  do j = js, je
594  do i = is, ie
595  k = kref(i,j)
596  rtot = rdry * ( 1.0_rp - qv(k,i,j) - qc(k,i,j) ) &
597  + rvap * qv(k,i,j)
598  cvtot = cvdry * ( 1.0_rp - qv(k,i,j) - qc(k,i,j) ) &
599  + cv_vapor * qv(k,i,j) &
600  + cv_water * qc(k,i,j)
601  cptot = cpdry * ( 1.0_rp - qv(k,i,j) - qc(k,i,j) ) &
602  + cp_vapor * qv(k,i,j) &
603  + cp_water * qc(k,i,j)
604  dens(k,i,j) = p00 / ( rtot * pott(k,i,j) ) * ( pres(k,i,j)/p00 )**(cvtot/cptot)
605  enddo
606  enddo
607  !$acc end kernels
608 
609  !--- from lowermost atmosphere to top of atmosphere
610  call atmos_hydrostatic_buildrho_atmos_3d( ka, ks, ke, ia, is, ie, ja, js, je, &
611  pott(:,:,:), qv(:,:,:), qc(:,:,:), & ! [IN]
612  dz(:,:,:), & ! [IN]
613  dens(:,:,:), & ! [INOUT]
614  temp(:,:,:), pres(:,:,:), & ! [OUT]
615  kref = kref ) ! [IN]
616  call atmos_hydrostatic_buildrho_atmos_rev_3d( ka, ks, ke, ia, is, ie, ja, js, je, &
617  pott(:,:,:), qv(:,:,:), qc(:,:,:), & ! [IN]
618  dz(:,:,:), & ! [IN]
619  dens(:,:,:), & ! [INOUT]
620  temp(:,:,:), pres(:,:,:), & ! [OUT]
621  kref = kref ) ! [IN]
622 
623 !!$ call ATMOS_HYDROSTATIC_buildrho_atmos_2D( IA, IS, IE, JA, JS, JE, &
624 !!$ pott_toa(:,:), qv_toa(:,:), qc_toa(:,:), & ! [IN]
625 !!$ dens(KE,:,:), & ! [IN]
626 !!$ pott(KE,:,:), qv(KE,:,:), qc(KE,:,:), & ! [IN]
627 !!$ dz(KE+1,:,:), KE+1, & ! [IN]
628 !!$ dens_toa(:,:), temp_toa(:,:), pres_toa(:,:) ) ! [OUT]
629 
630  ! density at TOA
631  !$acc kernels
632  dens( 1:ks-1,:,:) = 0.0_rp ! fill dummy
633 !!$ dens(KE+2:KA ,:,:) = 0.0_RP ! fill dummy
634  dens(ke+1:ka ,:,:) = 0.0_rp ! fill dummy
635  !$acc end kernels
636 
637  !$acc end data
638 
639  return
640  end subroutine atmos_hydrostatic_buildrho_real_3d
641 
642  !-----------------------------------------------------------------------------
644 !OCL SERIAL
645 !OCL NOSIMD
646  subroutine atmos_hydrostatic_buildrho_atmos_0d( &
647  pott_L2, qv_L2, qc_L2, &
648  dens_L1, pott_L1, qv_L1, qc_L1, &
649  dz, k, &
650  dens_L2, temp_L2, pres_L2, &
651  converged )
652  use scale_prc, only: &
653  prc_abort
654  use scale_atmos_hydrometeor, only: &
655  cv_vapor, &
656  cv_water, &
657  cp_vapor, &
658  cp_water
659  !$acc routine seq
660  implicit none
661  real(RP), intent(in) :: pott_L2
662  real(RP), intent(in) :: qv_L2
663  real(RP), intent(in) :: qc_L2
664  real(RP), intent(in) :: dens_L1
665  real(RP), intent(in) :: pott_L1
666  real(RP), intent(in) :: qv_L1
667  real(RP), intent(in) :: qc_L1
668  real(RP), intent(in) :: dz
669  integer, intent(in) :: k
670 
671  real(RP), intent(out) :: dens_L2
672  real(RP), intent(out) :: temp_L2
673  real(RP), intent(out) :: pres_L2
674  logical, intent(out) :: converged
675 
676  real(RP) :: Rtot_L1 , Rtot_L2
677  real(RP) :: CVtot_L1 , CVtot_L2
678  real(RP) :: CPtot_L1 , CPtot_L2
679  real(RP) :: CPovCV_L1, CPovCV_L2
680 
681  real(RP) :: pres_L1
682  real(RP) :: dens_s, dhyd, dgrd
683  integer :: ite
684  !---------------------------------------------------------------------------
685 
686  rtot_l1 = rdry * ( 1.0_rp - qv_l1 - qc_l1 ) &
687  + rvap * qv_l1
688  cvtot_l1 = cvdry * ( 1.0_rp - qv_l1 - qc_l1 ) &
689  + cv_vapor * qv_l1 &
690  + cv_water * qc_l1
691  cptot_l1 = cpdry * ( 1.0_rp - qv_l1 - qc_l1 ) &
692  + cp_vapor * qv_l1 &
693  + cp_water * qc_l1
694  cpovcv_l1 = cptot_l1 / cvtot_l1
695 
696  rtot_l2 = rdry * ( 1.0_rp - qv_l2 - qc_l2 ) &
697  + rvap * qv_l2
698  cvtot_l2 = cvdry * ( 1.0_rp - qv_l2 - qc_l2 ) &
699  + cv_vapor * qv_l2 &
700  + cv_water * qc_l2
701  cptot_l2 = cpdry * ( 1.0_rp - qv_l2 - qc_l2 ) &
702  + cp_vapor * qv_l2 &
703  + cp_water * qc_l2
704  cpovcv_l2 = cptot_l2 / cvtot_l2
705 
706  dens_s = 0.0_rp
707  dens_l2 = dens_l1 ! first guess
708 
709  pres_l1 = p00 * ( dens_l1 * rtot_l1 * pott_l1 / p00 )**cpovcv_l1
710 
711  converged = .false.
712  do ite = 1, itelim
713  if ( abs(dens_l2-dens_s) <= criteria ) then
714  converged = .true.
715  exit
716  endif
717 
718  dens_s = dens_l2
719  pres_l2 = p00 * ( dens_s * rtot_l2 * pott_l2 / p00 )**cpovcv_l2
720 
721  dhyd = + ( pres_l1 - pres_l2 ) / dz & ! dp/dz
722  - grav * 0.5_rp * ( dens_l1 + dens_s ) ! rho*g
723 
724  dgrd = - pres_l2 * cpovcv_l2 / dens_s / dz &
725  - grav * 0.5_rp
726 
727  dens_l2 = max( dens_s - dhyd/dgrd, dens_s * 0.1_rp )
728 
729  if ( dens_l2*0.0_rp /= 0.0_rp ) exit
730  enddo
731 
732  if ( converged ) then
733  pres_l2 = p00 * ( dens_l2 * rtot_l2 * pott_l2 / p00 )**cpovcv_l2
734  temp_l2 = pres_l2 / ( dens_l2 * rtot_l2 )
735 #ifndef _OPENACC
736  else
737  log_error("ATMOS_HYDROSTATIC_buildrho_atmos_0D",*) 'iteration not converged!', &
738  k,dens_l2,ite,dens_s,dhyd,dgrd
739  call prc_abort
740 #endif
741  endif
742 
743  return
744  end subroutine atmos_hydrostatic_buildrho_atmos_0d
745 
746  !-----------------------------------------------------------------------------
748 !OCL SERIAL
749  subroutine atmos_hydrostatic_buildrho_atmos_1d( &
750  KA, KS, KE, &
751  pott, qv, qc, &
752  dz, &
753 #ifdef _OPENACC
754  Rtot, CPovCV, &
755 #endif
756  dens, &
757  temp, pres, &
758  converged )
759  !$acc routine seq
760  use scale_const, only: &
761  undef => const_undef, &
762  eps => const_eps
763  use scale_prc, only: &
764  prc_abort
765  use scale_atmos_hydrometeor, only: &
766  cv_vapor, &
767  cv_water, &
768  cp_vapor, &
769  cp_water
770  implicit none
771  integer, intent(in) :: KA, KS, KE
772 
773  real(RP), intent(in) :: pott(KA)
774  real(RP), intent(in) :: qv (KA)
775  real(RP), intent(in) :: qc (KA)
776  real(RP), intent(in) :: dz (KA)
777 
778  real(RP), intent(inout) :: dens(KA)
779 
780  real(RP), intent(inout) :: temp(KA)
781  real(RP), intent(inout) :: pres(KA)
782  logical, intent(out) :: converged
783 
784 #ifdef _OPENACC
785  real(RP), intent(out) :: Rtot (KA)
786  real(RP), intent(out) :: CPovCV(KA)
787 #else
788  real(RP) :: Rtot (KA)
789  real(RP) :: CPovCV(KA)
790 #endif
791  real(RP) :: CVtot
792  real(RP) :: CPtot
793 
794  real(RP) :: dens_s, dhyd, dgrd
795  integer :: ite
796 
797  integer :: k
798  !---------------------------------------------------------------------------
799 
800  do k = ks, ke
801  rtot(k) = rdry * ( 1.0_rp - qv(k) - qc(k) ) &
802  + rvap * qv(k)
803  cvtot = cvdry * ( 1.0_rp - qv(k) - qc(k) ) &
804  + cv_vapor * qv(k) &
805  + cv_water * qc(k)
806  cptot = cpdry * ( 1.0_rp - qv(k) - qc(k) ) &
807  + cp_vapor * qv(k) &
808  + cp_water * qc(k)
809  cpovcv(k) = cptot / cvtot
810  enddo
811 
812  pres(ks) = p00 * ( dens(ks) * rtot(ks) * pott(ks) / p00 )**cpovcv(ks)
813 
814  do k = ks+1, ke
815 
816  dens_s = 0.0_rp
817  dens(k) = dens(k-1) ! first guess
818 
819  converged = .false.
820  do ite = 1, itelim
821  if ( abs(dens(k)-dens_s) <= criteria ) then
822  converged = .true.
823  exit
824  endif
825 
826  dens_s = dens(k)
827  pres(k) = p00 * ( dens_s * rtot(k) * pott(k) / p00 )**cpovcv(k)
828 
829  dhyd = + ( pres(k-1) - pres(k) ) / dz(k-1) & ! dpdz
830  - grav * 0.5_rp * ( dens(k-1) + dens_s ) ! rho*g
831 
832  dgrd = - cpovcv(k) * pres(k) / dens_s / dz(k-1) &
833  - grav * 0.5_rp
834 
835  dens(k) = max( dens_s - dhyd/dgrd, dens(k-1) * eps )
836 
837  if ( dens(k)*0.0_rp /= 0.0_rp ) exit
838  enddo
839 
840  if ( .NOT. converged ) then
841 #ifndef _OPENACC
842  log_error("ATMOS_HYDROSTATIC_buildrho_atmos_1D",*) 'iteration not converged!', &
843  k,dens(k),ite,dens_s,dhyd,dgrd
844  call prc_abort
845 #endif
846  exit
847  endif
848 
849  pres(k) = p00 * ( dens(k) * rtot(k) * pott(k) / p00 )**cpovcv(k)
850 
851  enddo
852 
853  if ( converged ) then
854  do k = ks, ke
855  temp(k) = pres(k) / ( dens(k) * rtot(k) )
856  enddo
857 
858  dens(ke+1:ka ) = dens(ke)
859  pres(ke+1:ka ) = pres(ke)
860  temp(ke+1:ka ) = temp(ke)
861  end if
862 
863  return
864  end subroutine atmos_hydrostatic_buildrho_atmos_1d
865 
866  !-----------------------------------------------------------------------------
868 !OCL SERIAL
869  subroutine atmos_hydrostatic_buildrho_atmos_rev_1d( &
870  KA, KS, KE, &
871  pott, qv, qc, &
872  dz, &
873 #ifdef _OPENACC
874  Rtot, CPovCV, &
875 #endif
876  dens, &
877  temp, pres, &
878  converged )
879  !$acc routine seq
880  use scale_prc, only: &
881  prc_abort
882  use scale_atmos_hydrometeor, only: &
883  cv_vapor, &
884  cv_water, &
885  cp_vapor, &
886  cp_water
887  implicit none
888  integer, intent(in) :: KA, KS, KE
889 
890  real(RP), intent(in) :: pott(KA)
891  real(RP), intent(in) :: qv (KA)
892  real(RP), intent(in) :: qc (KA)
893  real(RP), intent(in) :: dz (KA)
894 
895  real(RP), intent(inout) :: dens(KA)
896 
897  real(RP), intent(inout) :: temp(KA)
898  real(RP), intent(inout) :: pres(KA)
899  logical, intent(out) :: converged
900 
901 #ifdef _OPENACC
902  real(RP), intent(out) :: Rtot (KA)
903  real(RP), intent(out) :: CPovCV(KA)
904 #else
905  real(RP) :: Rtot (KA)
906  real(RP) :: CPovCV(KA)
907 #endif
908  real(RP) :: CVtot
909  real(RP) :: CPtot
910 
911  real(RP) :: dens_s, dhyd, dgrd
912  integer :: ite
913 
914  integer :: k
915  !---------------------------------------------------------------------------
916 
917  do k = ks, ke
918  rtot(k) = rdry * ( 1.0_rp - qv(k) - qc(k) ) &
919  + rvap * qv(k)
920  cvtot = cvdry * ( 1.0_rp - qv(k) - qc(k) ) &
921  + cv_vapor * qv(k) &
922  + cv_water * qc(k)
923  cptot = cpdry * ( 1.0_rp - qv(k) - qc(k) ) &
924  + cp_vapor * qv(k) &
925  + cp_water * qc(k)
926  cpovcv(k) = cptot / cvtot
927  enddo
928 
929  pres(ke) = p00 * ( dens(ke) * rtot(ke) * pott(ke) / p00 )**cpovcv(ke)
930 
931  do k = ke-1, ks, -1
932 
933  dens_s = 0.0_rp
934  dens(k) = dens(k+1) ! first guess
935 
936  converged = .false.
937  do ite = 1, itelim
938  if ( abs(dens(k)-dens_s) <= criteria ) then
939  converged = .true.
940  exit
941  endif
942 
943  dens_s = dens(k)
944  pres(k) = p00 * ( dens_s * rtot(k) * pott(k) / p00 )**cpovcv(k)
945 
946  dhyd = - ( pres(k+1) - pres(k) ) / dz(k) & ! dpdz
947  - grav * 0.5_rp * ( dens(k+1) + dens_s ) ! rho*g
948 
949  dgrd = + cpovcv(k) * pres(k) / dens_s / dz(k) &
950  - grav * 0.5_rp
951 
952  dens(k) = dens_s - dhyd/dgrd
953 
954  if ( dens(k)*0.0_rp /= 0.0_rp ) exit
955  enddo
956 
957  if ( .NOT. converged ) then
958 #ifndef _OPENACC
959  log_error("ATMOS_HYDROSTATIC_buildrho_atmos_rev_1D",*) 'iteration not converged!', &
960  k,dens(k),ite,dens_s,dhyd,dgrd
961  call prc_abort
962 #endif
963  exit
964  endif
965 
966  pres(k) = p00 * ( dens(k) * rtot(k) * pott(k) / p00 )**cpovcv(k)
967 
968  enddo
969 
970  if ( converged ) then
971  do k = ks, ke
972  temp(k) = pres(k) / ( dens(k) * rtot(k) )
973  enddo
974 
975 !!$ dens( 1:KS-1) = dens(KS)
976 !!$ pres( 1:KS-1) = pres(KS)
977 !!$ temp( 1:KS-1) = temp(KS)
978  end if
979 
980  return
981  end subroutine atmos_hydrostatic_buildrho_atmos_rev_1d
982 
983  !-----------------------------------------------------------------------------
985  subroutine atmos_hydrostatic_buildrho_atmos_2d( &
986  IA, IS, IE, JA, JS, JE, &
987  pott_L2, qv_L2, qc_L2, &
988  dens_L1, pott_L1, qv_L1, qc_L1, &
989  dz, k, &
990  dens_L2, temp_L2, pres_L2 )
991  use scale_prc, only: &
992  prc_abort
993  implicit none
994  integer, intent(in) :: IA, IS, IE
995  integer, intent(in) :: JA, JS, JE
996 
997  real(RP), intent(in) :: pott_L2(IA,JA)
998  real(RP), intent(in) :: qv_L2 (IA,JA)
999  real(RP), intent(in) :: qc_L2 (IA,JA)
1000  real(RP), intent(in) :: dens_L1(IA,JA)
1001  real(RP), intent(in) :: pott_L1(IA,JA)
1002  real(RP), intent(in) :: qv_L1 (IA,JA)
1003  real(RP), intent(in) :: qc_L1 (IA,JA)
1004  real(RP), intent(in) :: dz (IA,JA)
1005  integer, intent(in) :: k
1006 
1007  real(RP), intent(out) :: dens_L2(IA,JA)
1008  real(RP), intent(out) :: temp_L2(IA,JA)
1009  real(RP), intent(out) :: pres_L2(IA,JA)
1010 
1011  logical :: converged, flag
1012  integer :: i, j
1013  !---------------------------------------------------------------------------
1014 
1015  converged = .true.
1016 
1017  !$omp parallel do OMP_SCHEDULE_ collapse(2) private(flag) reduction(.and.:converged)
1018  !$acc kernels copyin(pott_L2, qv_l2, qc_l2, dens_l1, pott_L1, qv_L1, qc_L1, dz) copyout(dens_L2, temp_L2, pres_L2)
1019  !$acc loop independent collapse(2) &
1020  !$acc private(flag) reduction(.and.: converged)
1021  do j = js, je
1022  do i = is, ie
1023  call atmos_hydrostatic_buildrho_atmos_0d( &
1024  pott_l2(i,j), qv_l2(i,j), qc_l2(i,j), & ! [IN]
1025  dens_l1(i,j), pott_l1(i,j), qv_l1(i,j), qc_l1(i,j), & ! [IN]
1026  dz(i,j), k, & ! [IN]
1027  dens_l2(i,j), temp_l2(i,j), pres_l2(i,j), & ! [OUT]
1028  flag ) ! [OUT]
1029  converged = converged .and. flag
1030  enddo
1031  enddo
1032  !$acc end kernels
1033 
1034  if ( .not. converged ) then
1035  log_error("ATMOS_HYDROSTATIC_buildrho_atmos_2D",*) "not converged"
1036  call prc_abort
1037  end if
1038 
1039  return
1040  end subroutine atmos_hydrostatic_buildrho_atmos_2d
1041 
1042  !-----------------------------------------------------------------------------
1044  subroutine atmos_hydrostatic_buildrho_atmos_rev_2d( &
1045  IA, IS, IE, JA, JS, JE, &
1046  pott_L1, qv_L1, qc_L1, &
1047  dens_L2, pott_L2, qv_L2, qc_L2, &
1048  dz, k, &
1049  dens_L1, temp_L1, pres_L1 )
1050  use scale_prc, only: &
1051  prc_abort
1052  implicit none
1053  integer, intent(in) :: IA, IS, IE
1054  integer, intent(in) :: JA, JS, JE
1055 
1056  real(RP), intent(in) :: pott_L1(IA,JA)
1057  real(RP), intent(in) :: qv_L1 (IA,JA)
1058  real(RP), intent(in) :: qc_L1 (IA,JA)
1059  real(RP), intent(in) :: dens_L2(IA,JA)
1060  real(RP), intent(in) :: pott_L2(IA,JA)
1061  real(RP), intent(in) :: qv_L2 (IA,JA)
1062  real(RP), intent(in) :: qc_L2 (IA,JA)
1063  real(RP), intent(in) :: dz (IA,JA)
1064  integer, intent(in) :: k
1065 
1066  real(RP), intent(out) :: dens_L1(IA,JA)
1067  real(RP), intent(out) :: temp_L1(IA,JA)
1068  real(RP), intent(out) :: pres_L1(IA,JA)
1069 
1070  logical :: converged, flag
1071  integer :: i, j
1072  !---------------------------------------------------------------------------
1073 
1074  converged = .true.
1075 
1076  !$omp parallel do OMP_SCHEDULE_ collapse(2) private(flag) reduction(.and.:converged)
1077  !$acc kernels copyin(pott_L1, qv_L1, qc_L1, dens_L2, pott_L2, qv_L2, qc_L2, dz) copyout(dens_L1, temp_L1, pres_L1)
1078  !$acc loop independent collapse(2) &
1079  !$acc private(flag) reduction(.and.: converged)
1080  do j = js, je
1081  do i = is, ie
1082  call atmos_hydrostatic_buildrho_atmos_0d( &
1083  pott_l1(i,j), qv_l1(i,j), qc_l1(i,j), & ! [IN]
1084  dens_l2(i,j), pott_l2(i,j), qv_l2(i,j), qc_l2(i,j), & ! [IN]
1085  -dz(i,j), k, & ! [IN]
1086  dens_l1(i,j), temp_l1(i,j), pres_l1(i,j), & ! [OUT]
1087  flag ) ! [OUT]
1088  converged = converged .and. flag
1089  enddo
1090  enddo
1091  !$acc end kernels
1092 
1093  if ( .not. converged ) then
1094  log_error("ATMOS_HYDROSTATIC_buildrho_atmos_rev_2D",*) "not converged"
1095  call prc_abort
1096  end if
1097 
1098  return
1099  end subroutine atmos_hydrostatic_buildrho_atmos_rev_2d
1100 
1101  !-----------------------------------------------------------------------------
1103  subroutine atmos_hydrostatic_buildrho_atmos_3d( &
1104  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1105  pott, qv, qc, &
1106  dz, &
1107  dens, &
1108  temp, pres, &
1109  kref )
1110  use scale_prc, only: &
1111  prc_abort
1112  implicit none
1113  integer, intent(in) :: KA, KS, KE
1114  integer, intent(in) :: IA, IS, IE
1115  integer, intent(in) :: JA, JS, JE
1116 
1117  real(RP), intent(in) :: pott(KA,IA,JA)
1118  real(RP), intent(in) :: qv (KA,IA,JA)
1119  real(RP), intent(in) :: qc (KA,IA,JA)
1120  real(RP), intent(in) :: dz (KA,IA,JA)
1121 
1122  real(RP), intent(inout) :: dens(KA,IA,JA)
1123 
1124  real(RP), intent(inout) :: temp(KA,IA,JA)
1125  real(RP), intent(inout) :: pres(KA,IA,JA)
1126 
1127  integer, intent(in), optional, target :: kref(IA,JA)
1128 
1129  integer, pointer :: kref_(:,:)
1130 
1131 #ifdef _OPENACC
1132  real(RP) :: work1(KA,IA,JA)
1133  real(RP) :: work2(KA,IA,JA)
1134 #endif
1135 
1136  logical :: converged, flag
1137  integer :: i, j
1138  !---------------------------------------------------------------------------
1139 
1140  if ( present(kref) ) then
1141  kref_ => kref
1142  else
1143  allocate( kref_(ia,ja) )
1144  !$acc enter data create(kref_)
1145  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1146  !$acc kernels
1147  do j = js, je
1148  do i = is, ie
1149  kref_(i,j) = ks
1150  end do
1151  end do
1152  !$acc end kernels
1153  end if
1154 
1155  converged = .true.
1156 
1157  !$omp parallel do OMP_SCHEDULE_ collapse(2) private(flag) reduction(.and.:converged)
1158  !$acc kernels copyin(pott, qv, qc, dz, kref_) copy(dens) copyout(temp, pres)
1159  !$acc loop independent collapse(2) &
1160  !$acc private(flag) reduction(.and.: converged)
1161  do j = js, je
1162  do i = is, ie
1163  if ( kref_(i,j) < ke ) then
1164  call atmos_hydrostatic_buildrho_atmos_1d( ka, kref_(i,j), ke, &
1165  pott(:,i,j), qv(:,i,j), qc(:,i,j), & ! [IN]
1166  dz(:,i,j), & ! [IN]
1167 #ifdef _OPENACC
1168  work1(:,i,j), work2(:,i,j), & ! [WORK]
1169 #endif
1170  dens(:,i,j), & ! [INOUT]
1171  temp(:,i,j), pres(:,i,j), & ! [OUT]
1172  flag ) ! [OUT]
1173  converged = converged .and. flag
1174  end if
1175  enddo
1176  enddo
1177  !$acc end kernels
1178 
1179  if ( .not. converged ) then
1180  log_error("ATMOS_HYDROSTATIC_buildrho_atmos_3D",*) "not converged"
1181  call prc_abort
1182  end if
1183 
1184  if ( .not. present(kref) ) then
1185  !$acc exit data delete(kref_)
1186  deallocate( kref_ )
1187  end if
1188 
1189  return
1190  end subroutine atmos_hydrostatic_buildrho_atmos_3d
1191 
1192  !-----------------------------------------------------------------------------
1194  subroutine atmos_hydrostatic_buildrho_atmos_rev_3d( &
1195  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1196  pott, qv, qc, &
1197  dz, &
1198  dens, &
1199  temp, pres, &
1200  kref )
1201  use scale_prc, only: &
1202  prc_abort
1203  implicit none
1204 
1205  integer, intent(in) :: KA, KS, KE
1206  integer, intent(in) :: IA, IS, IE
1207  integer, intent(in) :: JA, JS, JE
1208  real(RP), intent(inout) :: dens(KA,IA,JA)
1209  real(RP), intent(inout) :: temp(KA,IA,JA)
1210  real(RP), intent(inout) :: pres(KA,IA,JA)
1211  real(RP), intent(in) :: pott(KA,IA,JA)
1212  real(RP), intent(in) :: qv (KA,IA,JA)
1213  real(RP), intent(in) :: qc (KA,IA,JA)
1214  real(RP), intent(in) :: dz (KA,IA,JA)
1215  integer, intent(in), optional, target :: kref(IA,JA)
1216 
1217  integer, pointer :: kref_(:,:)
1218  logical :: converged, flag
1219  integer :: i, j
1220 #ifdef _OPENACC
1221  real(RP) :: work1(KA,IA,JA)
1222  real(RP) :: work2(KA,IA,JA)
1223 #endif
1224  !---------------------------------------------------------------------------
1225 
1226  converged = .true.
1227 
1228  if ( present(kref) ) then
1229  kref_ => kref
1230 
1231  !$omp parallel do OMP_SCHEDULE_ collapse(2) private(flag) reduction(.and.:converged)
1232  !$acc kernels copyin(pott, qv, qc, dz, kref_) copy(dens) copyout(temp, pres)
1233  !$acc loop independent collapse(2) &
1234  !$acc private(flag) reduction(.and.: converged)
1235  do j = js, je
1236  do i = is, ie
1237  if ( kref_(i,j) > ks ) then
1238  call atmos_hydrostatic_buildrho_atmos_rev_1d( ka, ks, kref_(i,j), &
1239  pott(:,i,j), qv(:,i,j), qc(:,i,j), & ! [IN]
1240  dz(:,i,j), & ! [IN]
1241 #ifdef _OPENACC
1242  work1(:,i,j), work2(:,i,j), & ! [WORK]
1243 #endif
1244  dens(:,i,j), & ! [INOUT]
1245  temp(:,i,j), pres(:,i,j), & ! [OUT]
1246  flag ) ! [OUT]
1247  converged = converged .and. flag
1248  end if
1249  enddo
1250  enddo
1251  !$acc end kernels
1252 
1253  else
1254  !$omp parallel do OMP_SCHEDULE_ collapse(2) private(flag) reduction(.and.:converged)
1255  !$acc kernels copyin(pott, qv, qc, dz) copy(dens) copyout(temp, pres)
1256  !$acc loop independent collapse(2) &
1257  !$acc private(flag) reduction(.and.: converged)
1258  do j = js, je
1259  do i = is, ie
1260  call atmos_hydrostatic_buildrho_atmos_rev_1d( ka, ks, ke, &
1261  pott(:,i,j), qv(:,i,j), qc(:,i,j), & ! [IN]
1262  dz(:,i,j), & ! [IN]
1263 #ifdef _OPENACC
1264  work1(:,i,j), work2(:,i,j), & ! [WORK]
1265 #endif
1266  dens(:,i,j), & ! [INOUT]
1267  temp(:,i,j), pres(:,i,j), & ! [OUT]
1268  flag ) ! [OUT]
1269  converged = converged .and. flag
1270  enddo
1271  enddo
1272  !$acc end kernels
1273  end if
1274 
1275  if ( .not. converged ) then
1276  log_error("ATMOS_HYDROSTATIC_buildrho_atmos_rev_3D",*) "not converged"
1277  call prc_abort
1278  end if
1279 
1280  return
1281  end subroutine atmos_hydrostatic_buildrho_atmos_rev_3d
1282 
1283  !-----------------------------------------------------------------------------
1285 !OCL SERIAL
1286  subroutine atmos_hydrostatic_buildrho_bytemp_1d( &
1287  KA, KS, KE, &
1288  temp, qv, qc, &
1289  pres_sfc, temp_sfc, qv_sfc, qc_sfc, &
1290  cz, fz, &
1291 #ifdef _OPENACC
1292  dz, work1, work2, work3, &
1293 #endif
1294  dens, pott, pres, pott_sfc, &
1295  converged )
1296  !$acc routine seq
1297  use scale_prc, only: &
1298  prc_abort
1299  use scale_atmos_hydrometeor, only: &
1300  cv_vapor, &
1301  cv_water, &
1302  cp_vapor, &
1303  cp_water
1304  implicit none
1305  integer, intent(in) :: KA, KS, KE
1306 
1307  real(RP), intent(in) :: temp(KA)
1308  real(RP), intent(in) :: qv (KA)
1309  real(RP), intent(in) :: qc (KA)
1310  real(RP), intent(in) :: pres_sfc
1311  real(RP), intent(in) :: temp_sfc
1312  real(RP), intent(in) :: qv_sfc
1313  real(RP), intent(in) :: qc_sfc
1314  real(RP), intent(in) :: cz( KA)
1315  real(RP), intent(in) :: fz(0:KA)
1316 
1317  real(RP), intent(out) :: dens(KA)
1318  real(RP), intent(out) :: pott(KA)
1319  real(RP), intent(out) :: pres(KA)
1320  real(RP), intent(out) :: pott_sfc
1321  logical, intent(out) :: converged
1322 
1323 #ifdef _OPENACC
1324  real(RP), intent(out) :: dz(KA)
1325  real(RP), intent(out) :: work1(KA)
1326  real(RP), intent(out) :: work2(KA)
1327  real(RP), intent(out) :: work3(KA)
1328 #else
1329  real(RP) :: dz(KA)
1330 #endif
1331 
1332  real(RP) :: dens_sfc
1333 
1334  real(RP) :: Rtot_sfc
1335  real(RP) :: CVtot_sfc
1336  real(RP) :: CPtot_sfc
1337  real(RP) :: Rtot
1338  real(RP) :: CVtot
1339  real(RP) :: CPtot
1340 
1341  real(RP) :: RovCP_sfc
1342  real(RP) :: dens_s, dhyd, dgrd
1343  integer :: ite
1344 
1345  integer :: k
1346  !---------------------------------------------------------------------------
1347 
1348  !--- from surface to lowermost atmosphere
1349 
1350  rtot_sfc = rdry * ( 1.0_rp - qv_sfc - qc_sfc ) &
1351  + rvap * qv_sfc
1352  cvtot_sfc = cvdry * ( 1.0_rp - qv_sfc - qc_sfc ) &
1353  + cv_vapor * qv_sfc &
1354  + cv_water * qc_sfc
1355  cptot_sfc = cpdry * ( 1.0_rp - qv_sfc - qc_sfc ) &
1356  + cp_vapor * qv_sfc &
1357  + cp_water * qc_sfc
1358 
1359  rtot = rdry * ( 1.0_rp - qv(ks) - qc(ks) ) &
1360  + rvap * qv(ks)
1361  cvtot = cvdry * ( 1.0_rp - qv(ks) - qc(ks) ) &
1362  + cv_vapor * qv(ks) &
1363  + cv_water * qc(ks)
1364  cptot = cpdry * ( 1.0_rp - qv(ks) - qc(ks) ) &
1365  + cp_vapor * qv(ks) &
1366  + cp_water * qc(ks)
1367 
1368  ! density at surface
1369  rovcp_sfc = rtot_sfc / cptot_sfc
1370  dens_sfc = pres_sfc / ( rtot_sfc * temp_sfc )
1371  pott_sfc = temp_sfc * ( p00/pres_sfc )**rovcp_sfc
1372 
1373  ! make density at lowermost cell center
1374  dens_s = 0.0_rp
1375  dens(ks) = dens_sfc ! first guess
1376 
1377  dz(ks-1) = cz(ks) - fz(ks-1)
1378  do k = ks, ke-1
1379  dz(k) = cz(k+1) - cz(k)
1380  end do
1381 
1382  converged = .false.
1383  do ite = 1, itelim
1384  if ( abs(dens(ks)-dens_s) <= criteria ) then
1385  converged = .true.
1386  exit
1387  endif
1388 
1389  dens_s = dens(ks)
1390 
1391  dhyd = + ( pres_sfc - dens_s * rtot * temp(ks) ) / dz(ks-1) & ! dp/dz
1392  - grav * 0.5_rp * ( dens_sfc + dens_s ) ! rho*g
1393 
1394  dgrd = - rtot * temp(ks) / dz(ks-1) &
1395  - grav * 0.5_rp
1396 
1397  dens(ks) = dens_s - dhyd/dgrd
1398 
1399  if ( dens(ks)*0.0_rp /= 0.0_rp ) exit
1400  enddo
1401 
1402  if ( converged ) then
1403  !--- from lowermost atmosphere to top of atmosphere
1404  call atmos_hydrostatic_buildrho_bytemp_atmos_1d( ka, ks, ke, &
1405  temp(:), qv(:), qc(:), & ! [IN]
1406  dz(:), & ! [IN]
1407 #ifdef _OPENACC
1408  work1(:), work2(:), work3(:), & ! [WORK]
1409 #endif
1410  dens(:), & ! [INOUT]
1411  pott(:), pres(:), & ! [OUT]
1412  converged ) ! [OUT]
1413  end if
1414 
1415 #ifndef _OPENACC
1416  if ( .not. converged ) then
1417  log_error("ATMOS_HYDROSTATIC_buildrho_bytemp_1D",*) 'iteration not converged!', &
1418  dens(ks),ite,dens_s,dhyd,dgrd
1419  call prc_abort
1420  endif
1421 #endif
1422 
1423  return
1424  end subroutine atmos_hydrostatic_buildrho_bytemp_1d
1425 
1426  !-----------------------------------------------------------------------------
1428  subroutine atmos_hydrostatic_buildrho_bytemp_3d( &
1429  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1430  temp, qv, qc, &
1431  pres_sfc, temp_sfc, qv_sfc, qc_sfc, &
1432  cz, fz, &
1433  dens, pott, pres, &
1434  pott_sfc )
1435  use scale_prc, only: &
1436  prc_abort
1437  implicit none
1438  integer, intent(in) :: KA, KS, KE
1439  integer, intent(in) :: IA, IS, IE
1440  integer, intent(in) :: JA, JS, JE
1441 
1442  real(RP), intent(in) :: temp(KA,IA,JA)
1443  real(RP), intent(in) :: qv (KA,IA,JA)
1444  real(RP), intent(in) :: qc (KA,IA,JA)
1445  real(RP), intent(in) :: pres_sfc(IA,JA)
1446  real(RP), intent(in) :: temp_sfc(IA,JA)
1447  real(RP), intent(in) :: qv_sfc (IA,JA)
1448  real(RP), intent(in) :: qc_sfc (IA,JA)
1449  real(RP), intent(in) :: cz(KA,IA,JA)
1450  real(RP), intent(in) :: fz(KA,IA,JA)
1451 
1452  real(RP), intent(out) :: dens(KA,IA,JA)
1453  real(RP), intent(out) :: pott(KA,IA,JA)
1454  real(RP), intent(out) :: pres(KA,IA,JA)
1455  real(RP), intent(out) :: pott_sfc(IA,JA)
1456 
1457  logical :: converged, flag
1458 
1459 #ifdef _OPENACC
1460  real(RP) :: work1(KA,IA,JA)
1461  real(RP) :: work2(KA,IA,JA)
1462  real(RP) :: work3(KA,IA,JA)
1463  real(RP) :: work4(KA,IA,JA)
1464 #endif
1465 
1466  integer :: i, j
1467  !---------------------------------------------------------------------------
1468 
1469  converged = .true.
1470 
1471  !--- from surface to lowermost atmosphere
1472  !$omp parallel do OMP_SCHEDULE_ collapse(2) private(flag) reduction(.and.:converged)
1473  !$acc kernels copyin(temp, qv, qc, pres_sfc, temp_sfc, qv_sfc, qc_sfc, cz, fz) copyout(dens, pott, pres, pott_sfc)
1474  !$acc loop independent collapse(2) &
1475  !$acc private(flag) reduction(.and.: converged)
1476  do j = js, je
1477  do i = is, ie
1478  call atmos_hydrostatic_buildrho_bytemp_1d( ka, ks, ke, &
1479  temp(:,i,j), qv(:,i,j), qc(:,i,j), & ! [IN]
1480  pres_sfc(i,j), temp_sfc(i,j), qv_sfc(i,j), qc_sfc(i,j), & ! [IN]
1481  cz(:,i,j), fz(:,i,j), & ! [IN]
1482 #ifdef _OPENACC
1483  work1(:,i,j), work2(:,i,j), work3(:,i,j), work4(:,i,j), & ! [WORK]
1484 #endif
1485  dens(:,i,j), pott(:,i,j), pres(:,i,j), & ! [OUT]
1486  pott_sfc(i,j), & ! [OUT]
1487  flag ) ! [OUT]
1488  converged = converged .and. flag
1489  enddo
1490  enddo
1491  !$acc end kernels
1492 
1493  if ( .not. converged ) then
1494  log_error("ATMOS_HYDROSTATIC_buildrho_bytemp_3D",*) "not converged"
1495  call prc_abort
1496  end if
1497 
1498  return
1499  end subroutine atmos_hydrostatic_buildrho_bytemp_3d
1500 
1501  !-----------------------------------------------------------------------------
1503 !OCL SERIAL
1504  subroutine atmos_hydrostatic_buildrho_bytemp_atmos_1d( &
1505  KA, KS, KE, &
1506  temp, qv, qc, &
1507  dz, &
1508 #ifdef _OPENACC
1509  Rtot, &
1510  CVtot, &
1511  CPtot, &
1512 #endif
1513  dens, &
1514  pott, pres, &
1515  converged )
1516  use scale_prc, only: &
1517  prc_abort
1518  use scale_atmos_hydrometeor, only: &
1519  cv_vapor, &
1520  cv_water, &
1521  cp_vapor, &
1522  cp_water
1523  !$acc routine seq
1524  implicit none
1525  integer, intent(in) :: KA, KS, KE
1526 
1527  real(RP), intent(in) :: temp(KA)
1528  real(RP), intent(in) :: qv (KA)
1529  real(RP), intent(in) :: qc (KA)
1530  real(RP), intent(in) :: dz (KA)
1531 
1532  real(RP), intent(inout) :: dens(KA)
1533 
1534  real(RP), intent(out) :: pott(KA)
1535  real(RP), intent(out) :: pres(KA)
1536  logical, intent(out) :: converged
1537 
1538 #ifdef _OPENACC
1539  real(RP), intent(out) :: Rtot (KA)
1540  real(RP), intent(out) :: CVtot (KA)
1541  real(RP), intent(out) :: CPtot (KA)
1542 #else
1543  real(RP) :: Rtot (KA)
1544  real(RP) :: CVtot (KA)
1545  real(RP) :: CPtot (KA)
1546 #endif
1547 
1548  real(RP) :: RovCP
1549  real(RP) :: dens_s, dhyd, dgrd
1550  integer :: ite
1551 
1552  integer :: k
1553  !---------------------------------------------------------------------------
1554 
1555  do k = ks, ke
1556  rtot(k) = rdry * ( 1.0_rp - qv(k) - qc(k) ) &
1557  + rvap * qv(k)
1558  cvtot(k) = cvdry * ( 1.0_rp - qv(k) - qc(k) ) &
1559  + cv_vapor * qv(k) &
1560  + cv_water * qc(k)
1561  cptot(k) = cpdry * ( 1.0_rp - qv(k) - qc(k) ) &
1562  + cp_vapor * qv(k) &
1563  + cp_water * qc(k)
1564  enddo
1565 
1566  pres(ks) = dens(ks) * rtot(ks) * temp(ks)
1567 
1568  do k = ks+1, ke
1569 
1570  dens_s = 0.0_rp
1571  dens(k) = dens(k-1) ! first guess
1572 
1573  converged = .false.
1574  do ite = 1, itelim
1575  if ( abs(dens(k)-dens_s) <= criteria ) then
1576  converged = .true.
1577  exit
1578  endif
1579 
1580  dens_s = dens(k)
1581 
1582  dhyd = + ( pres(k-1) - dens_s * rtot(k) * temp(k) ) / dz(k-1) & ! dp/dz
1583  - grav * 0.5_rp * ( dens(k-1) + dens_s ) ! rho*g
1584 
1585  dgrd = - rtot(k) * temp(k) / dz(k-1) &
1586  - 0.5_rp * grav
1587 
1588  dens(k) = dens_s - dhyd/dgrd
1589 
1590  if ( dens(k)*0.0_rp /= 0.0_rp ) exit
1591  enddo
1592 
1593  if ( .NOT. converged ) then
1594 #ifndef _OPENACC
1595  log_error("ATMOS_HYDROSTATIC_buildrho_bytemp_atmos_1D",*) 'iteration not converged!', &
1596  k,dens(k),ite,dens_s,dhyd,dgrd
1597  call prc_abort
1598 #endif
1599  exit
1600  endif
1601 
1602  pres(k) = dens(k) * rtot(k) * temp(k)
1603 
1604  enddo
1605 
1606  if ( converged ) then
1607  do k = ks, ke
1608  rovcp = rtot(k) / cptot(k)
1609  pott(k) = temp(k) * ( p00 / pres(k) )**rovcp
1610  enddo
1611  end if
1612 
1613  return
1614  end subroutine atmos_hydrostatic_buildrho_bytemp_atmos_1d
1615 
1616  !-----------------------------------------------------------------------------
1618 !OCL SERIAL
1619  subroutine atmos_hydrostatic_buildrho_bytemp_atmos_rev_1d( &
1620  KA, KS, KE, &
1621  temp, qv, qc, &
1622  dz, &
1623 #ifdef _OPENACC
1624  Rtot, &
1625  CVtot, &
1626  CPtot, &
1627 #endif
1628  dens, &
1629  pott, pres, &
1630  converged )
1631  use scale_prc, only: &
1632  prc_abort
1633  use scale_atmos_hydrometeor, only: &
1634  cv_vapor, &
1635  cv_water, &
1636  cp_vapor, &
1637  cp_water
1638  !$acc routine seq
1639  implicit none
1640  integer, intent(in) :: KA, KS, KE
1641 
1642  real(RP), intent(in) :: temp(KA)
1643  real(RP), intent(in) :: qv (KA)
1644  real(RP), intent(in) :: qc (KA)
1645  real(RP), intent(in) :: dz (KA)
1646 
1647  real(RP), intent(inout) :: dens(KA)
1648 
1649  real(RP), intent(out) :: pott(KA)
1650  real(RP), intent(out) :: pres(KA)
1651  logical, intent(out) :: converged
1652 
1653 #ifdef _OPENACC
1654  real(RP), intent(out) :: Rtot (KA)
1655  real(RP), intent(out) :: CVtot (KA)
1656  real(RP), intent(out) :: CPtot (KA)
1657 #else
1658  real(RP) :: Rtot (KA)
1659  real(RP) :: CVtot (KA)
1660  real(RP) :: CPtot (KA)
1661 #endif
1662 
1663  real(RP) :: RovCP
1664  real(RP) :: dens_s, dhyd, dgrd
1665  integer :: ite
1666 
1667  integer :: k
1668  !---------------------------------------------------------------------------
1669 
1670  do k = ks, ke
1671  rtot(k) = rdry * ( 1.0_rp - qv(k) - qc(k) ) &
1672  + rvap * qv(k)
1673  cvtot(k) = cvdry * ( 1.0_rp - qv(k) - qc(k) ) &
1674  + cv_vapor * qv(k) &
1675  + cv_water * qc(k)
1676  cptot(k) = cpdry * ( 1.0_rp - qv(k) - qc(k) ) &
1677  + cp_vapor * qv(k) &
1678  + cp_water * qc(k)
1679  enddo
1680 
1681  pres(ke) = dens(ke) * rtot(ke) * temp(ke)
1682 
1683  do k = ke-1, ke, -1
1684 
1685  dens_s = 0.0_rp
1686  dens(k) = dens(k+1) ! first guess
1687 
1688  converged = .false.
1689  do ite = 1, itelim
1690  if ( abs(dens(k)-dens_s) <= criteria ) then
1691  converged = .true.
1692  exit
1693  endif
1694 
1695  dens_s = dens(k)
1696 
1697  dhyd = - ( pres(k+1) - dens_s * rtot(k) * temp(k) ) / dz(k) & ! dp/dz
1698  - grav * 0.5_rp * ( dens(k+1) + dens_s ) ! rho*g
1699 
1700  dgrd = + rtot(k) * temp(k) / dz(k) &
1701  - 0.5_rp * grav
1702 
1703  dens(k) = dens_s - dhyd/dgrd
1704 
1705  if ( dens(k)*0.0_rp /= 0.0_rp ) exit
1706  enddo
1707 
1708  if ( .NOT. converged ) then
1709 #ifdef _OPENACC
1710  log_error("ATMOS_HYDROSTATIC_buildrho_bytemp_atmos_rev_1D",*) 'iteration not converged!', &
1711  k,dens(k),ite,dens_s,dhyd,dgrd
1712  call prc_abort
1713 #endif
1714  exit
1715  endif
1716 
1717  pres(k) = dens(k) * rtot(k) * temp(k)
1718 
1719  enddo
1720 
1721  if ( converged ) then
1722  do k = ks, ke
1723  rovcp = rtot(k) / cptot(k)
1724  pott(k) = temp(k) * ( p00 / pres(k) )**rovcp
1725  enddo
1726  end if
1727 
1728  return
1729  end subroutine atmos_hydrostatic_buildrho_bytemp_atmos_rev_1d
1730 
1731  !-----------------------------------------------------------------------------
1733  subroutine atmos_hydrostatic_buildrho_bytemp_atmos_3d( &
1734  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1735  temp, qv, qc, &
1736  dz, &
1737  dens, &
1738  pott, pres )
1739  use scale_prc, only: &
1740  prc_abort
1741  implicit none
1742  integer, intent(in) :: KA, KS, KE
1743  integer, intent(in) :: IA, IS, IE
1744  integer, intent(in) :: JA, JS, JE
1745 
1746  real(RP), intent(in) :: temp(KA,IA,JA)
1747  real(RP), intent(in) :: qv (KA,IA,JA)
1748  real(RP), intent(in) :: qc (KA,IA,JA)
1749  real(RP), intent(in) :: dz (KA,IA,JA)
1750 
1751  real(RP), intent(inout) :: dens(KA,IA,JA)
1752  real(RP), intent(out) :: pott(KA,IA,JA)
1753  real(RP), intent(out) :: pres(KA,IA,JA)
1754 
1755  logical :: converged, flag
1756 
1757 #ifdef _OPENACC
1758  real(RP) :: work1(KA,IA,JA)
1759  real(RP) :: work2(KA,IA,JA)
1760  real(RP) :: work3(KA,IA,JA)
1761 #endif
1762  integer :: i, j
1763  !---------------------------------------------------------------------------
1764 
1765  converged = .true.
1766 
1767  !$omp parallel do OMP_SCHEDULE_ collapse(2) private(flag) reduction(.and.:converged)
1768  !$acc kernels copyin(temp, qv, qc, dz) copy(dens) copyout(pott, pres)
1769  !$acc loop independent collapse(2) &
1770  !$acc private(flag) reduction(.and.: converged)
1771  do j = js, je
1772  do i = is, ie
1773  call atmos_hydrostatic_buildrho_bytemp_atmos_1d( ka, ks, ke, &
1774  temp(:,i,j), qv(:,i,j), qc(:,i,j), & ! [IN]
1775  dz(:,i,j), & ! [IN]
1776 #ifdef _OPENACC
1777  work1(:,i,j), work2(:,i,j), work3(:,i,j), & ! [WORK]
1778 #endif
1779  dens(:,i,j), & ! [INOUT]
1780  pott(:,i,j), pres(:,i,j), & ! [OUT]
1781  flag ) ! [OUT]
1782  converged = converged .and. flag
1783  enddo
1784  enddo
1785  !$acc end kernels
1786 
1787  if ( .not. converged ) then
1788  log_error("ATMOS_HYDROSTATIC_buildrho_bytemp_atmos_3D",*) "iteration not converged"
1789  call prc_abort
1790  end if
1791 
1792  return
1793  end subroutine atmos_hydrostatic_buildrho_bytemp_atmos_3d
1794 
1795  !-----------------------------------------------------------------------------
1797 !OCL SERIAL
1798 !OCL NOSIMD
1799  subroutine atmos_hydrostatic_barometric_law_mslp_0d( &
1800  KA, KS, KE, &
1801  pres, temp, qv, &
1802  cz, &
1803  mslp )
1804  !$acc routine seq
1805  implicit none
1806  integer, intent(in) :: KA, KS, KE
1807 
1808  real(RP), intent(in) :: pres(KA)
1809  real(RP), intent(in) :: temp(KA)
1810  real(RP), intent(in) :: qv (KA)
1811  real(RP), intent(in) :: cz (KA)
1812 
1813  real(RP), intent(out) :: mslp
1814 
1815  ! work
1816  integer :: kref
1817  real(RP) :: vtemp
1818  !---------------------------------------------------------------------------
1819 
1820  kref = hydrostatic_barometric_law_mslp_kref + ks - 1
1821 
1822  ! virtual temperature
1823  vtemp = temp(kref) * (1.0_rp + epstvap * qv(kref) )
1824 
1825  ! barometric law assuming constant lapse rate
1826  mslp = pres(kref) * ( ( vtemp + laps * cz(kref) ) / vtemp ) ** ( grav / ( rdry * laps ) )
1827 
1828  return
1829  end subroutine atmos_hydrostatic_barometric_law_mslp_0d
1830 
1831  !-----------------------------------------------------------------------------
1833  subroutine atmos_hydrostatic_barometric_law_mslp_2d( &
1834  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1835  pres, temp, qv, &
1836  cz, &
1837  mslp )
1838  implicit none
1839  integer, intent(in) :: KA, KS, KE
1840  integer, intent(in) :: IA, IS, IE
1841  integer, intent(in) :: JA, JS, JE
1842 
1843  real(RP), intent(in) :: pres(KA,IA,JA)
1844  real(RP), intent(in) :: temp(KA,IA,JA)
1845  real(RP), intent(in) :: qv (KA,IA,JA)
1846  real(RP), intent(in) :: cz (KA,IA,JA)
1847 
1848  real(RP), intent(out) :: mslp(IA,JA)
1849 
1850  integer :: i, j
1851  !---------------------------------------------------------------------------
1852 
1853  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1854  !$acc kernels copyin(pres, temp, qv, cz) copyout(mslp)
1855  do j = js, je
1856  do i = is, ie
1857  call atmos_hydrostatic_barometric_law_mslp_0d( ka, ks, ke, &
1858  pres(:,i,j), temp(:,i,j), qv(:,i,j), & ! [IN]
1859  cz(:,i,j), & ! [IN]
1860  mslp(i,j) ) ! [OUT]
1861  enddo
1862  enddo
1863  !$acc end kernels
1864 
1865  return
1866  end subroutine atmos_hydrostatic_barometric_law_mslp_2d
1867 
1868  !-----------------------------------------------------------------------------
1870 !OCL SERIAL
1871 !OCL NOSIMD
1872  subroutine atmos_hydrostatic_barometric_law_pres_0d( &
1873  mslp, temp, &
1874  dz, &
1875  pres )
1876  !$acc routine seq
1877  implicit none
1878  real(RP), intent(in) :: mslp
1879  real(RP), intent(in) :: temp
1880  real(RP), intent(in) :: dz
1881 
1882  real(RP), intent(out) :: pres
1883 
1884  ! work
1885  real(RP) :: TM
1886  !---------------------------------------------------------------------------
1887 
1888  tm = temp + laps * dz * 0.5_rp ! column-mean air temperature
1889 
1890  ! barometric law
1891  pres = mslp / exp( grav * dz / ( rdry * tm ) )
1892 
1893  return
1894  end subroutine atmos_hydrostatic_barometric_law_pres_0d
1895 
1896  !-----------------------------------------------------------------------------
1898  subroutine atmos_hydrostatic_barometric_law_pres_2d( &
1899  IA, IS, IE, JA, JS, JE, &
1900  mslp, temp, &
1901  dz, &
1902  pres )
1903  implicit none
1904  integer, intent(in) :: IA, IS, IE
1905  integer, intent(in) :: JA, JS, JE
1906 
1907  real(RP), intent(in) :: mslp(IA,JA)
1908  real(RP), intent(in) :: temp(IA,JA)
1909  real(RP), intent(in) :: dz (IA,JA)
1910 
1911  real(RP), intent(out) :: pres(IA,JA)
1912 
1913  integer :: i, j
1914  !---------------------------------------------------------------------------
1915 
1916  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1917  !$acc kernels copyin(mslp, temp, dz) copyout(pres)
1918  do j = js, je
1919  do i = is, ie
1920  call atmos_hydrostatic_barometric_law_pres_0d( mslp(i,j), temp(i,j), dz(i,j), & ! [IN]
1921  pres(i,j) ) ! [OUT]
1922  enddo
1923  enddo
1924  !$acc end kernels
1925 
1926  return
1927  end subroutine atmos_hydrostatic_barometric_law_pres_2d
1928 
1929 end module scale_atmos_hydrostatic
scale_const::const_grav
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:49
scale_statistics
module Statistics
Definition: scale_statistics.F90:11
scale_const::const_lapsdry
real(rp), public const_lapsdry
dry adiabatic lapse rate [K/m]
Definition: scale_const.F90:63
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_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:350
scale_const::const_epstvap
real(rp), public const_epstvap
1 / epsilon - 1
Definition: scale_const.F90:76
scale_atmos_hydrometeor::cp_water
real(rp), public cp_water
CP for water [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:152
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_bytemp_1d
subroutine atmos_hydrostatic_buildrho_bytemp_1d(KA, KS, KE, temp, qv, qc, pres_sfc, temp_sfc, qv_sfc, qc_sfc, cz, fz, dens, pott, pres, pott_sfc, converged)
Build up density from surface (1D)
Definition: scale_atmos_hydrostatic.F90:1296
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_rvap
real(rp), parameter, public const_rvap
specific gas constant (water vapor) [J/kg/K]
Definition: scale_const.F90:68
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_real_3d
subroutine atmos_hydrostatic_buildrho_real_3d(KA, KS, KE, IA, IS, IE, JA, JS, JE, pott, qv, qc, cz, pres, dens, temp)
Build up density from surface (3D), not to reverse from TOA.
Definition: scale_atmos_hydrostatic.F90:511
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:35
scale_atmos_hydrometeor
module atmosphere / hydrometeor
Definition: scale_atmos_hydrometeor.F90:12
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_atmos_rev_1d
subroutine atmos_hydrostatic_buildrho_atmos_rev_1d(KA, KS, KE, pott, qv, qc, dz, dens, temp, pres, converged)
Build up density from upermost atmosphere (1D)
Definition: scale_atmos_hydrostatic.F90:879
scale_const::const_cpvap
real(rp), parameter, public const_cpvap
specific heat (water vapor, constant pressure) [J/kg/K]
Definition: scale_const.F90:69
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_atmos_hydrometeor::cv_vapor
real(rp), public cv_vapor
CV for vapor [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:149
scale_io
module STDIO
Definition: scale_io.F90:10
scale_const
module CONSTANT
Definition: scale_const.F90:11
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_bytemp_atmos_3d
subroutine atmos_hydrostatic_buildrho_bytemp_atmos_3d(KA, KS, KE, IA, IS, IE, JA, JS, JE, temp, qv, qc, dz, dens, pott, pres)
Build up density from lowermost atmosphere (3D)
Definition: scale_atmos_hydrostatic.F90:1739
scale_const::const_cvdry
real(rp), public const_cvdry
specific heat (dry air,constant volume) [J/kg/K]
Definition: scale_const.F90:61
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_bytemp_3d
subroutine atmos_hydrostatic_buildrho_bytemp_3d(KA, KS, KE, IA, IS, IE, JA, JS, JE, temp, qv, qc, pres_sfc, temp_sfc, qv_sfc, qc_sfc, cz, fz, dens, pott, pres, pott_sfc)
Build up density from surface (3D)
Definition: scale_atmos_hydrostatic.F90:1435
scale_const::const_cpdry
real(rp), public const_cpdry
specific heat (dry air,constant pressure) [J/kg/K]
Definition: scale_const.F90:60
scale_prof
module profiler
Definition: scale_prof.F90:11
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_atmos_rev_3d
subroutine atmos_hydrostatic_buildrho_atmos_rev_3d(KA, KS, KE, IA, IS, IE, JA, JS, JE, pott, qv, qc, dz, dens, temp, pres, kref)
Build up density from lowermost atmosphere (3D)
Definition: scale_atmos_hydrostatic.F90:1201
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_atmos_1d
subroutine atmos_hydrostatic_buildrho_atmos_1d(KA, KS, KE, pott, qv, qc, dz, dens, temp, pres, converged)
Build up density from lowermost atmosphere (1D)
Definition: scale_atmos_hydrostatic.F90:759
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_atmos_3d
subroutine atmos_hydrostatic_buildrho_atmos_3d(KA, KS, KE, IA, IS, IE, JA, JS, JE, pott, qv, qc, dz, dens, temp, pres, kref)
Build up density from lowermost atmosphere (3D)
Definition: scale_atmos_hydrostatic.F90:1110
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_atmos_hydrostatic
module atmosphere / hydrostatic barance
Definition: scale_atmos_hydrostatic.F90:12
scale_const::const_cvvap
real(rp), public const_cvvap
specific heat (water vapor, constant volume) [J/kg/K]
Definition: scale_const.F90:70
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_3d
subroutine atmos_hydrostatic_buildrho_3d(KA, KS, KE, IA, IS, IE, JA, JS, JE, pott, qv, qc, pres_sfc, pott_sfc, qv_sfc, qc_sfc, cz, fz, area, dens, temp, pres, temp_sfc)
Build up density from surface (3D)
Definition: scale_atmos_hydrostatic.F90:346
scale_const::const_rdry
real(rp), public const_rdry
specific gas constant (dry air) [J/kg/K]
Definition: scale_const.F90:59
scale_atmos_hydrostatic::atmos_hydrostatic_barometric_law_mslp_0d
subroutine atmos_hydrostatic_barometric_law_mslp_0d(KA, KS, KE, pres, temp, qv, cz, mslp)
Calculate mean sea-level pressure from barometric law (0D)
Definition: scale_atmos_hydrostatic.F90:1804
scale_atmos_hydrostatic::atmos_hydrostatic_setup
subroutine, public atmos_hydrostatic_setup
Setup.
Definition: scale_atmos_hydrostatic.F90:120
scale_const::const_laps
real(rp), public const_laps
lapse rate of ISA [K/m]
Definition: scale_const.F90:62
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_1d
subroutine atmos_hydrostatic_buildrho_1d(KA, KS, KE, pott, qv, qc, pres_sfc, pott_sfc, qv_sfc, qc_sfc, cz, fz, dens, temp, pres, temp_sfc, converged)
Build up density from surface (1D)
Definition: scale_atmos_hydrostatic.F90:218
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_bytemp_atmos_1d
subroutine atmos_hydrostatic_buildrho_bytemp_atmos_1d(KA, KS, KE, temp, qv, qc, dz, dens, pott, pres, converged)
Build up density from lowermost atmosphere (1D)
Definition: scale_atmos_hydrostatic.F90:1516
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_atmos_0d
subroutine atmos_hydrostatic_buildrho_atmos_0d(pott_L2, qv_L2, qc_L2, dens_L1, pott_L1, qv_L1, qc_L1, dz, k, dens_L2, temp_L2, pres_L2, converged)
Build up density (0D)
Definition: scale_atmos_hydrostatic.F90:652
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:43
scale_const::const_pre00
real(rp), public const_pre00
pressure reference [Pa]
Definition: scale_const.F90:97
scale_io::io_fid_conf
integer, public io_fid_conf
Config file ID.
Definition: scale_io.F90:57
scale_atmos_hydrometeor::cp_vapor
real(rp), public cp_vapor
CP for vapor [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:150
scale_atmos_hydrometeor::cv_water
real(rp), public cv_water
CV for water [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:151
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_atmos_rev_2d
subroutine atmos_hydrostatic_buildrho_atmos_rev_2d(IA, IS, IE, JA, JS, JE, pott_L1, qv_L1, qc_L1, dens_L2, pott_L2, qv_L2, qc_L2, dz, k, dens_L1, temp_L1, pres_L1)
Build up density (2D)
Definition: scale_atmos_hydrostatic.F90:1050
scale_atmos_hydrostatic::atmos_hydrostatic_buildrho_bytemp_atmos_rev_1d
subroutine atmos_hydrostatic_buildrho_bytemp_atmos_rev_1d(KA, KS, KE, temp, qv, qc, dz, dens, pott, pres, converged)
Build up density from upermost atmosphere (1D)
Definition: scale_atmos_hydrostatic.F90:1631
Generated by   doxygen 1.8.17