SCALE-RM
scale_atmos_saturation.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
20 
21  use scale_const, only: &
22  rvap => const_rvap, &
23  psat0 => const_psat0, &
24  epsvap => const_epsvap, &
25  tem00 => const_tem00
26  !-----------------------------------------------------------------------------
27  implicit none
28  private
29  !-----------------------------------------------------------------------------
30  !
31  !++ Public procedure
32  !
33  public :: atmos_saturation_setup
34 
35  public :: atmos_saturation_alpha
36 
37  public :: atmos_saturation_psat_all
38  public :: atmos_saturation_psat_liq
39  public :: atmos_saturation_psat_ice
40 
41  public :: atmos_saturation_psat2qsat_pres
42  public :: atmos_saturation_psat2qsat_dens
43 
44  public :: atmos_saturation_pres2qsat_all
45  public :: atmos_saturation_pres2qsat_liq
46  public :: atmos_saturation_pres2qsat_ice
47 
48  public :: atmos_saturation_dens2qsat_all
49  public :: atmos_saturation_dens2qsat_liq
50  public :: atmos_saturation_dens2qsat_ice
51 
52  public :: atmos_saturation_dalphadt
53 
54  public :: atmos_saturation_dqs_dtem_dens_liq
55  public :: atmos_saturation_dqs_dtem_dens_ice
56  public :: atmos_saturation_dqs_dtem_dpre_liq
57  public :: atmos_saturation_dqs_dtem_dpre_ice
58 
59  public :: atmos_saturation_tdew_liq
60 
61  public :: atmos_saturation_pote
62 
63  public :: atmos_saturation_moist_conversion_dens_liq
64  public :: atmos_saturation_moist_conversion_dens_all
65  public :: atmos_saturation_moist_conversion_pres_liq
66 
67  interface atmos_saturation_alpha
68  module procedure atmos_saturation_alpha_0d
69  module procedure atmos_saturation_alpha_1d
70  module procedure atmos_saturation_alpha_3d
71  end interface atmos_saturation_alpha
72 
73  interface atmos_saturation_psat_all
74  module procedure atmos_saturation_psat_all_0d
75  module procedure atmos_saturation_psat_all_1d
76  module procedure atmos_saturation_psat_all_2d
77  module procedure atmos_saturation_psat_all_3d
78  end interface atmos_saturation_psat_all
79  interface atmos_saturation_psat_liq
80  module procedure atmos_saturation_psat_liq_0d
81  module procedure atmos_saturation_psat_liq_1d
82  module procedure atmos_saturation_psat_liq_2d
83  module procedure atmos_saturation_psat_liq_3d
84  end interface atmos_saturation_psat_liq
85  interface atmos_saturation_psat_ice
86  module procedure atmos_saturation_psat_ice_0d
87  module procedure atmos_saturation_psat_ice_1d
88  module procedure atmos_saturation_psat_ice_2d
89  module procedure atmos_saturation_psat_ice_3d
90  end interface atmos_saturation_psat_ice
91 
92  interface atmos_saturation_psat2qsat_pres
93  module procedure atmos_saturation_psat2qsat_pres_0d
94  end interface atmos_saturation_psat2qsat_pres
95  interface atmos_saturation_psat2qsat_dens
96  module procedure atmos_saturation_psat2qsat_dens_0d
97  end interface atmos_saturation_psat2qsat_dens
98 
99  interface atmos_saturation_pres2qsat_all
100  module procedure atmos_saturation_pres2qsat_all_0d
101  module procedure atmos_saturation_pres2qsat_all_1d
102  module procedure atmos_saturation_pres2qsat_all_2d
103  module procedure atmos_saturation_pres2qsat_all_3d
104  end interface atmos_saturation_pres2qsat_all
105  interface atmos_saturation_pres2qsat_liq
106  module procedure atmos_saturation_pres2qsat_liq_0d
107  module procedure atmos_saturation_pres2qsat_liq_1d
108  module procedure atmos_saturation_pres2qsat_liq_3d
109  end interface atmos_saturation_pres2qsat_liq
110  interface atmos_saturation_pres2qsat_ice
111  module procedure atmos_saturation_pres2qsat_ice_0d
112  module procedure atmos_saturation_pres2qsat_ice_1d
113  module procedure atmos_saturation_pres2qsat_ice_3d
114  end interface atmos_saturation_pres2qsat_ice
115 
116  interface atmos_saturation_dens2qsat_all
117  module procedure atmos_saturation_dens2qsat_all_0d
118  module procedure atmos_saturation_dens2qsat_all_1d
119  module procedure atmos_saturation_dens2qsat_all_3d
120  end interface atmos_saturation_dens2qsat_all
121  interface atmos_saturation_dens2qsat_liq
122  module procedure atmos_saturation_dens2qsat_liq_0d
123  module procedure atmos_saturation_dens2qsat_liq_1d
124  module procedure atmos_saturation_dens2qsat_liq_3d
125  end interface atmos_saturation_dens2qsat_liq
126  interface atmos_saturation_dens2qsat_ice
127  module procedure atmos_saturation_dens2qsat_ice_0d
128  module procedure atmos_saturation_dens2qsat_ice_1d
129  module procedure atmos_saturation_dens2qsat_ice_3d
130  end interface atmos_saturation_dens2qsat_ice
131 
132  interface atmos_saturation_dalphadt
133  module procedure atmos_saturation_dalphadt_0d
134  module procedure atmos_saturation_dalphadt_1d
135  module procedure atmos_saturation_dalphadt_3d
136  end interface atmos_saturation_dalphadt
137 
138  interface atmos_saturation_dqs_dtem_dens_liq
139  module procedure atmos_saturation_dqs_dtem_dens_liq_0d
141  module procedure atmos_saturation_dqs_dtem_dens_liq_3d
142  end interface atmos_saturation_dqs_dtem_dens_liq
143  interface atmos_saturation_dqs_dtem_dens_ice
144  module procedure atmos_saturation_dqs_dtem_dens_ice_0d
146  module procedure atmos_saturation_dqs_dtem_dens_ice_3d
147  end interface atmos_saturation_dqs_dtem_dens_ice
148  interface atmos_saturation_dqs_dtem_dpre_liq
149  module procedure atmos_saturation_dqs_dtem_dpre_liq_0d
151  module procedure atmos_saturation_dqs_dtem_dpre_liq_3d
152  end interface atmos_saturation_dqs_dtem_dpre_liq
153  interface atmos_saturation_dqs_dtem_dpre_ice
154  module procedure atmos_saturation_dqs_dtem_dpre_ice_0d
156  module procedure atmos_saturation_dqs_dtem_dpre_ice_3d
157  end interface atmos_saturation_dqs_dtem_dpre_ice
158 
159  interface atmos_saturation_tdew_liq
160  module procedure atmos_saturation_tdew_liq_0d
161  module procedure atmos_saturation_tdew_liq_1d
162  module procedure atmos_saturation_tdew_liq_3d
163  end interface atmos_saturation_tdew_liq
164 
165  interface atmos_saturation_pote
166  module procedure atmos_saturation_pote_0d
167  module procedure atmos_saturation_pote_1d
168  module procedure atmos_saturation_pote_3d
169  end interface atmos_saturation_pote
170 
171  interface atmos_saturation_moist_conversion_dens_liq
172  module procedure atmos_saturation_moist_conversion_dens_liq_0d
173  end interface atmos_saturation_moist_conversion_dens_liq
174  interface atmos_saturation_moist_conversion_dens_all
176  end interface atmos_saturation_moist_conversion_dens_all
177  interface atmos_saturation_moist_conversion_pres_liq
179  end interface atmos_saturation_moist_conversion_pres_liq
180 
181  !-----------------------------------------------------------------------------
182  !
183  !++ Public parameters & variables
184  !
185  !-----------------------------------------------------------------------------
186  !
187  !++ Private procedure
188  !
189  !-----------------------------------------------------------------------------
190  !
191  !++ Private parameters & variables
192  !
193  real(RP), private, parameter :: TEM_MIN = 10.0_rp
194 
195  real(RP), private :: ATMOS_SATURATION_ULIMIT_TEMP = 273.15_rp
196  real(RP), private :: ATMOS_SATURATION_LLIMIT_TEMP = 233.15_rp
197 
198  real(RP), private :: RTEM00
199  real(RP), private :: dalphadT_const
200  real(RP), private :: psat_min_liq
201  real(RP), private :: psat_min_ice
202 
203  real(RP), private :: CPovR_liq
204  real(RP), private :: CPovR_ice
205  real(RP), private :: CVovR_liq
206  real(RP), private :: CVovR_ice
207  real(RP), private :: LovR_liq
208  real(RP), private :: LovR_ice
209 
210  !-----------------------------------------------------------------------------
211 contains
212  !-----------------------------------------------------------------------------
214  subroutine atmos_saturation_setup
215  use scale_prc, only: &
216  prc_abort
217  use scale_const, only: &
218  cpvap => const_cpvap, &
219  cvvap => const_cvvap, &
220  cl => const_cl, &
221  ci => const_ci, &
222  lhv00 => const_lhv00, &
223  lhs00 => const_lhs00, &
224  lhv0 => const_lhv0, &
225  lhs0 => const_lhs0, &
227  use scale_atmos_hydrometeor, only: &
229  implicit none
230 
231  namelist / param_atmos_saturation / &
232  atmos_saturation_ulimit_temp, &
233  atmos_saturation_llimit_temp
234 
235  integer :: ierr
236  !---------------------------------------------------------------------------
237 
238  log_newline
239  log_info("ATMOS_SATURATION_setup",*) 'Setup'
240 
241  !--- read namelist
242  rewind(io_fid_conf)
243  read(io_fid_conf,nml=param_atmos_saturation,iostat=ierr)
244  if( ierr < 0 ) then !--- missing
245  log_info("ATMOS_SATURATION_setup",*) 'Not found namelist. Default used.'
246  elseif( ierr > 0 ) then !--- fatal error
247  log_error("ATMOS_SATURATION_setup",*) 'Not appropriate names in namelist PARAM_ATMOS_SATURATION. Check!'
248  call prc_abort
249  endif
250  log_nml(param_atmos_saturation)
251 
252  rtem00 = 1.0_rp / tem00
253 
254  if ( const_thermodyn_type == 'EXACT' ) then
255 
256  cpovr_liq = ( cpvap - cl ) / rvap
257  cpovr_ice = ( cpvap - ci ) / rvap
258  cvovr_liq = ( cvvap - cl ) / rvap
259  cvovr_ice = ( cvvap - ci ) / rvap
260 
261  lovr_liq = lhv00 / rvap
262  lovr_ice = lhs00 / rvap
263 
264  elseif( const_thermodyn_type == 'SIMPLE' ) then
265 
266  cpovr_liq = 0.0_rp
267  cpovr_ice = 0.0_rp
268  cvovr_liq = 0.0_rp
269  cvovr_ice = 0.0_rp
270 
271  lovr_liq = lhv0 / rvap
272  lovr_ice = lhs0 / rvap
273 
274  endif
275 
276  dalphadt_const = 1.0_rp / ( atmos_saturation_ulimit_temp - atmos_saturation_llimit_temp )
277 
278  log_newline
279  log_info("ATMOS_SATURATION_setup",'(1x,A,F7.2,A,F7.2)') 'Temperature range for liquid/ice mixture : ', &
280  atmos_saturation_llimit_temp, ' - ', &
281  atmos_saturation_ulimit_temp
282 
284 
285  call atmos_saturation_psat_liq( tem_min, psat_min_liq ) ! [IN], [OUT]
286  call atmos_saturation_psat_ice( tem_min, psat_min_ice ) ! [IN], [OUT]
287 
288  return
289  end subroutine atmos_saturation_setup
290 
291  !-----------------------------------------------------------------------------
293  subroutine atmos_saturation_alpha_0d( &
294  temp, &
295  alpha )
296  implicit none
297 
298  real(RP), intent(in) :: temp
299  real(RP), intent(out) :: alpha
300  !---------------------------------------------------------------------------
301 
302  alpha = ( temp - atmos_saturation_llimit_temp ) &
303  / ( atmos_saturation_ulimit_temp - atmos_saturation_llimit_temp )
304 
305  alpha = min( max( alpha, 0.0_rp ), 1.0_rp )
306 
307  return
308  end subroutine atmos_saturation_alpha_0d
309 
310  !-----------------------------------------------------------------------------
312 !OCL SERIAL
313  subroutine atmos_saturation_alpha_1d( &
314  KA, KS, KE, &
315  temp, &
316  alpha )
317  implicit none
318  integer, intent(in) :: KA, KS, KE
319 
320  real(RP), intent(in) :: temp (KA)
321 
322  real(RP), intent(out) :: alpha(KA)
323 
324  integer :: k
325  !---------------------------------------------------------------------------
326 
327  do k = ks, ke
328  call atmos_saturation_alpha_0d( temp(k), & ! [IN]
329  alpha(k) ) ! [OUT]
330  enddo
331 
332  return
333  end subroutine atmos_saturation_alpha_1d
334 
335  !-----------------------------------------------------------------------------
337  subroutine atmos_saturation_alpha_3d( &
338  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
339  temp, &
340  alpha )
341  implicit none
342  integer, intent(in) :: KA, KS, KE
343  integer, intent(in) :: IA, IS, IE
344  integer, intent(in) :: JA, JS, JE
345 
346  real(RP), intent(in) :: temp (KA,IA,JA)
347  real(RP), intent(out) :: alpha(KA,IA,JA)
348 
349  integer :: k, i, j
350  !---------------------------------------------------------------------------
351 
352  !$omp parallel do OMP_SCHEDULE_ collapse(2)
353  do j = js, je
354  do i = is, ie
355  do k = ks, ke
356  call atmos_saturation_alpha_0d( temp(k,i,j), & ! [IN]
357  alpha(k,i,j) ) ! [OUT]
358 
359  enddo
360  enddo
361  enddo
362 
363  return
364  end subroutine atmos_saturation_alpha_3d
365 
366  !-----------------------------------------------------------------------------
368  subroutine atmos_saturation_psat_all_0d( &
369  temp, &
370  psat )
371  implicit none
372 
373  real(RP), intent(in) :: temp
374  real(RP), intent(out) :: psat
375 
376  real(RP) :: alpha, psatl, psati
377  !---------------------------------------------------------------------------
378 
379  call atmos_saturation_alpha ( temp, alpha )
380  call atmos_saturation_psat_liq( temp, psatl )
381  call atmos_saturation_psat_ice( temp, psati )
382 
383  psat = psatl * ( alpha ) &
384  + psati * ( 1.0_rp - alpha )
385 
386  return
387  end subroutine atmos_saturation_psat_all_0d
388 
389  !-----------------------------------------------------------------------------
391 !OCL SERIAL
392  subroutine atmos_saturation_psat_all_1d( &
393  KA, KS, KE, &
394  temp, &
395  psat )
396  implicit none
397  integer, intent(in) :: KA, KS, KE
398 
399  real(RP), intent(in) :: temp(KA)
400  real(RP), intent(out) :: psat(KA)
401 
402  integer :: k
403  !---------------------------------------------------------------------------
404 
405  do k = ks, ke
406  call atmos_saturation_psat_all_0d( temp(k), & ! [IN]
407  psat(k) ) ! [OUT]
408  enddo
409 
410  return
411  end subroutine atmos_saturation_psat_all_1d
412 
413  !-----------------------------------------------------------------------------
415  subroutine atmos_saturation_psat_all_2d( &
416  IA, IS, IE, &
417  JA, JS, JE, &
418  temp, &
419  psat )
420  implicit none
421  integer, intent(in) :: IA, IS, IE
422  integer, intent(in) :: JA, JS, JE
423 
424  real(RP), intent(in) :: temp(IA,JA)
425 
426  real(RP), intent(out) :: psat(IA,JA)
427 
428  integer :: i, j
429  !---------------------------------------------------------------------------
430 
431  !$omp parallel do OMP_SCHEDULE_ collapse(2)
432  do j = js, je
433  do i = is, ie
434  call atmos_saturation_psat_all_0d( temp(i,j), & ! [IN]
435  psat(i,j) ) ! [OUT]
436  enddo
437  enddo
438 
439  return
440  end subroutine atmos_saturation_psat_all_2d
441 
442  !-----------------------------------------------------------------------------
444  subroutine atmos_saturation_psat_all_3d( &
445  KA, KS, KE, &
446  IA, IS, IE, &
447  JA, JS, JE, &
448  temp, &
449  psat )
450  implicit none
451 
452  integer, intent(in) :: KA, KS, KE
453  integer, intent(in) :: IA, IS, IE
454  integer, intent(in) :: JA, JS, JE
455 
456  real(RP), intent(in) :: temp(KA,IA,JA)
457  real(RP), intent(out) :: psat(KA,IA,JA)
458 
459  integer :: k, i, j
460  !---------------------------------------------------------------------------
461 
462  !$omp parallel do OMP_SCHEDULE_ collapse(2)
463  do j = js, je
464  do i = is, ie
465  do k = ks, ke
466  call atmos_saturation_psat_all_0d( temp(k,i,j), & ! [IN]
467  psat(k,i,j) ) ! [OUT]
468  enddo
469  enddo
470  enddo
471 
472  return
473  end subroutine atmos_saturation_psat_all_3d
474 
475  !-----------------------------------------------------------------------------
477  subroutine atmos_saturation_psat_liq_0d( &
478  temp, &
479  psat )
480  implicit none
481 
482  real(RP), intent(in) :: temp
483  real(RP), intent(out) :: psat
484  !---------------------------------------------------------------------------
485 
486  psat = psat0 * ( temp * rtem00 )**cpovr_liq &
487  * exp( lovr_liq * ( rtem00 - 1.0_rp/temp ) )
488 
489  return
490  end subroutine atmos_saturation_psat_liq_0d
491 
492  !-----------------------------------------------------------------------------
494 !OCL SERIAL
495  subroutine atmos_saturation_psat_liq_1d( &
496  KA, KS, KE, &
497  temp, &
498  psat )
499  implicit none
500  integer, intent(in) :: KA, KS, KE
501 
502  real(RP), intent(in) :: temp(KA)
503  real(RP), intent(out) :: psat(KA)
504 
505  integer :: k
506  !---------------------------------------------------------------------------
507 
508  do k = ks, ke
509  call atmos_saturation_psat_liq_0d( temp(k), &
510  psat(k) )
511  enddo
512 
513  return
514  end subroutine atmos_saturation_psat_liq_1d
515 
516  !-----------------------------------------------------------------------------
518  subroutine atmos_saturation_psat_liq_2d( &
519  IA, IS, IE, &
520  JA, JS, JE, &
521  temp, &
522  psat )
523  implicit none
524 
525  integer, intent(in) :: IA, IS, IE
526  integer, intent(in) :: JA, JS, JE
527 
528  real(RP), intent(in) :: temp(IA,JA)
529  real(RP), intent(out) :: psat(IA,JA)
530 
531  integer :: i, j
532  !---------------------------------------------------------------------------
533 
534  !$omp parallel do OMP_SCHEDULE_ collapse(2)
535  do j = js, je
536  do i = is, ie
537  call atmos_saturation_psat_liq_0d( temp(i,j), & ! [IN]
538  psat(i,j) ) ! [OUT]
539  enddo
540  enddo
541 
542  return
543  end subroutine atmos_saturation_psat_liq_2d
544 
545  !-----------------------------------------------------------------------------
547  subroutine atmos_saturation_psat_liq_3d( &
548  KA, KS, KE, &
549  IA, IS, IE, &
550  JA, JS, JE, &
551  temp, &
552  psat )
553  implicit none
554 
555  integer, intent(in) :: KA, KS, KE
556  integer, intent(in) :: IA, IS, IE
557  integer, intent(in) :: JA, JS, JE
558 
559  real(RP), intent(in) :: temp(KA,IA,JA)
560  real(RP), intent(out) :: psat(KA,IA,JA)
561 
562  integer :: k, i, j
563  !---------------------------------------------------------------------------
564 
565  !$omp parallel do OMP_SCHEDULE_ collapse(2)
566  do j = js, je
567  do i = is, ie
568  do k = ks, ke
569  call atmos_saturation_psat_liq_0d( temp(k,i,j), & ! [IN]
570  psat(k,i,j) ) ! [OUT]
571  enddo
572  enddo
573  enddo
574 
575  return
576  end subroutine atmos_saturation_psat_liq_3d
577 
578  !-----------------------------------------------------------------------------
580  subroutine atmos_saturation_psat_ice_0d( &
581  temp, &
582  psat )
583  implicit none
584 
585  real(RP), intent(in) :: temp
586 
587  real(RP), intent(out) :: psat
588  !---------------------------------------------------------------------------
589 
590  psat = psat0 * ( temp * rtem00 )**cpovr_ice &
591  * exp( lovr_ice * ( rtem00 - 1.0_rp/temp ) )
592 
593  return
594  end subroutine atmos_saturation_psat_ice_0d
595 
596  !-----------------------------------------------------------------------------
598 !OCL SERIAL
599  subroutine atmos_saturation_psat_ice_1d( &
600  KA, KS, KE, &
601  temp, &
602  psat )
603  implicit none
604  integer, intent(in) :: KA, KS, KE
605 
606  real(RP), intent(in) :: temp(KA)
607  real(RP), intent(out) :: psat(KA)
608 
609  integer :: k
610  !---------------------------------------------------------------------------
611 
612  do k = ks, ke
613  call atmos_saturation_psat_ice_0d( temp(k), & ! [IN]
614  psat(k) ) ! [OUT]
615  enddo
616 
617  return
618  end subroutine atmos_saturation_psat_ice_1d
619 
620  !-----------------------------------------------------------------------------
622  subroutine atmos_saturation_psat_ice_2d( &
623  IA, IS, IE, &
624  JA, JS, JE, &
625  temp, &
626  psat )
627  implicit none
628 
629  integer, intent(in) :: IA, IS, IE
630  integer, intent(in) :: JA, JS, JE
631 
632  real(RP), intent(in) :: temp(IA,JA)
633  real(RP), intent(out) :: psat(IA,JA)
634 
635  integer :: i, j
636  !---------------------------------------------------------------------------
637 
638  !$omp parallel do OMP_SCHEDULE_ collapse(2)
639  do j = js, je
640  do i = is, ie
641  call atmos_saturation_psat_ice_0d( temp(i,j), & ! [IN]
642  psat(i,j) ) ! [OUT]
643  enddo
644  enddo
645 
646  return
647  end subroutine atmos_saturation_psat_ice_2d
648 
649  !-----------------------------------------------------------------------------
651  subroutine atmos_saturation_psat_ice_3d( &
652  KA, KS, KE, &
653  IA, IS, IE, &
654  JA, JS, JE, &
655  temp, &
656  psat )
657  implicit none
658 
659  integer, intent(in) :: KA, KS, KE
660  integer, intent(in) :: IA, IS, IE
661  integer, intent(in) :: JA, JS, JE
662 
663  real(RP), intent(in) :: temp(KA,IA,JA)
664  real(RP), intent(out) :: psat(KA,IA,JA)
665 
666  integer :: k, i, j
667  !---------------------------------------------------------------------------
668 
669  !$omp parallel do OMP_SCHEDULE_ collapse(2)
670  do j = js, je
671  do i = is, ie
672  do k = ks, ke
673  call atmos_saturation_psat_ice_0d( temp(k,i,j), & ! [IN]
674  psat(k,i,j) ) ! [OUT]
675  enddo
676  enddo
677  enddo
678 
679  return
680  end subroutine atmos_saturation_psat_ice_3d
681 
682 
683  !-----------------------------------------------------------------------------
685  subroutine atmos_saturation_psat2qsat_pres_0d( &
686  psat, pres, qdry, &
687  qsat )
688  implicit none
689 
690  real(RP), intent(in) :: psat
691  real(RP), intent(in) :: pres
692  real(RP), intent(in) :: qdry
693 
694  real(RP), intent(out) :: qsat
695 
696  ! ! qdry is assumed to be 1 - qsat
697  ! qsat = EPSvap * psat / ( pres - ( 1.0_RP-EPSvap ) * psat )
698 
699  qsat = epsvap * qdry * psat / ( pres - psat )
700 
701  return
702  end subroutine atmos_saturation_psat2qsat_pres_0d
703  !-----------------------------------------------------------------------------
705  subroutine atmos_saturation_pres2qsat_all_0d( &
706  temp, pres, qdry, &
707  qsat )
708  implicit none
709 
710  real(RP), intent(in) :: temp
711  real(RP), intent(in) :: pres
712  real(RP), intent(in) :: qdry
713 
714  real(RP), intent(out) :: qsat
715 
716  real(RP) :: psat
717  !---------------------------------------------------------------------------
718 
719  call atmos_saturation_psat_all_0d( temp, psat ) ! [IN], [OUT]
720  call atmos_saturation_psat2qsat_pres_0d( psat, pres, qdry, & ! [IN]
721  qsat ) ! [OUT]
722 
723  return
724  end subroutine atmos_saturation_pres2qsat_all_0d
725 
726  !-----------------------------------------------------------------------------
728 !OCL SERIAL
729  subroutine atmos_saturation_pres2qsat_all_1d( &
730  KA, KS, KE, &
731  temp, pres, qdry, &
732  qsat )
733  implicit none
734  integer, intent(in) :: KA, KS, KE
735 
736  real(RP), intent(in) :: temp(KA)
737  real(RP), intent(in) :: pres(KA)
738  real(RP), intent(in) :: qdry(KA)
739 
740  real(RP), intent(out) :: qsat(KA)
741 
742  integer :: k
743  !---------------------------------------------------------------------------
744 
745  do k = ks, ke
746  call atmos_saturation_pres2qsat_all_0d( temp(k), pres(k), qdry(k), & ! [IN]
747  qsat(k) ) ! [OUT]
748  enddo
749 
750  return
751  end subroutine atmos_saturation_pres2qsat_all_1d
752 
753  !-----------------------------------------------------------------------------
755  subroutine atmos_saturation_pres2qsat_all_2d( &
756  IA, IS, IE, JA, JS, JE, &
757  temp, pres, qdry, &
758  qsat )
759  implicit none
760  integer, intent(in) :: IA, IS, IE
761  integer, intent(in) :: JA, JS, JE
762 
763  real(RP), intent(in) :: temp(IA,JA)
764  real(RP), intent(in) :: pres(IA,JA)
765  real(RP), intent(in) :: qdry(IA,JA)
766 
767  real(RP), intent(out) :: qsat(IA,JA)
768 
769  integer :: i, j
770  !---------------------------------------------------------------------------
771 
772  !$omp parallel do OMP_SCHEDULE_ collapse(2)
773  do j = js, je
774  do i = is, ie
775  call atmos_saturation_pres2qsat_all_0d( temp(i,j), pres(i,j), qdry(i,j), &
776  qsat(i,j) )
777  enddo
778  enddo
779 
780  return
781  end subroutine atmos_saturation_pres2qsat_all_2d
782 
783  !-----------------------------------------------------------------------------
785  subroutine atmos_saturation_pres2qsat_all_3d( &
786  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
787  temp, pres, qdry, &
788  qsat )
789  implicit none
790  integer, intent(in) :: KA, KS, KE
791  integer, intent(in) :: IA, IS, IE
792  integer, intent(in) :: JA, JS, JE
793 
794  real(RP), intent(in) :: temp(KA,IA,JA)
795  real(RP), intent(in) :: pres(KA,IA,JA)
796  real(RP), intent(in) :: qdry(KA,IA,JA)
797 
798  real(RP), intent(out) :: qsat(KA,IA,JA)
799 
800  integer :: k, i, j
801  !---------------------------------------------------------------------------
802 
803  !$omp parallel do OMP_SCHEDULE_ collapse(2)
804  do j = js, je
805  do i = is, ie
806  do k = ks, ke
807  call atmos_saturation_pres2qsat_all_0d( temp(k,i,j), pres(k,i,j), qdry(k,i,j), & ! [IN]
808  qsat(k,i,j) ) ! [OUT]
809  enddo
810  enddo
811  enddo
812 
813  return
814  end subroutine atmos_saturation_pres2qsat_all_3d
815 
816  !-----------------------------------------------------------------------------
818  subroutine atmos_saturation_pres2qsat_liq_0d( &
819  temp, pres, qdry, &
820  qsat )
821  implicit none
822 
823  real(RP), intent(in) :: temp
824  real(RP), intent(in) :: pres
825  real(RP), intent(in) :: qdry
826 
827  real(RP), intent(out) :: qsat
828 
829  real(RP) :: psat
830  !---------------------------------------------------------------------------
831 
832  call atmos_saturation_psat_liq( temp, psat ) ! [IN], [OUT]
833  call atmos_saturation_psat2qsat_pres( psat, pres, qdry, & ! [IN]
834  qsat ) ! [OUT]
835 
836  return
837  end subroutine atmos_saturation_pres2qsat_liq_0d
838 
839  !-----------------------------------------------------------------------------
841 !OCL SERIAL
842  subroutine atmos_saturation_pres2qsat_liq_1d( &
843  KA, KS, KE, &
844  temp, pres, qdry, &
845  qsat )
846  implicit none
847  integer, intent(in) :: KA, KS, KE
848 
849  real(RP), intent(in) :: temp(KA)
850  real(RP), intent(in) :: pres(KA)
851  real(RP), intent(in) :: qdry(KA)
852 
853  real(RP), intent(out) :: qsat(KA)
854 
855  integer :: k
856  !---------------------------------------------------------------------------
857 
858  do k = ks, ke
859  call atmos_saturation_pres2qsat_liq_0d( temp(k), pres(k), qdry(k), & ! [IN]
860  qsat(k) ) ! [OUT]
861  enddo
862 
863  return
864  end subroutine atmos_saturation_pres2qsat_liq_1d
865 
866  !-----------------------------------------------------------------------------
868  subroutine atmos_saturation_pres2qsat_liq_3d( &
869  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
870  temp, pres, qdry, &
871  qsat )
872  implicit none
873  integer, intent(in) :: KA, KS, KE
874  integer, intent(in) :: IA, IS, IE
875  integer, intent(in) :: JA, JS, JE
876 
877  real(RP), intent(in) :: temp(KA,IA,JA)
878  real(RP), intent(in) :: pres(KA,IA,JA)
879  real(RP), intent(in) :: qdry(KA,IA,JA)
880 
881  real(RP), intent(out) :: qsat(KA,IA,JA)
882 
883  integer :: k, i, j
884  !---------------------------------------------------------------------------
885 
886  !$omp parallel do OMP_SCHEDULE_ collapse(2)
887  do j = js, je
888  do i = is, ie
889  do k = ks, ke
890  call atmos_saturation_pres2qsat_liq_0d( temp(k,i,j), pres(k,i,j), qdry(k,i,j), & ! [IN]
891  qsat(k,i,j) ) ! [OUT]
892  enddo
893  enddo
894  enddo
895 
896  return
897  end subroutine atmos_saturation_pres2qsat_liq_3d
898 
899  !-----------------------------------------------------------------------------
901  subroutine atmos_saturation_pres2qsat_ice_0d( &
902  temp, pres, qdry, &
903  qsat )
904  implicit none
905 
906  real(RP), intent(in) :: temp
907  real(RP), intent(in) :: pres
908  real(RP), intent(in) :: qdry
909 
910  real(RP), intent(out) :: qsat
911 
912  real(RP) :: psat
913  !---------------------------------------------------------------------------
914 
915  call atmos_saturation_psat_ice( temp, psat ) ! [IN], [OUT]
916  call atmos_saturation_psat2qsat_pres( psat, pres, qdry, & ! [IN]
917  qsat ) ! [OUT]
918 
919  return
920  end subroutine atmos_saturation_pres2qsat_ice_0d
921 
922  !-----------------------------------------------------------------------------
924 !OCL SERIAL
925  subroutine atmos_saturation_pres2qsat_ice_1d( &
926  KA, KS, KE, &
927  temp, pres, qdry, &
928  qsat )
929  implicit none
930  integer, intent(in) :: KA, KS, KE
931 
932  real(RP), intent(in) :: temp(KA)
933  real(RP), intent(in) :: pres(KA)
934  real(RP), intent(in) :: qdry(KA)
935 
936  real(RP), intent(out) :: qsat(KA)
937 
938  integer :: k
939  !---------------------------------------------------------------------------
940 
941  do k = ks, ke
942  call atmos_saturation_pres2qsat_ice_0d( temp(k), pres(k), qdry(k), & ! [IN]
943  qsat(k) ) ! [OUT]
944  enddo
945 
946  return
947  end subroutine atmos_saturation_pres2qsat_ice_1d
948 
949  !-----------------------------------------------------------------------------
951  subroutine atmos_saturation_pres2qsat_ice_3d( &
952  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
953  temp, pres, qdry, &
954  qsat )
955  implicit none
956  integer, intent(in) :: KA, KS, KE
957  integer, intent(in) :: IA, IS, IE
958  integer, intent(in) :: JA, JS, JE
959 
960  real(RP), intent(in) :: temp(KA,IA,JA)
961  real(RP), intent(in) :: pres(KA,IA,JA)
962  real(RP), intent(in) :: qdry(KA,IA,JA)
963 
964  real(RP), intent(out) :: qsat(KA,IA,JA)
965 
966  integer :: k, i, j
967  !---------------------------------------------------------------------------
968 
969  !$omp parallel do OMP_SCHEDULE_ collapse(2)
970  do j = js, je
971  do i = is, ie
972  do k = ks, ke
973  call atmos_saturation_pres2qsat_ice_0d( temp(k,i,j), pres(k,i,j), qdry(k,i,j), & ! [IN]
974  qsat(k,i,j) ) ! [OUT]
975  enddo
976  enddo
977  enddo
978 
979  return
980  end subroutine atmos_saturation_pres2qsat_ice_3d
981 
982  !-----------------------------------------------------------------------------
984  subroutine atmos_saturation_psat2qsat_dens_0d( &
985  psat, temp, dens, &
986  qsat )
987  implicit none
988  real(RP), intent(in) :: psat
989  real(RP), intent(in) :: temp
990  real(RP), intent(in) :: dens
991 
992  real(RP), intent(out) :: qsat
993 
994  qsat = psat / ( dens * rvap * temp )
995 
996  return
997  end subroutine atmos_saturation_psat2qsat_dens_0d
998 
999  !-----------------------------------------------------------------------------
1001  subroutine atmos_saturation_dens2qsat_all_0d( &
1002  temp, &
1003  dens, &
1004  qsat )
1005  implicit none
1006  real(RP), intent(in) :: temp
1007  real(RP), intent(in) :: dens
1008 
1009  real(RP), intent(out) :: qsat
1010 
1011  real(RP) :: psat
1012  !---------------------------------------------------------------------------
1013 
1014  call atmos_saturation_psat_all( temp, psat ) ! [IN], [OUT]
1015  call atmos_saturation_psat2qsat_dens_0d( psat, temp, dens, & ! [IN]
1016  qsat ) ! [OUT]
1017 
1018  return
1019  end subroutine atmos_saturation_dens2qsat_all_0d
1020 
1021  !-----------------------------------------------------------------------------
1023 !OCL SERIAL
1024  subroutine atmos_saturation_dens2qsat_all_1d( &
1025  KA, KS, KE, &
1026  temp, dens, &
1027  qsat )
1028  implicit none
1029  integer, intent(in) :: KA, KS, KE
1030 
1031  real(RP), intent(in) :: temp(KA)
1032  real(RP), intent(in) :: dens(KA)
1033 
1034  real(RP), intent(out) :: qsat(KA)
1035 
1036  integer :: k
1037  !---------------------------------------------------------------------------
1038 
1039  do k = ks, ke
1040  call atmos_saturation_dens2qsat_all_0d( temp(k), dens(k), & ! [IN]
1041  qsat(k) ) ! [OUT]
1042  enddo
1043 
1044  return
1045  end subroutine atmos_saturation_dens2qsat_all_1d
1046 
1047  !-----------------------------------------------------------------------------
1049  subroutine atmos_saturation_dens2qsat_all_3d( &
1050  KA, KS, KE, &
1051  IA, IS, IE, &
1052  JA, JS, JE, &
1053  temp, &
1054  dens, &
1055  qsat )
1056  implicit none
1057  integer, intent(in) :: KA, KS, KE
1058  integer, intent(in) :: IA, IS, IE
1059  integer, intent(in) :: JA, JS, JE
1060 
1061  real(RP), intent(in) :: temp(KA,IA,JA)
1062  real(RP), intent(in) :: dens(KA,IA,JA)
1063 
1064  real(RP), intent(out) :: qsat(KA,IA,JA)
1065 
1066  integer :: k, i, j
1067  !---------------------------------------------------------------------------
1068 
1069  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1070  do j = js, je
1071  do i = is, ie
1072  do k = ks, ke
1073  call atmos_saturation_dens2qsat_all_0d( temp(k,i,j), dens(k,i,j), & ! [IN]
1074  qsat(k,i,j) ) ! [OUT]
1075  enddo
1076  enddo
1077  enddo
1078 
1079  return
1080  end subroutine atmos_saturation_dens2qsat_all_3d
1081 
1082  !-----------------------------------------------------------------------------
1084  subroutine atmos_saturation_dens2qsat_liq_0d( &
1085  temp, &
1086  dens, &
1087  qsat )
1088  implicit none
1089  real(RP), intent(in) :: temp
1090  real(RP), intent(in) :: dens
1091 
1092  real(RP), intent(out) :: qsat
1093 
1094  real(RP) :: psat
1095  !---------------------------------------------------------------------------
1096 
1097  call atmos_saturation_psat_liq_0d( temp, psat ) ! [IN], [OUT]
1098  call atmos_saturation_psat2qsat_dens_0d( psat, temp, dens, & ! [IN]
1099  qsat ) ! [OUT]
1100 
1101  return
1102  end subroutine atmos_saturation_dens2qsat_liq_0d
1103 
1104  !-----------------------------------------------------------------------------
1106 !OCL SERIAL
1107  subroutine atmos_saturation_dens2qsat_liq_1d( &
1108  KA, KS, KE, &
1109  temp, dens, &
1110  qsat )
1111  implicit none
1112  integer, intent(in) :: KA, KS, KE
1113 
1114  real(RP), intent(in) :: temp(KA)
1115  real(RP), intent(in) :: dens(KA)
1116 
1117  real(RP), intent(out) :: qsat(KA)
1118 
1119  integer :: k
1120  !---------------------------------------------------------------------------
1121 
1122  do k = ks, ke
1123  call atmos_saturation_dens2qsat_liq_0d( temp(k), dens(k), & ! [IN]
1124  qsat(k) ) ! [OUT]
1125  enddo
1126 
1127  return
1128  end subroutine atmos_saturation_dens2qsat_liq_1d
1129 
1130  !-----------------------------------------------------------------------------
1132  subroutine atmos_saturation_dens2qsat_liq_3d( &
1133  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1134  temp, dens, &
1135  qsat )
1136  implicit none
1137  integer, intent(in) :: KA, KS, KE
1138  integer, intent(in) :: IA, IS, IE
1139  integer, intent(in) :: JA, JS, JE
1140 
1141  real(RP), intent(in) :: temp(KA,IA,JA)
1142  real(RP), intent(in) :: dens(KA,IA,JA)
1143 
1144  real(RP), intent(out) :: qsat(KA,IA,JA)
1145 
1146  integer :: k, i, j
1147  !---------------------------------------------------------------------------
1148 
1149  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1150  do j = js, je
1151  do i = is, ie
1152  do k = ks, ke
1153  call atmos_saturation_dens2qsat_liq_0d( temp(k,i,j), dens(k,i,j), & ! [IN]
1154  qsat(k,i,j) ) ! [OUT]
1155  enddo
1156  enddo
1157  enddo
1158 
1159  return
1160  end subroutine atmos_saturation_dens2qsat_liq_3d
1161 
1162  !-----------------------------------------------------------------------------
1164  subroutine atmos_saturation_dens2qsat_ice_0d( &
1165  temp, &
1166  dens, &
1167  qsat )
1168  implicit none
1169  real(RP), intent(in) :: temp
1170  real(RP), intent(in) :: dens
1171 
1172  real(RP), intent(out) :: qsat
1173 
1174  real(RP) :: psat
1175  !---------------------------------------------------------------------------
1176 
1177  call atmos_saturation_psat_ice( temp, psat ) ! [IN], [OUT]
1178  call atmos_saturation_psat2qsat_dens( psat, temp, dens, & ! [IN]
1179  qsat ) ! [OUT]
1180  return
1181  end subroutine atmos_saturation_dens2qsat_ice_0d
1182 
1183  !-----------------------------------------------------------------------------
1185 !OCL SERIAL
1186  subroutine atmos_saturation_dens2qsat_ice_1d( &
1187  KA, KS, KE, &
1188  temp, dens, &
1189  qsat )
1190  implicit none
1191  integer, intent(in) :: KA, KS, KE
1192 
1193  real(RP), intent(in) :: temp(KA)
1194  real(RP), intent(in) :: dens(KA)
1195 
1196  real(RP), intent(out) :: qsat(KA)
1197 
1198  integer :: k
1199  !---------------------------------------------------------------------------
1200 
1201  do k = ks, ke
1202  call atmos_saturation_dens2qsat_ice_0d( temp(k), dens(k), & ! [IN]
1203  qsat(k) ) ! [OUT]
1204  enddo
1205 
1206  return
1207  end subroutine atmos_saturation_dens2qsat_ice_1d
1208 
1209  !-----------------------------------------------------------------------------
1211  subroutine atmos_saturation_dens2qsat_ice_3d( &
1212  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1213  temp, dens, &
1214  qsat )
1215  implicit none
1216  integer, intent(in) :: KA, KS, KE
1217  integer, intent(in) :: IA, IS, IE
1218  integer, intent(in) :: JA, JS, JE
1219 
1220  real(RP), intent(in) :: temp(KA,IA,JA)
1221  real(RP), intent(in) :: dens(KA,IA,JA)
1222 
1223  real(RP), intent(out) :: qsat(KA,IA,JA)
1224 
1225  integer :: k, i, j
1226  !---------------------------------------------------------------------------
1227 
1228  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1229  do j = js, je
1230  do i = is, ie
1231  do k = ks, ke
1232  call atmos_saturation_dens2qsat_ice_0d( temp(k,i,j), dens(k,i,j), & ! [IN]
1233  qsat(k,i,j) ) ! [OUT]
1234  enddo
1235  enddo
1236  enddo
1237 
1238  return
1239  end subroutine atmos_saturation_dens2qsat_ice_3d
1240 
1241  !-----------------------------------------------------------------------------
1243  subroutine atmos_saturation_dalphadt_0d( &
1244  temp, &
1245  dalpha_dT )
1246  implicit none
1247  real(RP), intent(in) :: temp
1248 
1249  real(RP), intent(out) :: dalpha_dT
1250 
1251  real(RP) :: lim1, lim2
1252  !---------------------------------------------------------------------------
1253 
1254  ! if Tup < temp, dalpha/dT = 0 (no slope)
1255  lim1 = 0.5_rp + sign( 0.5_rp, atmos_saturation_ulimit_temp - temp )
1256  ! if Tdn > temp, dalpha/dT = 0 (no slope)
1257  lim2 = 0.5_rp + sign( 0.5_rp, temp - atmos_saturation_llimit_temp )
1258 
1259  dalpha_dt = dalphadt_const * lim1 * lim2
1260 
1261  return
1262  end subroutine atmos_saturation_dalphadt_0d
1263 
1264  !-----------------------------------------------------------------------------
1266 !OCL SERIAL
1267  subroutine atmos_saturation_dalphadt_1d( &
1268  KA, KS, KE, &
1269  temp, &
1270  dalpha_dT )
1271  implicit none
1272  integer, intent(in) :: KA, KS, KE
1273 
1274  real(RP), intent(in) :: temp (KA)
1275 
1276  real(RP), intent(out) :: dalpha_dT(KA)
1277 
1278  integer :: k
1279  !---------------------------------------------------------------------------
1280 
1281  do k = ks, ke
1282  call atmos_saturation_dalphadt_0d( temp(k), dalpha_dt(k) ) ! [IN], [OUT]
1283  enddo
1284 
1285  return
1286  end subroutine atmos_saturation_dalphadt_1d
1287 
1288  !-----------------------------------------------------------------------------
1290  subroutine atmos_saturation_dalphadt_3d( &
1291  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1292  temp, &
1293  dalpha_dT )
1294  implicit none
1295  integer, intent(in) :: KA, KS, KE
1296  integer, intent(in) :: IA, IS, IE
1297  integer, intent(in) :: JA, JS, JE
1298 
1299  real(RP), intent(in) :: temp (KA,IA,JA)
1300 
1301  real(RP), intent(out) :: dalpha_dT(KA,IA,JA)
1302 
1303  integer :: k, i, j
1304  !---------------------------------------------------------------------------
1305 
1306  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1307  do j = js, je
1308  do i = is, ie
1309  do k = ks, ke
1310  call atmos_saturation_dalphadt_0d( temp(k,i,j), dalpha_dt(k,i,j) ) ! [IN], [OUT]
1311  enddo
1312  enddo
1313  enddo
1314 
1315  return
1316  end subroutine atmos_saturation_dalphadt_3d
1317 
1318  !-----------------------------------------------------------------------------
1319  ! (d qsw/d T)_{rho}: partial difference of qsat_water at constant density
1320  subroutine atmos_saturation_dqs_dtem_dens_liq_0d( &
1321  temp, dens, &
1322  dqsdtem, qsat )
1323  use scale_atmos_hydrometeor, only: &
1324  hydrometeor_lhv => atmos_hydrometeor_lhv
1325  implicit none
1326 
1327  real(RP), intent(in) :: temp
1328  real(RP), intent(in) :: dens
1329 
1330  real(RP), intent(out) :: dqsdtem
1331  real(RP), intent(out), optional :: qsat
1332 
1333  real(RP) :: LHV ! latent heat of vaporization [J/kg]
1334  real(RP) :: psat
1335  !---------------------------------------------------------------------------
1336 
1337  call hydrometeor_lhv( temp, lhv )
1338  call atmos_saturation_psat_liq( temp, psat ) ! [IN], [OUT]
1339 
1340  dqsdtem = psat / ( dens* rvap * temp**2 ) &
1341  * ( lhv / ( rvap * temp ) - 1.0_rp )
1342 
1343  if ( present(qsat) ) &
1344  call atmos_saturation_psat2qsat_dens( psat, temp, dens, & ! [IN]
1345  qsat ) ! [OUT]
1346 
1347  return
1348  end subroutine atmos_saturation_dqs_dtem_dens_liq_0d
1349 
1350  !-----------------------------------------------------------------------------
1351 !OCL SERIAL
1353  KA, KS, KE, &
1354  temp, dens, &
1355  dqsdtem )
1356  implicit none
1357  integer, intent(in) :: KA, KS, KE
1358 
1359  real(RP), intent(in) :: temp (KA)
1360  real(RP), intent(in) :: dens (KA)
1361 
1362  real(RP), intent(out) :: dqsdtem(KA)
1363 
1364  integer :: k
1365  !---------------------------------------------------------------------------
1366 
1367  do k = ks, ke
1368  call atmos_saturation_dqs_dtem_dens_liq_0d( temp(k), dens(k), & ! [IN]
1369  dqsdtem(k) ) ! [OUT]
1370  enddo
1371 
1372  return
1374 
1375  !-----------------------------------------------------------------------------
1376  ! (d qsw/d T)_{rho}: partial difference of qsat_water at constant density
1377  subroutine atmos_saturation_dqs_dtem_dens_liq_3d( &
1378  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1379  temp, dens, &
1380  dqsdtem )
1381  implicit none
1382  integer, intent(in) :: KA, KS, KE
1383  integer, intent(in) :: IA, IS, IE
1384  integer, intent(in) :: JA, JS, JE
1385 
1386  real(RP), intent(in) :: temp (KA,IA,JA)
1387  real(RP), intent(in) :: dens (KA,IA,JA)
1388 
1389  real(RP), intent(out) :: dqsdtem(KA,IA,JA)
1390 
1391  integer :: k, i, j
1392  !---------------------------------------------------------------------------
1393 
1394  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1395  do j = js, je
1396  do i = is, ie
1397  do k = ks, ke
1398  call atmos_saturation_dqs_dtem_dens_liq_0d( temp(k,i,j), dens(k,i,j), & ! [IN]
1399  dqsdtem(k,i,j) ) ! [OUT]
1400  enddo
1401  enddo
1402  enddo
1403 
1404  return
1405  end subroutine atmos_saturation_dqs_dtem_dens_liq_3d
1406 
1407  !-----------------------------------------------------------------------------
1408  ! (d qsi/d T)_{rho}: partial difference of qsat_ice at constant density
1409  subroutine atmos_saturation_dqs_dtem_dens_ice_0d( &
1410  temp, dens, &
1411  dqsdtem, qsat )
1412  use scale_atmos_hydrometeor, only: &
1413  hydrometeor_lhs => atmos_hydrometeor_lhs
1414  implicit none
1415 
1416  real(RP), intent(in) :: temp
1417  real(RP), intent(in) :: dens
1418 
1419  real(RP), intent(out) :: dqsdtem
1420  real(RP), intent(out), optional :: qsat
1421 
1422  real(RP) :: LHS ! latent heat of sublimation [J/kg]
1423  real(RP) :: psat
1424 
1425  !---------------------------------------------------------------------------
1426 
1427  call hydrometeor_lhs( temp, lhs ) ! [IN], [OUT]
1428  call atmos_saturation_psat_ice( temp, psat ) ! [IN], [OUT]
1429 
1430  dqsdtem = psat / ( dens * rvap * temp**2 ) &
1431  * ( lhs / ( rvap * temp ) - 1.0_rp )
1432 
1433  if ( present(qsat) ) &
1434  call atmos_saturation_psat2qsat_dens( psat, temp, dens, & ! [IN]
1435  qsat ) ! [OUT]
1436 
1437  return
1438  end subroutine atmos_saturation_dqs_dtem_dens_ice_0d
1439 
1440  !-----------------------------------------------------------------------------
1441 !OCL SERIAL
1443  KA, KS, KE, &
1444  temp, dens, &
1445  dqsdtem )
1446  implicit none
1447  integer, intent(in) :: KA, KS, KE
1448 
1449  real(RP), intent(in) :: temp (KA)
1450  real(RP), intent(in) :: dens (KA)
1451 
1452  real(RP), intent(out) :: dqsdtem(KA)
1453 
1454  integer :: k
1455  !---------------------------------------------------------------------------
1456 
1457  do k = ks, ke
1458  call atmos_saturation_dqs_dtem_dens_ice_0d( temp(k), dens(k), & ! [IN]
1459  dqsdtem(k) ) ! [OUT]
1460  enddo
1461 
1462  return
1464 
1465  !-----------------------------------------------------------------------------
1466  ! (d qsi/d T)_{rho}: partial difference of qsat_ice at constant density
1467  subroutine atmos_saturation_dqs_dtem_dens_ice_3d( &
1468  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1469  temp, dens, &
1470  dqsdtem )
1471  implicit none
1472  integer, intent(in) :: KA, KS, KE
1473  integer, intent(in) :: IA, IS, IE
1474  integer, intent(in) :: JA, JS, JE
1475 
1476  real(RP), intent(in) :: temp (KA,IA,JA)
1477  real(RP), intent(in) :: dens (KA,IA,JA)
1478 
1479  real(RP), intent(out) :: dqsdtem(KA,IA,JA)
1480 
1481  integer :: k, i, j
1482  !---------------------------------------------------------------------------
1483 
1484  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1485  do j = js, je
1486  do i = is, ie
1487  do k = ks, ke
1488  call atmos_saturation_dqs_dtem_dens_ice_0d( temp(k,i,j), dens(k,i,j), & ! [IN]
1489  dqsdtem(k,i,j) ) ! [OUT]
1490  enddo
1491  enddo
1492  enddo
1493 
1494  return
1495  end subroutine atmos_saturation_dqs_dtem_dens_ice_3d
1496 
1497  !-----------------------------------------------------------------------------
1498  ! (d qsi/d T)_{rho}: partial difference of qsat_all at constant density
1499  subroutine atmos_saturation_dqs_dtem_dens_all_0d( &
1500  temp, dens, &
1501  dqsat_dT, &
1502  qsat, qsat_liq, qsat_ice, &
1503  alpha )
1504  implicit none
1505 
1506  real(RP), intent(in) :: temp
1507  real(RP), intent(in) :: dens
1508 
1509  real(RP), intent(out) :: dqsat_dT
1510 
1511  real(RP), intent(out), optional :: qsat
1512  real(RP), intent(out), optional :: qsat_liq
1513  real(RP), intent(out), optional :: qsat_ice
1514  real(RP), intent(out), optional :: alpha
1515 
1516  real(RP) :: qsat_liq_, qsat_ice_, alpha_
1517  real(RP) :: dqsat_dT_liq, dqsat_dT_ice, dalpha_dT
1518 
1519  !---------------------------------------------------------------------------
1520 
1521  call atmos_saturation_dqs_dtem_dens_liq_0d( temp, dens, & ! [IN]
1522  dqsat_dt_liq, qsat_liq_ ) ! [OUT]
1523  call atmos_saturation_dqs_dtem_dens_ice_0d( temp, dens, & ! [IN]
1524  dqsat_dt_ice, qsat_ice_ ) ! [OUT]
1525  call atmos_saturation_alpha ( temp, alpha_ ) ! [IN], [OUT]
1526  call atmos_saturation_dalphadt( temp, dalpha_dt ) ! [IN], [OUT]
1527 
1528  dqsat_dt = qsat_liq_ * dalpha_dt + dqsat_dt_liq * ( alpha_ ) &
1529  - qsat_ice_ * dalpha_dt + dqsat_dt_ice * ( 1.0_rp-alpha_ )
1530 
1531  if ( present(qsat) ) qsat = qsat_liq_ * alpha_ + qsat_ice_ * ( 1.0_rp-alpha_ )
1532  if ( present(qsat_liq) ) qsat_liq = qsat_liq_
1533  if ( present(qsat_ice) ) qsat_ice = qsat_ice_
1534  if ( present(alpha ) ) alpha = alpha_
1535 
1536  return
1537  end subroutine atmos_saturation_dqs_dtem_dens_all_0d
1538 
1539  !-----------------------------------------------------------------------------
1540  ! (d qs/d T)_{p} and (d qs/d p)_{T}
1541  subroutine atmos_saturation_dqs_dtem_dpre_liq_0d( &
1542  temp, pres, qdry, &
1543  dqsat_dT, dqsat_dP, &
1544  qsat, psat )
1545  use scale_atmos_hydrometeor, only: &
1546  hydrometeor_lhv => atmos_hydrometeor_lhv
1547  implicit none
1548 
1549  real(RP), intent(in) :: temp
1550  real(RP), intent(in) :: pres
1551  real(RP), intent(in) :: qdry
1552 
1553  real(RP), intent(out) :: dqsat_dT
1554  real(RP), intent(out) :: dqsat_dP
1555 
1556  real(RP), intent(out), optional :: qsat
1557  real(RP), intent(out), optional :: psat
1558 
1559  real(RP) :: LHV ! latent heat of vaporization [J/kg]
1560  real(RP) :: psat_
1561  real(RP) :: den1, den2 ! denominator
1562 
1563  !---------------------------------------------------------------------------
1564 
1565  call hydrometeor_lhv( temp, lhv )
1566  call atmos_saturation_psat_liq( temp, psat_ ) ! [IN], [OUT]
1567 
1568  den1 = ( pres - (1.0_rp-epsvap) * psat_ )**2
1569  den2 = den1 * rvap * temp**2
1570 
1571  dqsat_dp = - epsvap * psat_ / den1
1572  dqsat_dt = epsvap * psat_ / den2 * lhv * pres
1573 
1574  if ( present(qsat) ) &
1575  call atmos_saturation_psat2qsat_pres( psat_, pres, qdry, & ! [IN]
1576  qsat ) ! [OUT]
1577  if ( present(psat) ) psat = psat_
1578 
1579  return
1580  end subroutine atmos_saturation_dqs_dtem_dpre_liq_0d
1581 
1582  !-----------------------------------------------------------------------------
1583 !OCL SERIAL
1585  KA, KS, KE, &
1586  temp, pres, qdry, &
1587  dqsat_dT, dqsat_dP )
1588  implicit none
1589  integer, intent(in) :: KA, KS, KE
1590 
1591  real(RP), intent(in) :: temp (KA)
1592  real(RP), intent(in) :: pres (KA)
1593  real(RP), intent(in) :: qdry (KA)
1594 
1595  real(RP), intent(out) :: dqsat_dT(KA)
1596  real(RP), intent(out) :: dqsat_dP(KA)
1597 
1598  integer :: k
1599  !---------------------------------------------------------------------------
1600 
1601  do k = ks, ke
1602  call atmos_saturation_dqs_dtem_dpre_liq_0d( temp(k), pres(k), qdry(k), & ! [IN]
1603  dqsat_dt(k), dqsat_dp(k) ) ! [OUT]
1604  enddo
1605 
1606  return
1608 
1609  !-----------------------------------------------------------------------------
1610  ! (d qs/d T)_{p} and (d qs/d p)_{T}
1611  subroutine atmos_saturation_dqs_dtem_dpre_liq_3d( &
1612  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1613  temp, pres, qdry, &
1614  dqsat_dT, dqsat_dP )
1615  implicit none
1616  integer, intent(in) :: KA, KS, KE
1617  integer, intent(in) :: IA, IS, IE
1618  integer, intent(in) :: JA, JS, JE
1619 
1620  real(RP), intent(in) :: temp (KA,IA,JA)
1621  real(RP), intent(in) :: pres (KA,IA,JA)
1622  real(RP), intent(in) :: qdry (KA,IA,JA)
1623 
1624  real(RP), intent(out) :: dqsat_dT(KA,IA,JA)
1625  real(RP), intent(out) :: dqsat_dP(KA,IA,JA)
1626 
1627  integer :: k, i, j
1628  !---------------------------------------------------------------------------
1629 
1630  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1631  do j = js, je
1632  do i = is, ie
1633  do k = ks, ke
1634  call atmos_saturation_dqs_dtem_dpre_liq_0d( temp(k,i,j), pres(k,i,j), qdry(k,i,j), & ! [IN]
1635  dqsat_dt(k,i,j), dqsat_dp(k,i,j) ) ! [OUT]
1636  enddo
1637  enddo
1638  enddo
1639 
1640  return
1641  end subroutine atmos_saturation_dqs_dtem_dpre_liq_3d
1642 
1643  !-----------------------------------------------------------------------------
1644  ! (d qsi/d T)_{p} and (d qs/d p)_{T}
1645  !-----------------------------------------------------------------------------
1646  subroutine atmos_saturation_dqs_dtem_dpre_ice_0d( &
1647  temp, pres, qdry, &
1648  dqsat_dT, dqsat_dP, &
1649  qsat )
1650  use scale_atmos_hydrometeor, only: &
1651  hydrometeor_lhs => atmos_hydrometeor_lhs
1652  implicit none
1653 
1654  real(RP), intent(in) :: temp
1655  real(RP), intent(in) :: pres
1656  real(RP), intent(in) :: qdry
1657 
1658  real(RP), intent(out) :: dqsat_dT
1659  real(RP), intent(out) :: dqsat_dP
1660 
1661  real(RP), intent(out), optional :: qsat
1662 
1663  real(RP) :: LHS ! latent heat of sublimation [J/kg]
1664  real(RP) :: psat
1665  real(RP) :: den1, den2 ! denominator
1666 
1667  !---------------------------------------------------------------------------
1668 
1669  call hydrometeor_lhs( temp, lhs ) ! [IN], [OUT]
1670  call atmos_saturation_psat_ice( temp, psat ) ! [IN], [OUT]
1671 
1672  den1 = ( pres - (1.0_rp-epsvap) * psat )**2
1673  den2 = den1 * rvap * temp**2
1674 
1675  dqsat_dp = - epsvap * psat / den1
1676  dqsat_dt = epsvap * psat / den2 * lhs * pres
1677 
1678  if ( present(qsat) ) &
1679  call atmos_saturation_psat2qsat_pres( psat, pres, qdry, & ! [IN]
1680  qsat ) ! [OUT]
1681  return
1682  end subroutine atmos_saturation_dqs_dtem_dpre_ice_0d
1683 
1684  !-----------------------------------------------------------------------------
1685 !OCL SERIAL
1687  KA, KS, KE, &
1688  temp, pres, qdry, &
1689  dqsat_dT, dqsat_dP )
1690  implicit none
1691  integer, intent(in) :: KA, KS, KE
1692 
1693  real(RP), intent(in) :: temp (KA)
1694  real(RP), intent(in) :: pres (KA)
1695  real(RP), intent(in) :: qdry (KA)
1696 
1697  real(RP), intent(out) :: dqsat_dT(KA)
1698  real(RP), intent(out) :: dqsat_dP(KA)
1699 
1700  integer :: k
1701  !---------------------------------------------------------------------------
1702 
1703  do k = ks, ke
1704  call atmos_saturation_dqs_dtem_dpre_ice( temp(k), pres(k), qdry(k), & ! [IN]
1705  dqsat_dt(k), dqsat_dp(k) ) ! [OUT]
1706  enddo
1707 
1708  return
1710 
1711  !-----------------------------------------------------------------------------
1712  ! (d qsi/d T)_{p} and (d qs/d p)_{T}
1713  !-----------------------------------------------------------------------------
1714  subroutine atmos_saturation_dqs_dtem_dpre_ice_3d( &
1715  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1716  temp, pres, qdry, &
1717  dqsat_dT, dqsat_dP )
1718  implicit none
1719  integer, intent(in) :: KA, KS, KE
1720  integer, intent(in) :: IA, IS, IE
1721  integer, intent(in) :: JA, JS, JE
1722 
1723  real(RP), intent(in) :: temp (KA,IA,JA)
1724  real(RP), intent(in) :: pres (KA,IA,JA)
1725  real(RP), intent(in) :: qdry (KA,IA,JA)
1726 
1727  real(RP), intent(out) :: dqsat_dT(KA,IA,JA)
1728  real(RP), intent(out) :: dqsat_dP(KA,IA,JA)
1729 
1730  integer :: k, i, j
1731  !---------------------------------------------------------------------------
1732 
1733  !$omp parallel do OMP_SCHEDULE_ collapse(2)
1734  do j = js, je
1735  do i = is, ie
1736  do k = ks, ke
1737  call atmos_saturation_dqs_dtem_dpre_ice( temp(k,i,j), pres(k,i,j), qdry(k,i,j), & ! [IN]
1738  dqsat_dt(k,i,j), dqsat_dp(k,i,j) ) ! [OUT]
1739  enddo
1740  enddo
1741  enddo
1742 
1743  return
1744  end subroutine atmos_saturation_dqs_dtem_dpre_ice_3d
1745 
1746  !-----------------------------------------------------------------------------
1747  ! (d qsi/d T)_{p} and (d qs/d p)_{T}
1748  !-----------------------------------------------------------------------------
1749  subroutine atmos_saturation_dqs_dtem_dpre_all_0d( &
1750  temp, pres, qdry, &
1751  dqsat_dT, dqsat_dP, &
1752  qsat, qsat_liq, qsat_ice, &
1753  alpha )
1754  implicit none
1755 
1756  real(RP), intent(in) :: temp
1757  real(RP), intent(in) :: pres
1758  real(RP), intent(in) :: qdry
1759 
1760  real(RP), intent(out) :: dqsat_dT
1761  real(RP), intent(out) :: dqsat_dP
1762 
1763  real(RP), intent(out), optional :: qsat
1764  real(RP), intent(out), optional :: qsat_liq
1765  real(RP), intent(out), optional :: qsat_ice
1766  real(RP), intent(out), optional :: alpha
1767 
1768  real(RP) :: qsat_liq_, qsat_ice_, alpha_
1769  real(RP) :: dqsat_dT_liq, dqsat_dT_ice
1770  real(RP) :: dqsat_dP_liq, dqsat_dP_ice
1771  real(RP) :: dalpha_dT
1772 
1773  !---------------------------------------------------------------------------
1774 
1775  call atmos_saturation_dqs_dtem_dpre_liq_0d( temp, pres, qdry, & ! [IN]
1776  dqsat_dt_liq, dqsat_dp_liq, & ! [OUT]
1777  qsat_liq_ ) ! [OUT]
1778  call atmos_saturation_dqs_dtem_dpre_ice_0d( temp, pres, qdry, & ! [IN]
1779  dqsat_dt_ice, dqsat_dp_ice, & ! [OUT]
1780  qsat_ice_ ) ! [OUT]
1781  call atmos_saturation_alpha ( temp, alpha_ ) ! [IN], [OUT]
1782  call atmos_saturation_dalphadt( temp, dalpha_dt ) ! [IN], [OUT]
1783 
1784  dqsat_dt = qsat_liq_ * dalpha_dt + dqsat_dt_liq * ( alpha_ ) &
1785  - qsat_ice_ * dalpha_dt + dqsat_dt_ice * ( 1.0_rp-alpha_ )
1786  dqsat_dp = dqsat_dp_liq * ( alpha_ ) &
1787  + dqsat_dp_ice * ( 1.0_rp-alpha_ )
1788 
1789  if ( present(qsat) ) qsat = qsat_liq * alpha_ + qsat_ice * ( 1.0_rp-alpha_ )
1790  if ( present(qsat_liq) ) qsat_liq = qsat_liq_
1791  if ( present(qsat_ice) ) qsat_ice = qsat_ice_
1792  if ( present(alpha ) ) alpha = alpha_
1793 
1794  return
1795  end subroutine atmos_saturation_dqs_dtem_dpre_all_0d
1796 
1797  !-----------------------------------------------------------------------------
1799  !-----------------------------------------------------------------------------
1800  subroutine atmos_saturation_tdew_liq_0d( &
1801  DENS, TEMP, QV, &
1802  Tdew, &
1803  converged )
1804  use scale_const, only: &
1805  undef => const_undef
1806  use scale_atmos_hydrometeor, only: &
1807  atmos_hydrometeor_lhv
1808  real(RP), intent(in) :: DENS
1809  real(RP), intent(in) :: TEMP
1810  real(RP), intent(in) :: QV
1811 
1812  real(RP), intent(out) :: Tdew
1813  logical, intent(out) :: converged
1814 
1815  real(RP), parameter :: A = 17.625_rp
1816  real(RP), parameter :: B = 243.04_rp
1817  real(RP), parameter :: C = 610.94_rp
1818  integer, parameter :: itelim = 100
1819  real(RP), parameter :: criteria = 0.1_rp**(2+rp/2)
1820 
1821  real(RP) :: lhv
1822  real(RP) :: pvap, psat
1823  real(RP) :: dpsat_dT
1824  real(RP) :: dTdew
1825 
1826  integer :: ite
1827  !---------------------------------------------------------------------------
1828 
1829  pvap = dens * qv * rvap * temp
1830 
1831  if ( pvap < psat_min_liq ) then
1832  converged = .true.
1833  tdew = undef
1834  return
1835  end if
1836 
1837  ! first guess is calculated by Alduchov and Eskridge (1996)
1838  ! See Lawrence (2005) BAMS
1839  tdew = b * log( pvap / c ) / ( a - log( pvap / c ) ) + tem00
1840  converged = .false.
1841  do ite = 1, itelim
1842 
1843  call atmos_saturation_psat_liq( tdew, psat ) ! [IN], [OUT]
1844  call atmos_hydrometeor_lhv( tdew, lhv )
1845 
1846  dpsat_dt = psat * lhv / ( rvap * tdew**2 )
1847  dtdew = ( psat - pvap ) / dpsat_dt
1848  if ( dtdew < criteria ) then
1849  converged = .true.
1850  exit
1851  end if
1852 
1853  tdew = tdew - dtdew
1854  end do
1855  if( .not. converged ) then
1856  log_warn("ATMOS_SATURATION_tdew_liq_0D",*) dens, temp, qv, pvap, tdew, dtdew, dpsat_dt
1857  endif
1858 
1859  return
1860  end subroutine atmos_saturation_tdew_liq_0d
1861 
1862 !OCL SERIAL
1863  subroutine atmos_saturation_tdew_liq_1d( &
1864  KA, KS, KE, &
1865  DENS, TEMP, QV, &
1866  Tdew, &
1867  converged )
1868  integer, intent(in) :: KA, KS, KE
1869 
1870  real(RP), intent(in) :: DENS(KA)
1871  real(RP), intent(in) :: TEMP(KA)
1872  real(RP), intent(in) :: QV (KA)
1873 
1874  real(RP), intent(out) :: Tdew(KA)
1875  logical, intent(out) :: converged
1876 
1877  integer :: k
1878  !---------------------------------------------------------------------------
1879 
1880  do k = ks, ke
1881  call atmos_saturation_tdew_liq_0d( dens(k), temp(k), qv(k), & ! [IN]
1882  tdew(k), converged ) ! [OUT]
1883  if ( .not. converged ) exit
1884  end do
1885 
1886  end subroutine atmos_saturation_tdew_liq_1d
1887 
1888  subroutine atmos_saturation_tdew_liq_3d( &
1889  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1890  DENS, TEMP, QV, &
1891  Tdew )
1892  use scale_prc, only: &
1893  prc_abort
1894  integer, intent(in) :: KA, KS, KE
1895  integer, intent(in) :: IA, IS, IE
1896  integer, intent(in) :: JA, JS, JE
1897 
1898  real(RP), intent(in) :: DENS(KA,IA,JA)
1899  real(RP), intent(in) :: TEMP(KA,IA,JA)
1900  real(RP), intent(in) :: QV (KA,IA,JA)
1901 
1902  real(RP), intent(out) :: Tdew(KA,IA,JA)
1903 
1904  logical :: converged, error
1905  integer :: k, i, j
1906  !---------------------------------------------------------------------------
1907 
1908  error = .false.
1909  !$omp parallel do OMP_SCHEDULE_ collapse(2) &
1910  !$omp private(converged)
1911  do j = js, je
1912  do i = is, ie
1913  do k = ks, ke
1914  call atmos_saturation_tdew_liq_0d( dens(k,i,j), temp(k,i,j), qv(k,i,j), & ! [IN]
1915  tdew(k,i,j), converged ) ! [OUT]
1916  if ( .not. converged ) then
1917  log_error("ATMOS_SATURATION_tdew_liq_3D",*) 'not converged! ', k,i,j
1918  error = .true.
1919  exit
1920  end if
1921  end do
1922  end do
1923  end do
1924 
1925  if ( error ) call prc_abort
1926 
1927  end subroutine atmos_saturation_tdew_liq_3d
1928 
1929  !-----------------------------------------------------------------------------
1936  !-----------------------------------------------------------------------------
1937  subroutine atmos_saturation_pote_0d( &
1938  DENS, POTT, TEMP, QV, &
1939  POTE )
1940  use scale_const, only: &
1941  rdry => const_rdry, &
1942  cpdry => const_cpdry
1943  use scale_atmos_hydrometeor, only: &
1944  atmos_hydrometeor_lhv
1945  real(RP), intent(in) :: DENS
1946  real(RP), intent(in) :: POTT
1947  real(RP), intent(in) :: TEMP
1948  real(RP), intent(in) :: QV
1949 
1950  real(RP), intent(out) :: POTE
1951 
1952  real(RP) :: TL
1953  real(RP) :: Pv
1954  real(RP) :: LHV
1955 
1956  pv = dens * qv * rvap * temp
1957  tl = 55.0_rp + 2840.0_rp / ( cpdry / rdry * log(temp) - log(pv) - 4.805_rp )
1958  call atmos_hydrometeor_lhv( temp, lhv ) ! [IN], [OUT]
1959 
1960  pote = pott * exp( lhv * qv / ( cpdry * temp ) &
1961  * 1.0784_rp * ( 1.0_rp + 0.810_rp * qv ) )
1962 
1963  return
1964  end subroutine atmos_saturation_pote_0d
1965 
1966 !OCL SERIAL
1967  subroutine atmos_saturation_pote_1d( &
1968  KA, KS, KE, &
1969  DENS, POTT, TEMP, QV, &
1970  POTE )
1971  integer, intent(in) :: KA, KS, KE
1972 
1973  real(RP), intent(in) :: DENS(KA)
1974  real(RP), intent(in) :: POTT(KA)
1975  real(RP), intent(in) :: TEMP(KA)
1976  real(RP), intent(in) :: QV (KA)
1977 
1978  real(RP), intent(out) :: POTE(KA)
1979 
1980  integer :: k
1981 
1982  do k = ks, ke
1983  call atmos_saturation_pote_0d( dens(k), pott(k), temp(k), qv(k), & ! [IN]
1984  pote(k) ) ! [OUT]
1985  end do
1986 
1987  return
1988  end subroutine atmos_saturation_pote_1d
1989 
1990  subroutine atmos_saturation_pote_3d( &
1991  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
1992  DENS, POTT, TEMP, QV, &
1993  POTE )
1994  integer, intent(in) :: KA, KS, KE
1995  integer, intent(in) :: IA, IS, IE
1996  integer, intent(in) :: JA, JS, JE
1997 
1998  real(RP), intent(in) :: DENS(KA,IA,JA)
1999  real(RP), intent(in) :: POTT(KA,IA,JA)
2000  real(RP), intent(in) :: TEMP(KA,IA,JA)
2001  real(RP), intent(in) :: QV (KA,IA,JA)
2002 
2003  real(RP), intent(out) :: POTE(KA,IA,JA)
2004 
2005  integer :: k, i, j
2006 
2007  !$omp parallel do
2008  do j = js, je
2009  do i = is, ie
2010  do k = ks, ke
2011  call atmos_saturation_pote_0d( &
2012  dens(k,i,j), pott(k,i,j), temp(k,i,j), qv(k,i,j), & ! [IN]
2013  pote(k,i,j) ) ! [OUT]
2014  end do
2015  end do
2016  end do
2017 
2018  return
2019  end subroutine atmos_saturation_pote_3d
2020 
2021  !-----------------------------------------------------------------------------
2023  !-----------------------------------------------------------------------------
2024 !OCL SERIAL
2025  subroutine atmos_saturation_moist_conversion_dens_liq_0d( &
2026  DENS, Emoist0, &
2027  TEMP, QV, QC, CPtot, CVtot, &
2028  converged )
2029  use scale_prc, only: &
2030  prc_abort
2031  use scale_atmos_hydrometeor, only: &
2032  cp_vapor, &
2033  cp_water, &
2034  cv_vapor, &
2035  cv_water, &
2036  lhv
2037  implicit none
2038 
2039  real(RP), intent(in) :: DENS
2040  real(RP), intent(in) :: Emoist0
2041 
2042  real(RP), intent(inout) :: TEMP
2043  real(RP), intent(inout) :: QV
2044  real(RP), intent(inout) :: QC
2045  real(RP), intent(inout) :: CPtot
2046  real(RP), intent(inout) :: CVtot
2047 
2048  logical, intent(out) :: converged
2049 
2050  ! working
2051  real(RP) :: QSUM
2052  real(RP) :: CVtot0
2053  real(RP) :: qsat
2054  real(RP) :: Emoist ! moist internal energy
2055 
2056  ! d(X)/dT
2057  real(RP) :: dqsatl_dT
2058  real(RP) :: dqc_dT
2059  real(RP) :: dCVtot_dT
2060  real(RP) :: dEmoist_dT
2061  real(RP) :: dtemp
2062 
2063  integer, parameter :: itelim = 100
2064  real(RP), parameter :: dtemp_criteria = 0.1_rp**(2+rp/2)
2065  integer :: ite
2066  !---------------------------------------------------------------------------
2067 
2068  qsum = qv + qc
2069 
2070  call atmos_saturation_dens2qsat_liq( temp, dens, & ! [IN]
2071  qsat ) ! [OUT]
2072 
2073  if ( qsum <= qsat ) then
2074  qv = qsum
2075  cptot = cptot + qc * ( cp_vapor - cp_water )
2076  cvtot = cvtot + qc * ( cv_vapor - cv_water )
2077  qc = 0.0_rp
2078  temp = ( emoist0 - lhv * qv ) / cvtot
2079  converged = .true.
2080  return
2081  end if
2082 
2083  cvtot0 = cvtot
2084 
2085  converged = .false.
2086  do ite = 1, itelim
2087 
2088  ! get qsat and dX/dT under the temp
2089  call atmos_saturation_dqs_dtem_dens_liq_0d( temp, dens, & ! [IN]
2090  dqsatl_dt, qsat ) ! [OUT]
2091 
2092  dqc_dt = - dqsatl_dt
2093 
2094  dcvtot_dt = dqsatl_dt * cv_vapor &
2095  + dqc_dt * cv_water
2096 
2097  demoist_dt = temp * dcvtot_dt + cvtot + dqsatl_dt * lhv
2098 
2099  ! diagnose quantities with the qsat
2100  qc = qsum - qsat
2101  cvtot = cvtot0 + ( cv_vapor - cv_water ) * ( qsat - qv )
2102  emoist = temp * cvtot + qsat * lhv
2103 
2104 
2105  ! update temp by the newtonian method
2106  dtemp = ( emoist - emoist0 ) / demoist_dt
2107  temp = temp - dtemp
2108 
2109  if ( abs(dtemp) < dtemp_criteria ) then
2110  converged = .true.
2111  exit
2112  endif
2113 
2114  if( temp*0.0_rp /= 0.0_rp ) exit
2115  enddo
2116 
2117  cptot = cptot + ( cp_vapor - cp_water ) * ( qsat - qv )
2118 
2119  qv = qsat
2120 
2121  return
2122  end subroutine atmos_saturation_moist_conversion_dens_liq_0d
2123 
2124  !-----------------------------------------------------------------------------
2126  !-----------------------------------------------------------------------------
2127 !OCL SERIAL
2129  DENS, Emoist0, &
2130  TEMP, QV, QC, QI, CPtot, CVtot, &
2131  converged )
2132  use scale_prc, only: &
2133  prc_abort
2134  use scale_atmos_hydrometeor, only: &
2135  cp_vapor, &
2136  cp_water, &
2137  cp_ice, &
2138  cv_vapor, &
2139  cv_water, &
2140  cv_ice, &
2141  lhv, &
2142  lhf
2143  implicit none
2144 
2145  real(RP), intent(in) :: DENS
2146  real(RP), intent(in) :: Emoist0
2147 
2148  real(RP), intent(inout) :: TEMP
2149  real(RP), intent(inout) :: QV
2150  real(RP), intent(inout) :: QC
2151  real(RP), intent(inout) :: QI
2152  real(RP), intent(inout) :: CPtot
2153  real(RP), intent(inout) :: CVtot
2154 
2155  logical, intent(out) :: converged
2156 
2157  ! working
2158  real(RP) :: QSUM
2159  real(RP) :: TEMP0
2160  real(RP) :: QV0, QC0, QI0
2161  real(RP) :: CVtot0
2162  real(RP) :: alpha
2163  real(RP) :: qsat
2164  real(RP) :: Emoist ! moist internal energy
2165 
2166  ! d(X)/dT
2167  real(RP) :: dalpha_dT
2168  real(RP) :: dqsat_dT
2169  real(RP) :: dqc_dT, dqi_dT
2170  real(RP) :: dCVtot_dT
2171  real(RP) :: dEmoist_dT
2172  real(RP) :: dtemp
2173 
2174  integer, parameter :: itelim = 100
2175  real(RP), parameter :: dtemp_criteria = 0.1_rp**(2+rp/2)
2176  integer :: ite
2177  !---------------------------------------------------------------------------
2178 
2179  temp0 = temp
2180  cvtot0 = cvtot
2181  qv0 = qv
2182  qc0 = qc
2183  qi0 = qi
2184  qsum = qv + qc + qi
2185 
2186  call atmos_saturation_dens2qsat_all( temp, dens, & ! [IN]
2187  qsat ) ! [OUT]
2188 
2189  if ( qsum <= qsat ) then
2190  qv = qsum
2191  cptot = cptot + qc * ( cp_vapor - cp_water ) + qi * ( cp_vapor - cp_ice )
2192  cvtot = cvtot + qc * ( cv_vapor - cv_water ) + qi * ( cv_vapor - cv_ice )
2193  qc = 0.0_rp
2194  qi = 0.0_rp
2195  temp = ( emoist0 - lhv * qv ) / cvtot
2196  converged = .true.
2197  return
2198  end if
2199 
2200  converged = .false.
2201  do ite = 1, itelim
2202 
2203  ! dX/dT
2204  call atmos_saturation_dqs_dtem_dens_all_0d( temp, dens, & ! [IN]
2205  dqsat_dt, & ! [OUT]
2206  qsat=qsat, alpha=alpha ) ! [OUT]
2207  call atmos_saturation_dalphadt( temp, dalpha_dt ) ! [IN], [OUT]
2208 
2209  dqc_dt = ( qsum - qv ) * dalpha_dt - dqsat_dt * ( alpha )
2210  dqi_dt = -( qsum - qv ) * dalpha_dt - dqsat_dt * ( 1.0_rp-alpha )
2211 
2212  dcvtot_dt = dqsat_dt * cv_vapor &
2213  + dqc_dt * cv_water &
2214  + dqi_dt * cv_ice
2215 
2216  demoist_dt = temp * dcvtot_dt + cvtot + dqsat_dt * lhv - dqi_dt * lhf
2217 
2218  ! Saturation
2219  qv = qsat
2220  qc = ( qsum - qsat ) * ( alpha )
2221  qi = ( qsum - qsat ) * ( 1.0_rp - alpha )
2222 
2223  cvtot = cvtot0 &
2224  + cv_vapor * ( qv - qv0 ) &
2225  + cv_water * ( qc - qc0 ) &
2226  + cv_ice * ( qi - qi0 )
2227 
2228  emoist = temp * cvtot + qv * lhv - qi * lhf
2229 
2230  ! update temp by the newtonian method
2231  dtemp = ( emoist - emoist0 ) / demoist_dt
2232  temp = temp - dtemp
2233 
2234  if ( abs(dtemp) < dtemp_criteria ) then
2235  converged = .true.
2236  exit
2237  endif
2238 
2239  if( temp*0.0_rp /= 0.0_rp ) exit
2240  enddo
2241 
2242  cptot = cptot &
2243  + cp_vapor * ( qv - qv0 ) &
2244  + cp_water * ( qc - qc0 ) &
2245  + cp_ice * ( qi - qi0 )
2246 
2247  return
2249 
2250  !-----------------------------------------------------------------------------
2252  !-----------------------------------------------------------------------------
2254  PRES, Entr, Qdry, &
2255  QV, QC, &
2256  Rtot, CPtot, &
2257  TEMP, &
2258  converged )
2259  use scale_const, only: &
2260  eps => const_eps, &
2261  lhv0 => const_lhv0
2262  use scale_atmos_hydrometeor, only: &
2263  atmos_hydrometeor_entr, &
2264  atmos_hydrometeor_entr2temp, &
2265  cp_vapor, &
2266  cp_water
2267 
2268  real(RP), intent(in) :: PRES
2269  real(RP), intent(in) :: ENTR
2270  real(RP), intent(in) :: Qdry
2271 
2272  real(RP), intent(inout) :: QV
2273  real(RP), intent(inout) :: QC
2274  real(RP), intent(inout) :: CPtot
2275  real(RP), intent(inout) :: Rtot
2276 
2277  real(RP), intent(out) :: TEMP
2278  logical, intent(out) :: converged
2279 
2280  real(RP), parameter :: TEMMIN = 0.1_rp
2281  real(RP), parameter :: criteria = 1.e-8_rp
2282  integer, parameter :: itelim = 100
2283  integer :: ite
2284 
2285  real(RP) :: Qsum
2286  real(RP) :: qsat, psat
2287 
2288  real(RP) :: TEMP_ite
2289  real(RP) :: QV_ite
2290  real(RP) :: ENTR_ite
2291  real(RP) :: Rtot_ite
2292  real(RP) :: CPtot_ite
2293 
2294  real(RP) :: dqsat_dT, dqsat_dP
2295 
2296  real(RP) :: TEMP_prev
2297  real(RP) :: dENTR_dT
2298  !---------------------------------------------------------------------------
2299 
2300  qsum = qv + qc
2301 
2302  ! all liquid evaporates
2303  rtot = rtot + qc * rvap
2304  cptot = cptot + qc * ( cp_vapor - cp_water )
2305  qv = qsum
2306  qc = 0.0_rp
2307 
2308  call atmos_hydrometeor_entr2temp( entr, pres, qsum, 0.0_rp, qdry, & ! [IN]
2309  rtot, cptot, & ! [IN]
2310  temp ) ! [OUT]
2311 
2312  call atmos_saturation_pres2qsat_liq( temp, pres, qdry, & ! [IN]
2313  qsat ) ! [OUT]
2314  if ( qsum <= qsat ) then
2315  ! unsaturated
2316 
2317  converged = .true.
2318 
2319  return
2320 
2321  else
2322  ! saturated
2323 
2324  temp_ite = temp
2325 
2326  do ite = 1, itelim
2327 
2328  call atmos_saturation_dqs_dtem_dpre_liq_0d( temp_ite, pres, qdry, & ! [IN]
2329  dqsat_dt, dqsat_dp, & ! [OUT]
2330  qsat=qsat, psat=psat ) ! [OUT]
2331 
2332  qv_ite = min( qsum, qsat )
2333 
2334  rtot_ite = rtot - ( qsum - qv_ite ) * rvap
2335  cptot_ite = cptot - ( qsum - qv_ite ) * ( cp_vapor - cp_water )
2336 
2337  dentr_dt = cptot_ite / temp_ite &
2338  + ( ( cp_vapor - cp_water ) * log( temp_ite / tem00 ) &
2339  - rvap * log( psat/ psat0 ) &
2340  + lhv0 / tem00 &
2341  ) * dqsat_dt
2342 
2343  call atmos_hydrometeor_entr( temp_ite, pres, & ! [IN]
2344  qv_ite, 0.0_rp, qdry, & ! [IN] ! QI = 0
2345  rtot_ite, cptot_ite, & ! [IN]
2346  entr_ite ) ! [OUT]
2347 
2348  temp_prev = temp_ite
2349  temp_ite = temp_ite - ( entr_ite - entr ) / max( dentr_dt, eps )
2350  temp_ite = max( temp_ite, temmin )
2351 
2352  if( abs(temp_ite-temp_prev) < criteria ) then
2353  converged = .true.
2354  exit
2355  end if
2356 
2357  enddo
2358 
2359  end if
2360 
2361  qv = qv_ite
2362  qc = qsum - qv
2363  cptot = cptot_ite
2364  rtot = rtot_ite
2365  temp = temp_ite
2366 
2367  return
2369 
2370 end module scale_atmos_saturation
scale_const::const_lhv0
real(rp), parameter, public const_lhv0
latent heat of vaporizaion at 0C [J/kg]
Definition: scale_const.F90:75
scale_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:342
scale_atmos_hydrometeor::cp_water
real(rp), public cp_water
CP for water [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:133
scale_const::const_lhs00
real(rp), public const_lhs00
latent heat of sublimation at 0K [J/kg]
Definition: scale_const.F90:78
scale_atmos_saturation::atmos_saturation_moist_conversion_pres_liq_0d
subroutine atmos_saturation_moist_conversion_pres_liq_0d(PRES, Entr, Qdry, QV, QC, Rtot, CPtot, TEMP, converged)
Iterative moist conversion for liquid water at constant pressure.
Definition: scale_atmos_saturation.F90:2259
scale_const::const_epsvap
real(rp), public const_epsvap
Rdry / Rvap.
Definition: scale_const.F90:69
scale_const::const_lhv00
real(rp), public const_lhv00
latent heat of vaporizaion at 0K [J/kg]
Definition: scale_const.F90:76
scale_precision
module PRECISION
Definition: scale_precision.F90:14
scale_const::const_rvap
real(rp), parameter, public const_rvap
specific gas constant (water vapor) [J/kg/K]
Definition: scale_const.F90:63
scale_atmos_saturation::atmos_saturation_dqs_dtem_dpre_ice_1d
subroutine atmos_saturation_dqs_dtem_dpre_ice_1d(KA, KS, KE, temp, pres, qdry, dqsat_dT, dqsat_dP)
Definition: scale_atmos_saturation.F90:1690
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:33
scale_atmos_hydrometeor
module atmosphere / hydrometeor
Definition: scale_atmos_hydrometeor.F90:12
scale_atmos_saturation::atmos_saturation_moist_conversion_dens_all_0d
subroutine atmos_saturation_moist_conversion_dens_all_0d(DENS, Emoist0, TEMP, QV, QC, QI, CPtot, CVtot, converged)
Iterative moist conversion (liquid/ice mixture) at constant density (volume)
Definition: scale_atmos_saturation.F90:2132
scale_const::const_thermodyn_type
character(len=h_short), public const_thermodyn_type
internal energy type
Definition: scale_const.F90:99
scale_const::const_cpvap
real(rp), parameter, public const_cpvap
specific heat (water vapor, constant pressure) [J/kg/K]
Definition: scale_const.F90:64
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:130
scale_io
module STDIO
Definition: scale_io.F90:10
scale_atmos_saturation::atmos_saturation_tdew_liq_1d
subroutine atmos_saturation_tdew_liq_1d(KA, KS, KE, DENS, TEMP, QV, Tdew, converged)
Definition: scale_atmos_saturation.F90:1868
scale_const
module CONSTANT
Definition: scale_const.F90:11
scale_atmos_hydrometeor::lhv
real(rp), public lhv
latent heat of vaporization for use [J/kg]
Definition: scale_atmos_hydrometeor.F90:126
scale_atmos_saturation::atmos_saturation_dqs_dtem_dens_ice_1d
subroutine atmos_saturation_dqs_dtem_dens_ice_1d(KA, KS, KE, temp, dens, dqsdtem)
Definition: scale_atmos_saturation.F90:1446
scale_const::const_cpdry
real(rp), public const_cpdry
specific heat (dry air,constant pressure) [J/kg/K]
Definition: scale_const.F90:56
scale_prof
module profiler
Definition: scale_prof.F90:11
scale_const::const_psat0
real(rp), parameter, public const_psat0
saturate pressure of water vapor at 0C [Pa]
Definition: scale_const.F90:81
scale_atmos_saturation::atmos_saturation_setup
subroutine, public atmos_saturation_setup
Setup.
Definition: scale_atmos_saturation.F90:215
scale_const::const_ci
real(rp), parameter, public const_ci
specific heat (ice) [J/kg/K]
Definition: scale_const.F90:67
scale_atmos_saturation::atmos_saturation_pote_0d
subroutine atmos_saturation_pote_0d(DENS, POTT, TEMP, QV, POTE)
calculate equivalent potential temperature Bolton, D., 1980: The computation of equivalent potential ...
Definition: scale_atmos_saturation.F90:1940
scale_const::const_tem00
real(rp), parameter, public const_tem00
temperature reference (0C) [K]
Definition: scale_const.F90:90
scale_const::const_cl
real(rp), parameter, public const_cl
specific heat (liquid water) [J/kg/K]
Definition: scale_const.F90:66
scale_const::const_cvvap
real(rp), public const_cvvap
specific heat (water vapor, constant volume) [J/kg/K]
Definition: scale_const.F90:65
scale_atmos_saturation::atmos_saturation_pote_1d
subroutine atmos_saturation_pote_1d(KA, KS, KE, DENS, POTT, TEMP, QV, POTE)
Definition: scale_atmos_saturation.F90:1971
scale_atmos_hydrometeor::lhf
real(rp), public lhf
latent heat of fusion for use [J/kg]
Definition: scale_atmos_hydrometeor.F90:128
scale_const::const_rdry
real(rp), public const_rdry
specific gas constant (dry air) [J/kg/K]
Definition: scale_const.F90:55
scale_atmos_saturation::atmos_saturation_alpha_0d
subroutine atmos_saturation_alpha_0d(temp, alpha)
calc liquid/ice separation factor (0D)
Definition: scale_atmos_saturation.F90:296
scale_atmos_saturation::atmos_saturation_dqs_dtem_dpre_liq_1d
subroutine atmos_saturation_dqs_dtem_dpre_liq_1d(KA, KS, KE, temp, pres, qdry, dqsat_dT, dqsat_dP)
Definition: scale_atmos_saturation.F90:1588
scale_const::const_lhs0
real(rp), parameter, public const_lhs0
latent heat of sublimation at 0C [J/kg]
Definition: scale_const.F90:77
scale_atmos_saturation::atmos_saturation_dqs_dtem_dens_liq_1d
subroutine atmos_saturation_dqs_dtem_dens_liq_1d(KA, KS, KE, temp, dens, dqsdtem)
Definition: scale_atmos_saturation.F90:1356
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:41
scale_atmos_hydrometeor::atmos_hydrometeor_setup
subroutine, public atmos_hydrometeor_setup
Setup.
Definition: scale_atmos_hydrometeor.F90:154
scale_atmos_saturation
module atmosphere / saturation
Definition: scale_atmos_saturation.F90:12
scale_io::io_fid_conf
integer, public io_fid_conf
Config file ID.
Definition: scale_io.F90:56
scale_atmos_hydrometeor::cp_vapor
real(rp), public cp_vapor
CP for vapor [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:131
scale_atmos_hydrometeor::cv_water
real(rp), public cv_water
CV for water [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:132
scale_atmos_hydrometeor::cv_ice
real(rp), public cv_ice
CV for ice [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:134
scale_atmos_hydrometeor::cp_ice
real(rp), public cp_ice
CP for ice [J/kg/K].
Definition: scale_atmos_hydrometeor.F90:135