CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
module_ccpp_suite_simulator.F90
1
4
6
7 use machine, only : kind_phys
8 implicit none
9
11
14 real(kind_phys), dimension(:), allocatable :: t
15 real(kind_phys), dimension(:), allocatable :: u
16 real(kind_phys), dimension(:), allocatable :: v
17 real(kind_phys), dimension(:), allocatable :: q
18 real(kind_phys), dimension(:), allocatable :: p
19 real(kind_phys), dimension(:), allocatable :: z
20 end type phys_tend_1d
21
24 real(kind_phys), dimension(:), allocatable :: time
25 real(kind_phys), dimension(:,:), allocatable :: t
26 real(kind_phys), dimension(:,:), allocatable :: u
27 real(kind_phys), dimension(:,:), allocatable :: v
28 real(kind_phys), dimension(:,:), allocatable :: q
29 real(kind_phys), dimension(:,:), allocatable :: p
30 real(kind_phys), dimension(:,:), allocatable :: z
31 end type phys_tend_2d
32
33 ! Type containing 3D (loc,lev,time) physics tendencies.
35 real(kind_phys), dimension(:), allocatable :: time
36 real(kind_phys), dimension(:), allocatable :: lon
37 real(kind_phys), dimension(:), allocatable :: lat
38 real(kind_phys), dimension(:,:,:), allocatable :: t
39 real(kind_phys), dimension(:,:,:), allocatable :: u
40 real(kind_phys), dimension(:,:,:), allocatable :: v
41 real(kind_phys), dimension(:,:,:), allocatable :: q
42 end type phys_tend_3d
43
46 real(kind_phys), dimension(:), allocatable :: time
47 real(kind_phys), dimension(:,:), allocatable :: lon
48 real(kind_phys), dimension(:,:), allocatable :: lat
49 real(kind_phys), dimension(:,:,:,:), allocatable :: t
50 real(kind_phys), dimension(:,:,:,:), allocatable :: u
51 real(kind_phys), dimension(:,:,:,:), allocatable :: v
52 real(kind_phys), dimension(:,:,:,:), allocatable :: q
53 end type phys_tend_4d
54
57 character(len=16) :: name
58 logical :: time_split = .false.
59 logical :: use_sim = .false.
60 integer :: order
61 type(phys_tend_1d) :: tend1d
62 type(phys_tend_2d) :: tend2d
63 type(phys_tend_3d) :: tend3d
64 type(phys_tend_4d) :: tend4d
65 character(len=16) :: active_name
66 integer :: iactive_scheme
67 logical :: active_tsp
68 integer :: nprg_active
69 contains
70 generic, public :: linterp => linterp_1d, linterp_2d
71 procedure, private :: linterp_1d
72 procedure, private :: linterp_2d
73 procedure, public :: find_nearest_loc_2d_1d
74 procedure, public :: cmp_time_wts
76
77contains
78
81 function linterp_1d(this, var_name, year, month, day, hour, min, sec) result(err_message)
82 class(base_physics_process), intent(inout) :: this
83 character(len=*), intent(in) :: var_name
84 integer, intent(in) :: year, month, day, hour, min, sec
85 character(len=128) :: err_message
86 integer :: ti(1), tf(1), ntime
87 real(kind_phys) :: w1, w2
88
89 ! Interpolation weights
90 call this%cmp_time_wts(year, month, day, hour, min, sec, w1, w2, ti, tf)
91
92 ntime = size(this%tend2d%T(1,:))
93
94 select case(var_name)
95 case("T")
96 if (tf(1) .le. ntime) then
97 this%tend1d%T = w1*this%tend2d%T(:,ti(1)) + w2*this%tend2d%T(:,tf(1))
98 else
99 this%tend1d%T = this%tend2d%T(:,1)
100 endif
101 case("u")
102 if (tf(1) .le. ntime) then
103 this%tend1d%u = w1*this%tend2d%u(:,ti(1)) + w2*this%tend2d%u(:,tf(1))
104 else
105 this%tend1d%u = this%tend2d%u(:,1)
106 endif
107 case("v")
108 if (tf(1) .le. ntime) then
109 this%tend1d%v = w1*this%tend2d%v(:,ti(1)) + w2*this%tend2d%v(:,tf(1))
110 else
111 this%tend1d%v = this%tend2d%v(:,1)
112 endif
113 case("q")
114 if (tf(1) .le. ntime) then
115 this%tend1d%q = w1*this%tend2d%q(:,ti(1)) + w2*this%tend2d%q(:,tf(1))
116 else
117 this%tend1d%q = this%tend2d%q(:,1)
118 endif
119 end select
120
121 end function linterp_1d
122
127 function linterp_2d(this, var_name, lon, lat, year, month, day, hour, min, sec) result(err_message)
128 class(base_physics_process), intent(inout) :: this
129 character(len=*), intent(in) :: var_name
130 integer, intent(in) :: year, month, day, hour, min, sec
131 real(kind_phys), intent(in) :: lon, lat
132 character(len=128) :: err_message
133 integer :: ti(1), tf(1), inearest
134 real(kind_phys) :: w1, w2
135
136 ! Interpolation weights (temporal)
137 call this%cmp_time_wts(year, month, day, hour, min, sec, w1, w2, ti, tf)
138
139 ! Grab data tendency closest to column [lon,lat]
140 inearest = this%find_nearest_loc_2d_1d(lon,lat)
141
142 select case(var_name)
143 case("T")
144 this%tend1d%T = w1*this%tend3d%T(inearest,:,ti(1)) + w2*this%tend3d%T(inearest,:,tf(1))
145 case("u")
146 this%tend1d%u = w1*this%tend3d%u(inearest,:,ti(1)) + w2*this%tend3d%u(inearest,:,tf(1))
147 case("v")
148 this%tend1d%v = w1*this%tend3d%v(inearest,:,ti(1)) + w2*this%tend3d%v(inearest,:,tf(1))
149 case("q")
150 this%tend1d%q = w1*this%tend3d%q(inearest,:,ti(1)) + w2*this%tend3d%q(inearest,:,tf(1))
151 end select
152 end function linterp_2d
153
156 pure function find_nearest_loc_2d_1d(this, lon, lat)
157 class(base_physics_process), intent(in) :: this
158 real(kind_phys), intent(in) :: lon, lat
159 integer :: find_nearest_loc_2d_1d
160
161 find_nearest_loc_2d_1d = 1
162 end function find_nearest_loc_2d_1d
163
166 subroutine cmp_time_wts(this, year, month, day, hour, minute, sec, w1, w2, ti, tf)
167 ! Inputs
168 class(base_physics_process), intent(in) :: this
169 integer, intent(in) :: year, month, day, hour, minute, sec
170 ! Outputs
171 integer,intent(out) :: ti(1), tf(1)
172 real(kind_phys),intent(out) :: w1, w2
173 ! Locals
174 real(kind_phys) :: hrofday
175
176 hrofday = hour*3600. + minute*60. + sec
177 ti = max(hour,1)
178 tf = min(ti + 1,24)
179 w1 = ((hour+1)*3600 - hrofday)/3600
180 w2 = 1 - w1
181
182 end subroutine cmp_time_wts
183
185 subroutine sim_lwrad( year, month, day, hour, min, sec, process)
186 type(base_physics_process), intent(inout) :: process
187 integer, intent(in) :: year, month, day, hour, min, sec
188 character(len=128) :: errmsg
189
190 if (allocated(process%tend2d%T)) then
191 errmsg = process%linterp("T", year,month,day,hour,min,sec)
192 endif
193
194 end subroutine sim_lwrad
195
197 subroutine sim_swrad( year, month, day, hour, min, sec, process)
198 type(base_physics_process), intent(inout) :: process
199 integer, intent(in) :: year, month, day, hour, min, sec
200 character(len=128) :: errmsg
201
202 if (allocated(process%tend2d%T)) then
203 errmsg = process%linterp("T", year,month,day,hour,min,sec)
204 endif
205
206 end subroutine sim_swrad
207
209 subroutine sim_gwd( year, month, day, hour, min, sec, process)
210 type(base_physics_process), intent(inout) :: process
211 integer, intent(in) :: year, month, day, hour, min, sec
212 character(len=128) :: errmsg
213
214 if (allocated(process%tend2d%T)) then
215 errmsg = process%linterp("T", year,month,day,hour,min,sec)
216 endif
217 if (allocated(process%tend2d%u)) then
218 errmsg = process%linterp("u", year,month,day,hour,min,sec)
219 endif
220 if (allocated(process%tend2d%v)) then
221 errmsg = process%linterp("v", year,month,day,hour,min,sec)
222 endif
223
224 end subroutine sim_gwd
225
227 subroutine sim_pbl( year, month, day, hour, min, sec, process)
228 type(base_physics_process), intent(inout) :: process
229 integer, intent(in) :: year, month, day, hour, min, sec
230 character(len=128) :: errmsg
231
232 if (allocated(process%tend2d%T)) then
233 errmsg = process%linterp("T", year,month,day,hour,min,sec)
234 endif
235 if (allocated(process%tend2d%u)) then
236 errmsg = process%linterp("u", year,month,day,hour,min,sec)
237 endif
238 if (allocated(process%tend2d%v)) then
239 errmsg = process%linterp("v", year,month,day,hour,min,sec)
240 endif
241 if (allocated(process%tend2d%q)) then
242 errmsg = process%linterp("q", year,month,day,hour,min,sec)
243 endif
244
245 end subroutine sim_pbl
246
248 subroutine sim_dcnv( year, month, day, hour, min, sec, process)
249 type(base_physics_process), intent(inout) :: process
250 integer, intent(in) :: year, month, day, hour, min, sec
251 character(len=128) :: errmsg
252
253 if (allocated(process%tend2d%T)) then
254 errmsg = process%linterp("T", year,month,day,hour,min,sec)
255 endif
256 if (allocated(process%tend2d%u)) then
257 errmsg = process%linterp("u", year,month,day,hour,min,sec)
258 endif
259 if (allocated(process%tend2d%v)) then
260 errmsg = process%linterp("v", year,month,day,hour,min,sec)
261 endif
262 if (allocated(process%tend2d%q)) then
263 errmsg = process%linterp("q", year,month,day,hour,min,sec)
264 endif
265
266 end subroutine sim_dcnv
267
269 subroutine sim_scnv( year, month, day, hour, min, sec, process)
270 type(base_physics_process), intent(inout) :: process
271 integer, intent(in) :: year, month, day, hour, min, sec
272 character(len=128) :: errmsg
273
274 if (allocated(process%tend2d%T)) then
275 errmsg = process%linterp("T", year,month,day,hour,min,sec)
276 endif
277 if (allocated(process%tend2d%u)) then
278 errmsg = process%linterp("u", year,month,day,hour,min,sec)
279 endif
280 if (allocated(process%tend2d%v)) then
281 errmsg = process%linterp("v", year,month,day,hour,min,sec)
282 endif
283 if (allocated(process%tend2d%q)) then
284 errmsg = process%linterp("q", year,month,day,hour,min,sec)
285 endif
286
287 end subroutine sim_scnv
288
290 subroutine sim_cldmp( year, month, day, hour, min, sec, process)
291 type(base_physics_process), intent(inout) :: process
292 integer, intent(in) :: year, month, day, hour, min, sec
293 character(len=128) :: errmsg
294
295 if (allocated(process%tend2d%T)) then
296 errmsg = process%linterp("T", year,month,day,hour,min,sec)
297 endif
298 if (allocated(process%tend2d%q)) then
299 errmsg = process%linterp("q", year,month,day,hour,min,sec)
300 endif
301 end subroutine sim_cldmp
302
This type contains the meta information and data for each physics process.
Type containing 1D (time) physics tendencies.
Type containing 2D (lev,time) physics tendencies.
Type containing 4D (lon,lat,lev,time) physics tendencies.