ComIn 0.5.1
ICON Community Interface
Loading...
Searching...
No Matches
calc_water_column_plugin.F90
Go to the documentation of this file.
1! --------------------------------------------------------------------
2!> Example plugin for the ICON Community Interface (ComIn)
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! --------------------------------------------------------------------
28
30 & comin_var_get, &
31 & t_comin_var_descriptor, t_comin_var_handle, &
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, &
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
52 !> working precision (will be compared to ComIn's and ICON's)
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
73 !> access descriptive data structures
74 TYPE(t_comin_descrdata_domain) :: p_patch
75 ! access global setup information
76 TYPE(t_comin_descrdata_global) :: p_global
77 !
78 LOGICAL :: lqg = .true. ! qg is not always available
79
80 CHARACTER(LEN=120) :: text
81
82 PUBLIC :: comin_main
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%get_n_dom()
111 CALL comin_print_info(text)
112
113 !> check plugin id
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
118 !> add requests for additional ICON variables
119
120 ! request host model to add local variables
121
122 DO jg = 1, p_global%get_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)
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%get_n_dom()))
157
158 domain_loop: DO jg = 1, p_global%get_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)
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)
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(:,:,:) :: hhl
239
240 REAL(wp), POINTER, DIMENSION(:,:,:) :: rho_3d => null()
241
242 REAL(wp), POINTER, DIMENSION(:,:,:) :: lwp_3d => null()
243 REAL(wp), POINTER, DIMENSION(:,:,:) :: iwp_3d => null()
244 REAL(wp), POINTER, DIMENSION(:,:,:) :: twc_3d => null()
245
246 REAL(wp), POINTER, DIMENSION(:,:,:) :: qv_3d => null()
247 REAL(wp), POINTER, DIMENSION(:,:,:) :: qc_3d => null()
248 REAL(wp), POINTER, DIMENSION(:,:,:) :: qi_3d => null()
249 REAL(wp), POINTER, DIMENSION(:,:,:) :: qr_3d => null()
250 REAL(wp), POINTER, DIMENSION(:,:,:) :: qs_3d => null()
251 REAL(wp), POINTER, DIMENSION(:,:,:) :: qg_3d => null()
252 ! conversion factor kg / kg => kg / m2
253 REAL(wp), POINTER, DIMENSION(:,:,:) :: conv => null()
254
255 CHARACTER(LEN=:), ALLOCATABLE :: units
256 CHARACTER(LEN=200) :: text = ''
257
258 ! get current domain:
260 ! get current patch description data
261 p_patch = comin_descrdata_get_domain(jg)
262
263 WRITE (text,'(a,i4)') '- calculate water columns before physics. dom: ', jg
264 CALL comin_print_info(text)
265
266 CALL pvpp(jg)%qv%to_3d(qv_3d)
267 CALL pvpp(jg)%qi%to_3d(qi_3d)
268 CALL pvpp(jg)%qr%to_3d(qr_3d)
269 CALL pvpp(jg)%qs%to_3d(qs_3d)
270 CALL pvpp(jg)%qc%to_3d(qc_3d)
271 IF (lqg) CALL pvpp(jg)%qg%to_3d(qg_3d)
272
273 CALL pvpp(jg)%rho%to_3d(rho_3d)
274
275 CALL pvpp(jg)%lwp%to_3d(lwp_3d)
276 CALL pvpp(jg)%iwp%to_3d(iwp_3d)
277 CALL pvpp(jg)%twc%to_3d(twc_3d)
278
279 ! conversion factor 1/kg to 1/m2
280 ALLOCATE(conv(SIZE(qv_3d,1),SIZE(qv_3d,2),SIZE(qv_3d,3)))
281 conv(:,:,:) = 0._wp
282 hhl => p_patch%cells%get_hhl()
283 DO jk=1,SIZE(qv_3d,2)
284 ! convert 1/kg to 1/m2
285 conv(:,jk,:) = rho_3d(:,jk,:) * &
286 (hhl(:,jk,:) - hhl(:,jk+1,:))
287 END DO
288
289 ! calculate liquid water part / ice_water_path and total water column
290 lwp_3d(:,:,:) = 0._wp
291 iwp_3d(:,:,:) = 0._wp
292 twc_3d(:,:,:) = 0._wp
293 DO jk=1,SIZE(qv_3d,2)
294 lwp_3d(:,:,1) = lwp_3d(:,:,1) &
295 + (qr_3d(:,jk,:) + qc_3d(:,jk,:)) * conv(:,jk,:)
296 iwp_3d(:,:,1) = iwp_3d(:,:,1) &
297 + (qi_3d(:,jk,:)+qs_3d(:,jk,:)) * conv(:,jk,:)
298 END DO
299 IF (lqg) THEN
300 DO jk=1,SIZE(qv_3d,2)
301 iwp_3d(:,:,1) = iwp_3d(:,:,1) + qg_3d(:,jk,:) * conv(:,jk,:)
302 END DO
303 END IF
304 twc_3d(:,:,1) = iwp_3d(:,:,1) + lwp_3d(:,:,1)
305 DO jk=1,SIZE(qv_3d,2)
306 twc_3d(:,:,1) = twc_3d(:,:,1) + qv_3d(:,jk,:) * conv(:,jk,:)
307 END DO
308
309 WRITE (text,'(a,i4)') '- results of plugin calc_water_column. dom: ', jg
310
311 lwp_d = t_comin_var_descriptor( id = jg, name = "lwp" )
312 CALL comin_metadata_get(lwp_d, "units", units)
313 write (text,fmt='(A,A5,A,F8.4)') '- maximum liquid water path (',trim(units),'): ', maxval(lwp_3d(:,:,1))
314 CALL message(text, idom=jg,lall=.true.)
315 iwp_d = t_comin_var_descriptor( id = jg, name = "iwp" )
316 CALL comin_metadata_get(iwp_d, "units", units)
317 write (text,fmt='(A,A5,A,F8.4)') '- maximum ice water path (',trim(units),'): ', maxval(iwp_3d(:,:,1))
318 CALL message(text, idom=jg,lall=.true.)
319 twc_d = t_comin_var_descriptor( id = jg, name = "twc" )
320 write (text,fmt='(A,A5,A,F8.4)') '- maximum total water column (',trim(units),'): ', maxval(twc_3d(:,:,1))
321 CALL message(text, idom=jg,lall=.true.)
322
323 ! CLEANUP
324 DEALLOCATE(conv); NULLIFY(conv)
325 NULLIFY(lwp_3d, iwp_3d, twc_3d)
326 NULLIFY(qv_3d,qc_3d,qi_3d,qr_3d,qs_3d,qg_3d, rho_3d)
327
328 END SUBROUTINE calc_water_column_diagfct
329
330 ! --------------------------------------------------------------------
331 ! ComIn callback function.
332 ! --------------------------------------------------------------------
333
334 SUBROUTINE calc_water_column_destructor() BIND(C)
335
336 IMPLICIT NONE
337
338 CALL comin_print_info(' - destructor.')
339
340 END SUBROUTINE calc_water_column_destructor
341
342 !---------------------------------------------------------------------
343 !---------------------------------------------------------------------
344 !---------------------------------------------------------------------
345 ! PRIVATE ROUTINES
346 !---------------------------------------------------------------------
347 !---------------------------------------------------------------------
348 !---------------------------------------------------------------------
349
350 SUBROUTINE comin_var_request_add_wrapper(descriptor, lmode_exclusive, zaxis_id &
351 , tracer, restart, units)
352
353 IMPLICIT NONE
354
355 ! I/O
356 TYPE(t_comin_var_descriptor), INTENT(IN) :: descriptor
357 LOGICAL, OPTIONAL, INTENT(IN) :: lmode_exclusive
358 INTEGER, OPTIONAL, INTENT(IN) :: zaxis_id
359 LOGICAL, OPTIONAL, INTENT(IN) :: tracer
360 LOGICAL, OPTIONAL, INTENT(IN) :: restart
361 CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: units
362 ! LOCAL
363 LOGICAL :: lexclusive
364
365 IF (PRESENT(lmode_exclusive)) THEN
366 lexclusive = lmode_exclusive
367 ELSE
368 lexclusive = .false.
369 END IF
370
371 CALL comin_var_request_add(descriptor, lexclusive)
372
373 IF (PRESENT(zaxis_id)) THEN
374 CALL comin_metadata_set(descriptor, "zaxis_id", zaxis_id)
375 END IF
376 IF (PRESENT(tracer)) THEN
377 CALL comin_metadata_set(descriptor, "tracer", tracer)
378 END IF
379 IF (PRESENT(restart)) THEN
380 CALL comin_metadata_set(descriptor, "restart", restart)
381 END IF
382 IF (PRESENT(units)) THEN
383 CALL comin_metadata_set(descriptor, "units", trim(units))
384 END IF
385
386 END SUBROUTINE comin_var_request_add_wrapper
387
388 ! --------------------------------------------------------------------------------------------------
389
390 SUBROUTINE check_variable(var, name, dim)
391
392 IMPLICIT NONE
393
394 ! I/O
395 TYPE(t_comin_var_handle), INTENT(IN) :: var
396 CHARACTER(LEN=*), INTENT(IN) :: name
397 INTEGER, INTENT(IN), OPTIONAL :: dim
398 ! LOCAL
399 INTEGER :: idim
400
401 IF (PRESENT(dim)) THEN
402 idim = dim
403 ELSE
404 idim = 3
405 END IF
406 SELECT CASE (idim)
407 CASE (2)
408 CALL check_dimensions_2d(var, name)
409 CASE (3)
410 CALL check_dimensions_3d(var, name)
411 CASE DEFAULT
412 CALL comin_plugin_finish (pluginname//' unsupported dimension', 'check variable')
413 END SELECT
414
415 END SUBROUTINE check_variable
416
417 ! --------------------------------------------------------------------------------------------------
418
419 SUBROUTINE check_dimensions_3d(var, name)
420
421 IMPLICIT NONE
422
423 TYPE(t_comin_var_handle), INTENT(IN) :: var
424 CHARACTER(LEN=*), INTENT(IN) :: name
425
426 IF (any(var%dim_semantics() /= [comin_dim_semantics_nproma, &
427 & comin_dim_semantics_level, &
428 & comin_dim_semantics_block, &
429 & comin_dim_semantics_unused, &
430 & comin_dim_semantics_unused] )) &
431 & CALL comin_plugin_finish(pluginname//': '//trim(name) &
432 & , "Dimension check failed!")
433
434 END SUBROUTINE check_dimensions_3d
435
436 ! --------------------------------------------------------------------------------------------------
437
438 ! --------------------------------------------------------------------------------------------------
439
440 SUBROUTINE check_dimensions_2d(var, name)
441
442 IMPLICIT NONE
443
444 TYPE(t_comin_var_handle), INTENT(IN) :: var
445 CHARACTER(LEN=*), INTENT(IN) :: name
446
447 IF (any(var%dim_semantics() /= [comin_dim_semantics_nproma, &
448 & comin_dim_semantics_block, &
449 & comin_dim_semantics_unused, &
450 & comin_dim_semantics_unused, &
451 & comin_dim_semantics_unused] )) &
452 & CALL comin_plugin_finish(pluginname//': '//trim(name) &
453 , "Dimension check failed!")
454
455 END SUBROUTINE check_dimensions_2d
456
457 ! --------------------------------------------------------------------------------------------------
458
459 SUBROUTINE message(text, idom, lall)
460
461 IMPLICIT NONE
462
463 CHARACTER(LEN=*) :: text
464 INTEGER, OPTIONAL :: idom
465 LOGICAL, OPTIONAL :: lall
466
467 ! LOCAL
468 LOGICAL :: la
469
470 IF (PRESENT(lall)) THEN
471 la = lall
472 ELSE
473 la = .false.
474 END IF
475
476 IF (comin_parallel_get_host_mpi_rank() == 0 .OR. la) THEN
477 IF (PRESENT(idom)) THEN
478 WRITE (0,fmt='(A,I2,A)') pluginname//' domain: ',idom,' '//trim(text)
479 ELSE
480 WRITE (0,fmt='(A,A)') pluginname//' ',trim(text)
481 END IF
482 END IF
483
484 END SUBROUTINE message
485
487! --------------------------------------------------------------------------------------------------
488! --------------------------------------------------------------------------------------------------
t_comin_error_code comin_error_get()
Get the current ComIn error code.
void comin_plugin_finish(const char *routine, const char *text)
const T * comin_metadata_get(t_comin_var_descriptor descriptor, const char *key)
void comin_error_set_errors_return(bool errors_return)
void comin_setup_get_version(unsigned int *major, unsigned int *minor, unsigned int *patch)
int comin_current_get_domain_id()
void comin_var_request_add(t_comin_var_descriptor var_desc, bool lmodexclusive)
void comin_callback_register(t_comin_entry_point entry_point, t_comin_callback_function fct_ptr)
Example plugin for the ICON Community Interface (ComIn)
subroutine, public calc_water_column_diagfct()
subroutine, public calc_water_column_destructor()
type(t_plugin_vars), dimension(:), allocatable pvpp
subroutine, public calc_water_column_constructor()