Line data Source code
1 : !-------------------------------------------------------------------------------
2 : ! $Id$
3 : !===============================================================================
4 : module parameters_model
5 :
6 : ! Description:
7 : ! Contains model parameters that are determined at run time rather than
8 : ! compile time.
9 : !
10 : ! References:
11 : ! None
12 : !-------------------------------------------------------------------------------
13 :
14 : use clubb_precision, only: &
15 : core_rknd
16 :
17 : implicit none
18 :
19 : private ! Default scope
20 :
21 : integer, parameter :: &
22 : sp = selected_real_kind(6) ! 32-bit floating point number
23 :
24 : ! Maximum magnitude of PDF parameter 'mixt_frac'.
25 : real( kind = core_rknd ), public :: mixt_frac_max_mag
26 :
27 : !$omp threadprivate(mixt_frac_max_mag)
28 :
29 : ! Model parameters and constraints setup in the namelists
30 : real( kind = core_rknd ), public :: &
31 : T0 = 300._core_rknd, & ! Reference temperature (usually 300) [K]
32 : ts_nudge = 0._core_rknd ! Timescale of u/v nudging [s]
33 :
34 : #ifdef GFDL
35 : real( kind = core_rknd ), public :: & ! h1g, 2010-06-15
36 : cloud_frac_min ! minimum cloud fraction for droplet #
37 : !$omp threadprivate( cloud_frac_min )
38 : #endif
39 :
40 :
41 : !$omp threadprivate(T0, ts_nudge)
42 :
43 : real( kind = core_rknd), public :: &
44 : rtm_min = epsilon( rtm_min ), & ! Value below which rtm will be nudged [kg/kg]
45 : rtm_nudge_max_altitude = 10000._core_rknd ! Highest altitude at which to nudge rtm [m]
46 : !$omp threadprivate( rtm_min, rtm_nudge_max_altitude )
47 :
48 : integer, public :: &
49 : sclr_dim = 0, & ! Number of passive scalars
50 : edsclr_dim = 0, & ! Number of eddy-diff. passive scalars
51 : hydromet_dim = 0 ! Number of hydrometeor species
52 :
53 : !$omp threadprivate(sclr_dim, edsclr_dim, hydromet_dim)
54 :
55 : real( kind = core_rknd ), dimension(:), allocatable, public :: &
56 : sclr_tol ! Threshold(s) on the passive scalars [units vary]
57 :
58 : !$acc declare create(sclr_tol)
59 : !$omp threadprivate(sclr_tol)
60 :
61 : real( kind = sp ), public :: PosInf
62 :
63 : !$omp threadprivate(PosInf)
64 :
65 : public :: setup_parameters_model
66 :
67 : contains
68 :
69 : !-------------------------------------------------------------------------------
70 1536 : subroutine setup_parameters_model &
71 : ( T0_in, ts_nudge_in, Skw_max_mag, &
72 : hydromet_dim_in, &
73 1536 : sclr_dim_in, sclr_tol_in, edsclr_dim_in &
74 : #ifdef GFDL
75 : , cloud_frac_min_in & ! hlg, 2010-6-15
76 : #endif
77 :
78 : )
79 :
80 : ! Description:
81 : ! Sets parameters to their initial values
82 : !
83 : ! References:
84 : ! None
85 : !-------------------------------------------------------------------------------
86 :
87 : use clubb_precision, only: &
88 : core_rknd ! Variable(s)
89 :
90 : implicit none
91 :
92 : ! External
93 : intrinsic :: sqrt, allocated, transfer
94 :
95 : ! Constants
96 : integer(kind=4), parameter :: nanbits = 2139095040
97 :
98 : ! Input Variables
99 : real( kind = core_rknd ), intent(in) :: &
100 : T0_in, & ! Ref. temperature [K]
101 : ts_nudge_in, & ! Timescale for u/v nudging [s]
102 : Skw_max_mag ! Maximum allowable magnitude of Skewness [-]
103 :
104 : #ifdef GFDL
105 : real( kind = core_rknd ), intent(in) :: cloud_frac_min_in ! h1g, 2010-06-15
106 : #endif
107 :
108 :
109 : integer, intent(in) :: &
110 : hydromet_dim_in, & ! Number of hydrometeor species
111 : sclr_dim_in, & ! Number of passive scalars
112 : edsclr_dim_in ! Number of eddy-diff. passive scalars
113 :
114 : real( kind = core_rknd ), intent(in), dimension(sclr_dim_in) :: &
115 : sclr_tol_in ! Threshold on passive scalars
116 :
117 : ! --- Begin Code ---
118 :
119 : ! Formula from subroutine pdf_closure, where sigma_sqd_w = 0.4 and Skw =
120 : ! Skw_max_mag in this formula. Note that this is constant, but can't appear
121 : ! with a Fortran parameter attribute, so we define it here.
122 : mixt_frac_max_mag = 1.0_core_rknd &
123 : - ( 0.5_core_rknd * ( 1.0_core_rknd - Skw_max_mag / &
124 : sqrt( 4.0_core_rknd * ( 1.0_core_rknd - 0.4_core_rknd )**3 &
125 1536 : + Skw_max_mag**2 ) ) ) ! Known magic number
126 :
127 1536 : T0 = T0_in
128 1536 : ts_nudge = ts_nudge_in
129 :
130 1536 : hydromet_dim = hydromet_dim_in
131 1536 : sclr_dim = sclr_dim_in
132 1536 : edsclr_dim = edsclr_dim_in
133 :
134 : ! In a tuning run, this array has the potential to be allocated already
135 1536 : if ( .not. allocated( sclr_tol ) ) then
136 4608 : allocate( sclr_tol(1:max(1,sclr_dim)) )
137 : else
138 0 : deallocate( sclr_tol )
139 0 : allocate( sclr_tol(1:max(1,sclr_dim)) )
140 : end if
141 :
142 : ! Set to 0 before setting to sclr_tol_in because the above allocates sclr_tol
143 : ! with a minimum dimension of 1, while it may be the case that sclr_dim=0.
144 : ! This ensures that if we are not using scalars that sclr_tol will be
145 : ! initialized to 0.0
146 3072 : sclr_tol = 0.0_core_rknd
147 1536 : sclr_tol(1:sclr_dim) = sclr_tol_in(1:sclr_dim)
148 :
149 1536 : PosInf = transfer( nanbits, PosInf )
150 :
151 : #ifdef GFDL
152 : cloud_frac_min = cloud_frac_min_in ! h1g, 2010-06-15
153 : #endif
154 :
155 1536 : return
156 : end subroutine setup_parameters_model
157 : !-------------------------------------------------------------------------------
158 :
159 : end module parameters_model
|