ICON Community Interface 0.4.0
Loading...
Searching...
No Matches
calc_water_column_plugin.F90
Go to the documentation of this file.
1! --------------------------------------------------------------------
3!
4! This example illustrates a simple diagnostic application, i.e.
5! the calculations of the liquid water path(lwp), the ice water path (iwp)
6! and the total water column (twc).
7! For this
8! - the variables lwp, iwp and twc need to be added to the ICON variables
9! including the definition of metadata, as the units
10! - the humidity tracers from ICON need to be accessed and it needs to be
11! checked, if qg exists or not, which depends on the microphysics scheme
12! chosen in ICON.
13! - for the calculation itself, the ICON prognostic variable rho is
14! required and
15! - the descriptive data (p_patch%cell%hhl) is used
16!
17! The calculation is performed for each active domain.
18!
19! @authors 11/2023 :: ICON Community Interface <comin@icon-model.org>
20!
21! SPDX-License-Identifier: BSD-3-Clause
22!
23! Please see the file LICENSE in the root of the source tree for this code.
24! Where software is supplied by third parties, it is indicated in the
25! headers of the routines.
26! --------------------------------------------------------------------
27MODULE calc_water_column_plugin
28
29 USE comin_plugin_interface, ONLY : comin_callback_register, &
30 & comin_var_get, &
31 & t_comin_var_descriptor, t_comin_var_handle, &
32 & comin_var_request_add, &
33 & comin_descrdata_get_domain, t_comin_descrdata_domain, &
34 & comin_descrdata_get_global, t_comin_descrdata_global, &
35 & t_comin_setup_version_info, comin_setup_get_version, &
36 & ep_secondary_constructor, ep_destructor, ep_atm_physics_before, &
37 & comin_flag_read, comin_flag_write, comin_zaxis_2d, &
38 & comin_parallel_get_host_mpi_rank, comin_current_get_domain_id, &
39 & t_comin_plugin_info, comin_current_get_plugin_info, &
40 & comin_plugin_finish, comin_metadata_set, &
41 & comin_metadata_get, &
42 & comin_error_set_errors_return, comin_success, comin_error_get, &
43 & comin_print_info, &
44 & comin_dim_semantics_nproma, comin_dim_semantics_level,&
45 & comin_dim_semantics_block, comin_dim_semantics_unused
46
47 IMPLICIT NONE
48 PRIVATE
49
50 CHARACTER(LEN=*), PARAMETER :: pluginname = "calc_water_column_plugin"
51
53 INTEGER, PARAMETER :: wp = selected_real_kind(12,307)
54 TYPE(t_comin_setup_version_info) :: version
55
56 TYPE :: t_plugin_vars
57 ! ICON data
58 TYPE(t_comin_var_handle) :: rho
59 TYPE(t_comin_var_handle) :: qv
60 TYPE(t_comin_var_handle) :: qc
61 TYPE(t_comin_var_handle) :: qi
62 TYPE(t_comin_var_handle) :: qr
63 TYPE(t_comin_var_handle) :: qs
64 TYPE(t_comin_var_handle) :: qg
65 ! variables added to ICON by this plugin
66 TYPE(t_comin_var_handle) :: lwp ! liquid water path (qc+qr)
67 TYPE(t_comin_var_handle) :: iwp ! ice water path (qi+qs+qg)
68 TYPE(t_comin_var_handle) :: twc ! (full) water path (qv+qc+qr+qi+qs+qg)
69 END type t_plugin_vars
70 ! pugin variables per patch
71 TYPE(t_plugin_vars), DIMENSION(:), ALLOCATABLE :: pvpp
72
74 TYPE(t_comin_descrdata_domain), POINTER :: p_patch
75 ! access global setup information
76 TYPE(t_comin_descrdata_global), POINTER :: p_global
77 !
78 LOGICAL :: lqg = .true. ! qg is not always available
79
80 CHARACTER(LEN=120) :: text
81
82 PUBLIC :: comin_main
83 PUBLIC :: calc_water_column_constructor
84 PUBLIC :: calc_water_column_diagfct
85 PUBLIC :: calc_water_column_destructor
86
87CONTAINS
88
89 ! --------------------------------------------------------------------
90 ! ComIn primary constructor.
91 ! --------------------------------------------------------------------
92 SUBROUTINE comin_main() BIND(C)
93 !
94 IMPLICIT NONE
95 !
96 TYPE(t_comin_plugin_info) :: this_plugin
97 TYPE(t_comin_var_descriptor) :: lwp_d, iwp_d, twc_d
98 INTEGER :: jg
99
100 CALL comin_print_info('- setup')
101
102 version = comin_setup_get_version()
103 IF (version%version_no_major > 1) THEN
104 CALL comin_plugin_finish('comin_main ('//pluginname//')', "incompatible version!")
105 END IF
106
107 ! get descriptive data structures
108 p_global => comin_descrdata_get_global()
109
110 WRITE (text,'(a,i4)') '- number of domains: ', p_global%n_dom
111 CALL comin_print_info(text)
112
114 CALL comin_current_get_plugin_info(this_plugin)
115 WRITE (text,'(a,a,a,i4)') "- plugin ", this_plugin%name, " has id: ", this_plugin%id
116 CALL comin_print_info(text)
117
119
120 ! request host model to add local variables
121
122 DO jg = 1, p_global%n_dom
123 ! liquid water path (lwp) for first domain
124 lwp_d = t_comin_var_descriptor( id = jg, name = "lwp" )
125 CALL comin_var_request_add_wrapper(lwp_d, lmode_exclusive=.false. &
126 , zaxis_id = comin_zaxis_2d, tracer =.false., restart=.false., units='kg/m2')
127
128 ! ice water path (iwp) for first domain
129 iwp_d = t_comin_var_descriptor( id = jg, name = "iwp" )
130 CALL comin_var_request_add_wrapper(iwp_d, lmode_exclusive=.false. &
131 , zaxis_id = comin_zaxis_2d, tracer =.false., restart=.false., units='kg/m2')
132
133 ! total water column (twc) for first domain
134 twc_d = t_comin_var_descriptor( id = jg, name = "twc" )
135 CALL comin_var_request_add_wrapper(twc_d, lmode_exclusive=.false. &
136 , zaxis_id = comin_zaxis_2d, tracer =.false., restart=.false., units='kg/m2')
137 END DO
138
139 ! register callbacks
140 CALL comin_callback_register(ep_atm_physics_before, calc_water_column_diagfct)
141 CALL comin_callback_register(ep_secondary_constructor, calc_water_column_constructor)
142 CALL comin_callback_register(ep_destructor, calc_water_column_destructor)
143
144 END SUBROUTINE comin_main
145
146 ! --------------------------------------------------------------------
147 ! ComIn secondary constructor.
148 ! --------------------------------------------------------------------
149 SUBROUTINE calc_water_column_constructor() BIND(C)
150
151 IMPLICIT NONE
152
153 TYPE(t_comin_var_descriptor) :: var_desc
154 INTEGER :: jg
155
156 ALLOCATE(pvpp(p_global%n_dom))
157
158 domain_loop: DO jg = 1, p_global%n_dom
159 WRITE (text,'(a,i4)') ' - get required meteorological data. dom: ', jg
160 CALL comin_print_info(text)
161
162 ! 1 - GET ICON VARIABLES:
163 ! 1.1 - A RHO
164 var_desc = t_comin_var_descriptor('rho', jg)
165 CALL comin_var_get([ep_atm_physics_before], &
166 & var_desc, comin_flag_read, pvpp(jg)%rho)
167 CALL check_variable(pvpp(jg)%rho, 'rho')
168
169 WRITE (text,'(a,i4)') ' - get humidity tracer. dom: ',jg
170 CALL comin_print_info(text)
171 ! 1.2 qv
172 var_desc = t_comin_var_descriptor('qv', jg)
173 CALL comin_var_get([ep_atm_physics_before], &
174 & var_desc, comin_flag_read, pvpp(jg)%qv)
175 CALL check_variable(pvpp(jg)%qv, 'qv')
176 ! 1.3 qc
177 var_desc = t_comin_var_descriptor('qc', jg)
178 CALL comin_var_get([ep_atm_physics_before], &
179 & var_desc, comin_flag_read, pvpp(jg)%qc)
180 CALL check_variable(pvpp(jg)%qc,'qc')
181 ! 1.4 qi
182 var_desc = t_comin_var_descriptor('qi', jg)
183 CALL comin_var_get([ep_atm_physics_before], &
184 & var_desc, comin_flag_read, pvpp(jg)%qi)
185 CALL check_variable(pvpp(jg)%qi, 'qi')
186 ! 1.5 qr
187 var_desc = t_comin_var_descriptor('qr', jg)
188 CALL comin_var_get([ep_atm_physics_before], &
189 & var_desc, comin_flag_read, pvpp(jg)%qr)
190 CALL check_variable(pvpp(jg)%qr, 'qr')
191 ! 1.6 qs
192 var_desc = t_comin_var_descriptor('qs', jg)
193 CALL comin_var_get([ep_atm_physics_before], &
194 & var_desc, comin_flag_read, pvpp(jg)%qs)
195 CALL check_variable(pvpp(jg)%qs, 'qs')
196 ! 1.7 qg
197 var_desc = t_comin_var_descriptor('qg', jg)
198 CALL comin_error_set_errors_return(.true.)
199 CALL comin_var_get([ep_atm_physics_before], &
200 & var_desc, comin_flag_read, pvpp(jg)%qg)
201 lqg = (comin_error_get() == comin_success)
202 CALL comin_error_set_errors_return(.false.)
203 IF (lqg) CALL check_variable(pvpp(jg)%qg, 'qg')
204
205 WRITE (text,'(a,i4)') ' - get plugin variables - requested to be created by ICON. dom: ', jg
206 CALL comin_print_info(text)
207
208 ! 2 - GET OWN VARIABLE requested in comin_main
209 CALL comin_var_get([ep_atm_physics_before], &
210 & t_comin_var_descriptor(name='lwp', id=jg),&
211 & comin_flag_write, pvpp(jg)%lwp)
212 CALL check_variable(pvpp(jg)%lwp, 'lwp',dim=2)
213
214 CALL comin_var_get([ep_atm_physics_before], &
215 & t_comin_var_descriptor(name='iwp', id=jg), &
216 & comin_flag_write, pvpp(jg)%iwp)
217 CALL check_variable(pvpp(jg)%iwp, 'iwp',dim=2)
218
219 CALL comin_var_get([ep_atm_physics_before], &
220 & t_comin_var_descriptor(name='twc', id=jg), &
221 & comin_flag_write, pvpp(jg)%twc)
222 CALL check_variable(pvpp(jg)%twc, 'twc',dim=2)
223
224 END DO domain_loop
225
226 END SUBROUTINE calc_water_column_constructor
227
228 ! --------------------------------------------------------------------
229 ! ComIn callback function.
230 ! --------------------------------------------------------------------
231 SUBROUTINE calc_water_column_diagfct() BIND(C)
232
233 IMPLICIT NONE
234
235 TYPE(t_comin_var_descriptor) :: lwp_d, iwp_d, twc_d
236 INTEGER :: jk, jg
237
238 REAL(wp), POINTER, DIMENSION(:,:,:) :: rho_3d => null()
239
240 REAL(wp), POINTER, DIMENSION(:,:,:) :: lwp_3d => null()
241 REAL(wp), POINTER, DIMENSION(:,:,:) :: iwp_3d => null()
242 REAL(wp), POINTER, DIMENSION(:,:,:) :: twc_3d => null()
243
244 REAL(wp), POINTER, DIMENSION(:,:,:) :: qv_3d => null()
245 REAL(wp), POINTER, DIMENSION(:,:,:) :: qc_3d => null()
246 REAL(wp), POINTER, DIMENSION(:,:,:) :: qi_3d => null()
247 REAL(wp), POINTER, DIMENSION(:,:,:) :: qr_3d => null()
248 REAL(wp), POINTER, DIMENSION(:,:,:) :: qs_3d => null()
249 REAL(wp), POINTER, DIMENSION(:,:,:) :: qg_3d => null()
250 ! conversion factor kg / kg => kg / m2
251 REAL(wp), POINTER, DIMENSION(:,:,:) :: conv => null()
252
253 CHARACTER(LEN=:), ALLOCATABLE :: units
254 CHARACTER(LEN=200) :: text = ''
255
256 ! get current domain:
257 jg = comin_current_get_domain_id()
258 ! get current patch description data
259 p_patch => comin_descrdata_get_domain(jg)
260
261 WRITE (text,'(a,i4)') '- calculate water columns before physics. dom: ', jg
262 CALL comin_print_info(text)
263
264 CALL pvpp(jg)%qv%to_3d(qv_3d)
265 CALL pvpp(jg)%qi%to_3d(qi_3d)
266 CALL pvpp(jg)%qr%to_3d(qr_3d)
267 CALL pvpp(jg)%qs%to_3d(qs_3d)
268 CALL pvpp(jg)%qc%to_3d(qc_3d)
269 IF (lqg) CALL pvpp(jg)%qg%to_3d(qg_3d)
270
271 CALL pvpp(jg)%rho%to_3d(rho_3d)
272
273 CALL pvpp(jg)%lwp%to_3d(lwp_3d)
274 CALL pvpp(jg)%iwp%to_3d(iwp_3d)
275 CALL pvpp(jg)%twc%to_3d(twc_3d)
276
277 ! conversion factor 1/kg to 1/m2
278 ALLOCATE(conv(SIZE(qv_3d,1),SIZE(qv_3d,2),SIZE(qv_3d,3)))
279 conv(:,:,:) = 0._wp
280 DO jk=1,SIZE(qv_3d,2)
281 ! convert 1/kg to 1/m2
282 conv(:,jk,:) = rho_3d(:,jk,:) * &
283 (p_patch%cells%hhl(:,jk,:) - p_patch%cells%hhl(:,jk+1,:))
284 END DO
285
286 ! calculate liquid water part / ice_water_path and total water column
287 lwp_3d(:,:,:) = 0._wp
288 iwp_3d(:,:,:) = 0._wp
289 twc_3d(:,:,:) = 0._wp
290 DO jk=1,SIZE(qv_3d,2)
291 lwp_3d(:,:,1) = lwp_3d(:,:,1) &
292 + (qr_3d(:,jk,:) + qc_3d(:,jk,:)) * conv(:,jk,:)
293 iwp_3d(:,:,1) = iwp_3d(:,:,1) &
294 + (qi_3d(:,jk,:)+qs_3d(:,jk,:)) * conv(:,jk,:)
295 END DO
296 IF (lqg) THEN
297 DO jk=1,SIZE(qv_3d,2)
298 iwp_3d(:,:,1) = iwp_3d(:,:,1) + qg_3d(:,jk,:) * conv(:,jk,:)
299 END DO
300 END IF
301 twc_3d(:,:,1) = iwp_3d(:,:,1) + lwp_3d(:,:,1)
302 DO jk=1,SIZE(qv_3d,2)
303 twc_3d(:,:,1) = twc_3d(:,:,1) + qv_3d(:,jk,:) * conv(:,jk,:)
304 END DO
305
306 WRITE (text,'(a,i4)') '- results of plugin calc_water_column. dom: ', jg
307
308 lwp_d = t_comin_var_descriptor( id = jg, name = "lwp" )
309 CALL comin_metadata_get(lwp_d, "units", units)
310 write (text,fmt='(A,A5,A,F8.4)') '- maximum liquid water path (',trim(units),'): ', maxval(lwp_3d(:,:,1))
311 CALL message(text, idom=jg,lall=.true.)
312 iwp_d = t_comin_var_descriptor( id = jg, name = "iwp" )
313 CALL comin_metadata_get(iwp_d, "units", units)
314 write (text,fmt='(A,A5,A,F8.4)') '- maximum ice water path (',trim(units),'): ', maxval(iwp_3d(:,:,1))
315 CALL message(text, idom=jg,lall=.true.)
316 twc_d = t_comin_var_descriptor( id = jg, name = "twc" )
317 write (text,fmt='(A,A5,A,F8.4)') '- maximum total water column (',trim(units),'): ', maxval(twc_3d(:,:,1))
318 CALL message(text, idom=jg,lall=.true.)
319
320 ! CLEANUP
321 DEALLOCATE(conv); NULLIFY(conv)
322 NULLIFY(lwp_3d, iwp_3d, twc_3d)
323 NULLIFY(qv_3d,qc_3d,qi_3d,qr_3d,qs_3d,qg_3d, rho_3d)
324 NULLIFY(p_patch)
325
326 END SUBROUTINE calc_water_column_diagfct
327
328 ! --------------------------------------------------------------------
329 ! ComIn callback function.
330 ! --------------------------------------------------------------------
331
332 SUBROUTINE calc_water_column_destructor() BIND(C)
333
334 IMPLICIT NONE
335
336 CALL comin_print_info(' - destructor.')
337
338 END SUBROUTINE calc_water_column_destructor
339
340 !---------------------------------------------------------------------
341 !---------------------------------------------------------------------
342 !---------------------------------------------------------------------
343 ! PRIVATE ROUTINES
344 !---------------------------------------------------------------------
345 !---------------------------------------------------------------------
346 !---------------------------------------------------------------------
347
348 SUBROUTINE comin_var_request_add_wrapper(descriptor, lmode_exclusive, zaxis_id &
349 , tracer, restart, units)
350
351 IMPLICIT NONE
352
353 ! I/O
354 TYPE(t_comin_var_descriptor), INTENT(IN) :: descriptor
355 LOGICAL, OPTIONAL, INTENT(IN) :: lmode_exclusive
356 INTEGER, OPTIONAL, INTENT(IN) :: zaxis_id
357 LOGICAL, OPTIONAL, INTENT(IN) :: tracer
358 LOGICAL, OPTIONAL, INTENT(IN) :: restart
359 CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: units
360 ! LOCAL
361 LOGICAL :: lexclusive
362
363 IF (PRESENT(lmode_exclusive)) THEN
364 lexclusive = lmode_exclusive
365 ELSE
366 lexclusive = .false.
367 END IF
368
369 CALL comin_var_request_add(descriptor, lexclusive)
370
371 IF (PRESENT(zaxis_id)) THEN
372 CALL comin_metadata_set(descriptor, "zaxis_id", zaxis_id)
373 END IF
374 IF (PRESENT(tracer)) THEN
375 CALL comin_metadata_set(descriptor, "tracer", tracer)
376 END IF
377 IF (PRESENT(restart)) THEN
378 CALL comin_metadata_set(descriptor, "restart", restart)
379 END IF
380 IF (PRESENT(units)) THEN
381 CALL comin_metadata_set(descriptor, "units", trim(units))
382 END IF
383
384 END SUBROUTINE comin_var_request_add_wrapper
385
386 ! --------------------------------------------------------------------------------------------------
387
388 SUBROUTINE check_variable(var, name, dim)
389
390 IMPLICIT NONE
391
392 ! I/O
393 TYPE(t_comin_var_handle), INTENT(IN) :: var
394 CHARACTER(LEN=*), INTENT(IN) :: name
395 INTEGER, INTENT(IN), OPTIONAL :: dim
396 ! LOCAL
397 INTEGER :: idim
398
399 IF (PRESENT(dim)) THEN
400 idim = dim
401 ELSE
402 idim = 3
403 END IF
404 SELECT CASE (idim)
405 CASE (2)
406 CALL check_dimensions_2d(var, name)
407 CASE (3)
408 CALL check_dimensions_3d(var, name)
409 CASE DEFAULT
410 CALL comin_plugin_finish (pluginname//' unsupported dimension', 'check variable')
411 END SELECT
412
413 END SUBROUTINE check_variable
414
415 ! --------------------------------------------------------------------------------------------------
416
417 SUBROUTINE check_dimensions_3d(var, name)
418
419 IMPLICIT NONE
420
421 TYPE(t_comin_var_handle), INTENT(IN) :: var
422 CHARACTER(LEN=*), INTENT(IN) :: name
423
424 IF (any(var%dim_semantics() /= [comin_dim_semantics_nproma, &
425 & comin_dim_semantics_level, &
426 & comin_dim_semantics_block, &
427 & comin_dim_semantics_unused, &
428 & comin_dim_semantics_unused] )) &
429 & CALL comin_plugin_finish(pluginname//': '//trim(name) &
430 & , "Dimension check failed!")
431
432 END SUBROUTINE check_dimensions_3d
433
434 ! --------------------------------------------------------------------------------------------------
435
436 ! --------------------------------------------------------------------------------------------------
437
438 SUBROUTINE check_dimensions_2d(var, name)
439
440 IMPLICIT NONE
441
442 TYPE(t_comin_var_handle), INTENT(IN) :: var
443 CHARACTER(LEN=*), INTENT(IN) :: name
444
445 IF (any(var%dim_semantics() /= [comin_dim_semantics_nproma, &
446 & comin_dim_semantics_block, &
447 & comin_dim_semantics_unused, &
448 & comin_dim_semantics_unused, &
449 & comin_dim_semantics_unused] )) &
450 & CALL comin_plugin_finish(pluginname//': '//trim(name) &
451 , "Dimension check failed!")
452
453 END SUBROUTINE check_dimensions_2d
454
455 ! --------------------------------------------------------------------------------------------------
456
457 SUBROUTINE message(text, idom, lall)
458
459 IMPLICIT NONE
460
461 CHARACTER(LEN=*) :: text
462 INTEGER, OPTIONAL :: idom
463 LOGICAL, OPTIONAL :: lall
464
465 ! LOCAL
466 LOGICAL :: la
467
468 IF (PRESENT(lall)) THEN
469 la = lall
470 ELSE
471 la = .false.
472 END IF
473
474 IF (comin_parallel_get_host_mpi_rank() == 0 .OR. la) THEN
475 IF (PRESENT(idom)) THEN
476 WRITE (0,fmt='(A,I2,A)') pluginname//' domain: ',idom,' '//trim(text)
477 ELSE
478 WRITE (0,fmt='(A,A)') pluginname//' ',trim(text)
479 END IF
480 END IF
481
482 END SUBROUTINE message
483
484END MODULE calc_water_column_plugin
485! --------------------------------------------------------------------------------------------------
486! --------------------------------------------------------------------------------------------------
calc_water_column_plugin
Example plugin for the ICON Community Interface (ComIn)
Definition calc_water_column_plugin.F90:27
calc_water_column_plugin::pvpp
type(t_plugin_vars), dimension(:), allocatable pvpp
Definition calc_water_column_plugin.F90:71
calc_water_column_plugin::comin_main
subroutine, public comin_main()
Definition calc_water_column_plugin.F90:93
calc_water_column_plugin::calc_water_column_constructor
subroutine, public calc_water_column_constructor()
Definition calc_water_column_plugin.F90:150
calc_water_column_plugin::calc_water_column_destructor
subroutine, public calc_water_column_destructor()
Definition calc_water_column_plugin.F90:333
calc_water_column_plugin::calc_water_column_diagfct
subroutine, public calc_water_column_diagfct()
Definition calc_water_column_plugin.F90:232
comin_plugin_interface
Definition comin_plugin_interface.F90:25
Generated by doxygen 1.11.0