CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_ccpp_suite_sim_pre.F90
1
3
4! ########################################################################################
5!
6! Description: Interstitial CCPP suite to couple UFS physics to ccpp_suite_simulator.
7!
8! Contains:
9! - load_ccpp_suite_sim(): read and load data into type used by ccpp_suite_simulator.
10! called once during model initialization
11! - GFS_ccpp_suite_sim_pre_run(): prepare GFS diagnostic physics tendencies for
12! ccpp_suite_simulator.
13!
14! ########################################################################################
16 use machine, only: kind_phys
18 use netcdf
19 implicit none
20 public gfs_ccpp_suite_sim_pre_run, load_ccpp_suite_sim
21contains
22
26 subroutine gfs_ccpp_suite_sim_pre_run(do_ccpp_suite_sim, dtend, ntqv, dtidx, dtp, &
27 index_of_process_dcnv, index_of_process_longwave, index_of_process_shortwave, &
28 index_of_process_scnv, index_of_process_orographic_gwd, index_of_process_pbl, &
29 index_of_process_mp, index_of_temperature, index_of_x_wind, index_of_y_wind, &
30 physics_process, iactive_T, iactive_u, iactive_v, iactive_q, active_phys_tend, &
31 errmsg, errflg)
32
33 ! Inputs
34 logical, intent(in) :: do_ccpp_suite_sim
35 integer, intent(in) :: ntqv, index_of_process_dcnv, index_of_process_longwave, &
36 index_of_process_shortwave, index_of_process_scnv, &
37 index_of_process_orographic_gwd, index_of_process_pbl, index_of_process_mp, &
38 index_of_temperature, index_of_x_wind, index_of_y_wind
39 integer, intent(in), dimension(:,:) :: dtidx
40 real(kind_phys), intent(in) :: dtp
41 real(kind_phys), intent(in), dimension(:,:,:), optional :: dtend
42 type(base_physics_process),intent(in) :: physics_process(:)
43 integer, intent(in) :: iactive_t, iactive_u, iactive_v, iactive_q
44
45 ! Outputs
46 real(kind_phys), intent(out) :: active_phys_tend(:,:,:)
47 character(len=*),intent(out) :: errmsg
48 integer, intent(out) :: errflg
49
50 ! Locals
51 integer :: idtend, iactive
52
53 ! Initialize CCPP error handling variables
54 errmsg = ''
55 errflg = 0
56
57 if (.not. do_ccpp_suite_sim) return
58
59 ! Get tendency for "active" process.
60
61 ! ######################################################################################
62 ! DJS2023: For the UFS and SCM, the physics tendencies are stored in a multi-dimensional
63 ! array, CCPP standard_name = cumulative_change_of_state_variables.
64 ! These are not the instantaneous physics tendencies that are applied to the state by
65 ! the physics suites. Not all suites output physics tendencies...
66 ! Rather these are intended for diagnostic puposes and are accumulated over some
67 ! interval.
68 ! In the UFS/SCM this is controlled by the diagnostic bucket interval, namelist option
69 ! "fhzero". For this to work, you need to clear the diagnostic buckets after each
70 ! physics timestep when running in the UFS/SCM.
71 !
72 ! In the SCM this is done by adding the following runtime options:
73 ! --n_itt_out 1 --n_itt_diag 1
74 !
75 ! ######################################################################################
76 if (physics_process(1)%active_name == "LWRAD") iactive = index_of_process_longwave
77 if (physics_process(1)%active_name == "SWRAD") iactive = index_of_process_shortwave
78 if (physics_process(1)%active_name == "PBL") iactive = index_of_process_pbl
79 if (physics_process(1)%active_name == "GWD") iactive = index_of_process_orographic_gwd
80 if (physics_process(1)%active_name == "SCNV") iactive = index_of_process_scnv
81 if (physics_process(1)%active_name == "DCNV") iactive = index_of_process_dcnv
82 if (physics_process(1)%active_name == "cldMP") iactive = index_of_process_mp
83
84 ! Heat
85 idtend = dtidx(index_of_temperature,iactive)
86 if (idtend >= 1) then
87 active_phys_tend(:,:,iactive_t) = dtend(:,:,idtend)/dtp
88 endif
89
90 ! u-wind
91 idtend = dtidx(index_of_x_wind,iactive)
92 if (idtend >= 1) then
93 active_phys_tend(:,:,iactive_u) = dtend(:,:,idtend)/dtp
94 endif
95
96 ! v-wind
97 idtend = dtidx(index_of_y_wind,iactive)
98 if (idtend >= 1) then
99 active_phys_tend(:,:,iactive_v) = dtend(:,:,idtend)/dtp
100 endif
101
102 ! Moisture
103 idtend = dtidx(100+ntqv,iactive)
104 if (idtend >= 1) then
105 active_phys_tend(:,:,iactive_q) = dtend(:,:,idtend)/dtp
106 endif
107
108 end subroutine gfs_ccpp_suite_sim_pre_run
109
110 ! ######################################################################################
112 subroutine load_ccpp_suite_sim(nlunit, nml_file, physics_process, iactive_T, &
113 iactive_u, iactive_v, iactive_q, errmsg, errflg)
114
115 ! Inputs
116 integer, intent (in) :: nlunit
117 character(len=*), intent (in) :: nml_file
118
119 ! Outputs
120 type(base_physics_process),intent(inout),allocatable :: physics_process(:)
121 integer, intent(inout) :: iactive_t, iactive_u, iactive_v, iactive_q
122 integer, intent(out) :: errflg
123 character(len=256), intent(out) :: errmsg
124
125 ! Local variables
126 integer :: ncid, dimid, varid, status, ios, iprc, nlev_data, ntime_data
127 character(len=256) :: suite_sim_file
128 logical :: exists, do_ccpp_suite_sim
129 integer :: nprc_sim
130
131 ! For each process there is a corresponding namelist entry, which is constructed as
132 ! follows:
133 ! {use_suite_sim[0(no)/1(yes)], time_split[0(no)/1(yes)], order[1:nPhysProcess]}
134 integer, dimension(3) :: &
135 prc_lwrad_cfg = (/0,0,0/), &
136 prc_swrad_cfg = (/0,0,0/), &
137 prc_pbl_cfg = (/0,0,0/), &
138 prc_gwd_cfg = (/0,0,0/), &
139 prc_scnv_cfg = (/0,0,0/), &
140 prc_dcnv_cfg = (/0,0,0/), &
141 prc_cldmp_cfg = (/0,0,0/)
142
143 ! Namelist
144 namelist / ccpp_suite_sim_nml / do_ccpp_suite_sim, suite_sim_file, nprc_sim, &
145 prc_lwrad_cfg, prc_swrad_cfg, prc_pbl_cfg, prc_gwd_cfg, prc_scnv_cfg, &
146 prc_dcnv_cfg, prc_cldmp_cfg
147
148 errmsg = ''
149 errflg = 0
150
151 ! Read in namelist
152 inquire (file = trim(nml_file), exist = exists)
153 if (.not. exists) then
154 errmsg = 'CCPP suite simulator namelist file: '//trim(nml_file)//' does not exist.'
155 errflg = 1
156 return
157 else
158 open (unit = nlunit, file = nml_file, action = 'read', status = 'old', iostat = ios)
159 endif
160 rewind(nlunit)
161 read (nlunit, nml = ccpp_suite_sim_nml, iostat=status)
162 close (nlunit)
163
164 ! Only proceed if suite simulator requested.
165 if (prc_swrad_cfg(1) == 1 .or. prc_lwrad_cfg(1) == 1 .or. prc_pbl_cfg(1) == 1 .or. &
166 prc_gwd_cfg(1) == 1 .or. prc_scnv_cfg(1) == 1 .or. prc_dcnv_cfg(1) == 1 .or. &
167 prc_cldmp_cfg(1) == 1 ) then
168 else
169 return
170 endif
171
172 ! Check that input data file exists.
173 inquire (file = trim(suite_sim_file), exist = exists)
174 if (.not. exists) then
175 errmsg = 'CCPP suite simulator file: '//trim(suite_sim_file)//' does not exist'
176 errflg = 1
177 return
178 endif
179
180 !
181 ! Read data file...
182 !
183
184 ! Open file
185 status = nf90_open(trim(suite_sim_file), nf90_nowrite, ncid)
186 if (status /= nf90_noerr) then
187 errmsg = 'Error reading in CCPP suite simulator file: '//trim(suite_sim_file)
188 errflg = 1
189 return
190 endif
191
192 ! Metadata (dimensions)
193 status = nf90_inq_dimid(ncid, 'time', dimid)
194 if (status == nf90_noerr) then
195 status = nf90_inquire_dimension(ncid, dimid, len = ntime_data)
196 else
197 errmsg = 'CCPP suite simulator file: '//trim(suite_sim_file)//' does not contain [time] dimension'
198 errflg = 1
199 return
200 endif
201
202 status = nf90_inq_dimid(ncid, 'lev', dimid)
203 if (status == nf90_noerr) then
204 status = nf90_inquire_dimension(ncid, dimid, len = nlev_data)
205 else
206 errmsg = 'CCPP suite simulator file: '//trim(suite_sim_file)//' does not contain [lev] dimension'
207 errflg = 1
208 return
209 endif
210
211 ! Allocate space and read in data
212 allocate(physics_process(nprc_sim))
213 physics_process(1)%active_name = ''
214 physics_process(1)%iactive_scheme = 0
215 physics_process(1)%active_tsp = .false.
216 do iprc = 1,nprc_sim
217 allocate(physics_process(iprc)%tend1d%T( nlev_data ))
218 allocate(physics_process(iprc)%tend1d%u( nlev_data ))
219 allocate(physics_process(iprc)%tend1d%v( nlev_data ))
220 allocate(physics_process(iprc)%tend1d%q( nlev_data ))
221 allocate(physics_process(iprc)%tend2d%time( ntime_data))
222 allocate(physics_process(iprc)%tend2d%T( nlev_data, ntime_data))
223 allocate(physics_process(iprc)%tend2d%u( nlev_data, ntime_data))
224 allocate(physics_process(iprc)%tend2d%v( nlev_data, ntime_data))
225 allocate(physics_process(iprc)%tend2d%q( nlev_data, ntime_data))
226
227 ! Temporal info
228 status = nf90_inq_varid(ncid, 'times', varid)
229 if (status == nf90_noerr) then
230 status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%time)
231 else
232 errmsg = 'SCM data tendency file: '//trim(suite_sim_file)//' does not contain times variable'
233 errflg = 1
234 return
235 endif
236
237 if (iprc == prc_swrad_cfg(3)) then
238 ! Metadata
239 physics_process(iprc)%order = iprc
240 physics_process(iprc)%name = "SWRAD"
241 if (prc_swrad_cfg(1) == 1) then
242 physics_process(iprc)%use_sim = .true.
243 else
244 physics_process(1)%nprg_active = 1
245 iactive_t = 1
246 endif
247 if (prc_swrad_cfg(2) == 1) then
248 physics_process(iprc)%time_split = .true.
249 endif
250
251 ! Data
252 status = nf90_inq_varid(ncid, 'dT_dt_swrad', varid)
253 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%T)
254 endif
255
256 if (iprc == prc_lwrad_cfg(3)) then
257 ! Metadata
258 physics_process(iprc)%order = iprc
259 physics_process(iprc)%name = "LWRAD"
260 if (prc_lwrad_cfg(1) == 1) then
261 physics_process(iprc)%use_sim = .true.
262 else
263 physics_process(1)%nprg_active = 1
264 iactive_t = 1
265 endif
266 if (prc_lwrad_cfg(2) == 1) then
267 physics_process(iprc)%time_split = .true.
268 endif
269
270 ! Data
271 status = nf90_inq_varid(ncid, 'dT_dt_lwrad', varid)
272 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%T)
273 endif
274
275 if (iprc == prc_gwd_cfg(3)) then
276 ! Metadata
277 physics_process(iprc)%order = iprc
278 physics_process(iprc)%name = "GWD"
279 if (prc_gwd_cfg(1) == 1) then
280 physics_process(iprc)%use_sim = .true.
281 else
282 physics_process(1)%nprg_active = 3
283 iactive_t = 1
284 iactive_u = 2
285 iactive_v = 3
286 endif
287 if (prc_gwd_cfg(2) == 1) then
288 physics_process(iprc)%time_split = .true.
289 endif
290
291 ! Data
292 status = nf90_inq_varid(ncid, 'dT_dt_cgwd', varid)
293 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%T)
294 status = nf90_inq_varid(ncid, 'du_dt_cgwd', varid)
295 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%u)
296 status = nf90_inq_varid(ncid, 'dv_dt_cgwd', varid)
297 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%v)
298 endif
299
300 if (iprc == prc_pbl_cfg(3)) then
301 ! Metadata
302 physics_process(iprc)%order = iprc
303 physics_process(iprc)%name = "PBL"
304 if (prc_pbl_cfg(1) == 1) then
305 physics_process(iprc)%use_sim = .true.
306 else
307 physics_process(1)%nprg_active = 4
308 iactive_t = 1
309 iactive_u = 2
310 iactive_v = 3
311 iactive_q = 4
312 endif
313 if (prc_pbl_cfg(2) == 1) then
314 physics_process(iprc)%time_split = .true.
315 endif
316
317 ! Data
318 status = nf90_inq_varid(ncid, 'dT_dt_pbl', varid)
319 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%T)
320 status = nf90_inq_varid(ncid, 'dq_dt_pbl', varid)
321 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%q)
322 status = nf90_inq_varid(ncid, 'du_dt_pbl', varid)
323 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%u)
324 status = nf90_inq_varid(ncid, 'dv_dt_pbl', varid)
325 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%v)
326 endif
327
328 if (iprc == prc_scnv_cfg(3)) then
329 ! Metadata
330 physics_process(iprc)%order = iprc
331 physics_process(iprc)%name = "SCNV"
332 if (prc_scnv_cfg(1) == 1) then
333 physics_process(iprc)%use_sim = .true.
334 else
335 physics_process(1)%nprg_active = 4
336 iactive_t = 1
337 iactive_u = 2
338 iactive_v = 3
339 iactive_q = 4
340 endif
341 if (prc_scnv_cfg(2) == 1) then
342 physics_process(iprc)%time_split = .true.
343 endif
344
345 ! Data
346 status = nf90_inq_varid(ncid, 'dT_dt_shalconv', varid)
347 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%T)
348 status = nf90_inq_varid(ncid, 'du_dt_shalconv', varid)
349 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%u)
350 status = nf90_inq_varid(ncid, 'dv_dt_shalconv', varid)
351 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%v)
352 status = nf90_inq_varid(ncid, 'dq_dt_shalconv', varid)
353 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%q)
354 endif
355
356 if (iprc == prc_dcnv_cfg(3)) then
357 ! Metadata
358 physics_process(iprc)%order = iprc
359 physics_process(iprc)%name = "DCNV"
360 if (prc_dcnv_cfg(1) == 1) then
361 physics_process(iprc)%use_sim = .true.
362 else
363 physics_process(1)%nprg_active = 4
364 iactive_t = 1
365 iactive_u = 2
366 iactive_v = 3
367 iactive_q = 4
368 endif
369 if (prc_dcnv_cfg(2) == 1) then
370 physics_process(iprc)%time_split = .true.
371 endif
372 ! Data
373 status = nf90_inq_varid(ncid, 'dT_dt_deepconv', varid)
374 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%T)
375 status = nf90_inq_varid(ncid, 'du_dt_deepconv', varid)
376 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%u)
377 status = nf90_inq_varid(ncid, 'dv_dt_deepconv', varid)
378 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%v)
379 status = nf90_inq_varid(ncid, 'dq_dt_deepconv', varid)
380 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%q)
381 endif
382
383 if (iprc == prc_cldmp_cfg(3)) then
384 ! Metadata
385 physics_process(iprc)%order = iprc
386 physics_process(iprc)%name = "cldMP"
387 if (prc_cldmp_cfg(1) == 1) then
388 physics_process(iprc)%use_sim = .true.
389 else
390 physics_process(1)%nprg_active = 2
391 iactive_t = 1
392 iactive_q = 2
393 endif
394 if (prc_cldmp_cfg(2) == 1) then
395 physics_process(iprc)%time_split = .true.
396 endif
397
398 ! Data
399 status = nf90_inq_varid(ncid, 'dT_dt_micro', varid)
400 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%T)
401 status = nf90_inq_varid(ncid, 'dq_dt_micro', varid)
402 if (status == nf90_noerr) status = nf90_get_var( ncid, varid, physics_process(iprc)%tend2d%q)
403 endif
404
405 ! Which process-suite is "active"? Is process time-split?
406 if (.not. physics_process(iprc)%use_sim) then
407 physics_process(1)%iactive_scheme = iprc
408 physics_process(1)%active_name = physics_process(iprc)%name
409 if (physics_process(iprc)%time_split) then
410 physics_process(1)%active_tsp = .true.
411 endif
412 endif
413
414 enddo
415
416 if (physics_process(1)%iactive_scheme == 0) then
417 errflg = 1
418 errmsg = "ERROR: No active suite set for CCPP suite simulator"
419 return
420 endif
421
422 print*, "-----------------------------------"
423 print*, "--- Using CCPP suite simulator ---"
424 print*, "-----------------------------------"
425 do iprc = 1,nprc_sim
426 if (physics_process(iprc)%use_sim) then
427 print*," simulate_suite: ", trim(physics_process(iprc)%name)
428 print*," order: ", physics_process(iprc)%order
429 print*," time_split: ", physics_process(iprc)%time_split
430 else
431 print*, " active_suite: ", trim(physics_process(1)%active_name)
432 print*, " order: ", physics_process(physics_process(1)%iactive_scheme)%order
433 print*, " time_split : ", physics_process(1)%active_tsp
434 endif
435 enddo
436 print*, "-----------------------------------"
437 print*, "-----------------------------------"
438
439 end subroutine load_ccpp_suite_sim
440
441end module gfs_ccpp_suite_sim_pre
This type contains the meta information and data for each physics process.