source: palm/trunk/SOURCE/diagnostic_output_quantities_mod.f90 @ 4027

Last change on this file since 4027 was 3998, checked in by suehring, 5 years ago

Bugfix in output module for diagnostic quantities

  • Property svn:keywords set to Id
File size: 11.9 KB
Line 
1!> @file diagnostic_output_quantities_mod.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2019 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! ------------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: diagnostic_output_quantities_mod.f90 3998 2019-05-23 13:38:11Z schwenkel $
27! Bugfix in gathering all output strings
28!
29! 3995 2019-05-22 18:59:54Z suehring
30! Avoid compiler warnings about unused variable and fix string operation which
31! is not allowed with PGI compiler
32!
33! 3994 2019-05-22 18:08:09Z suehring
34! Initial revision
35!
36!
37! @author Farah Kanani-Suehring
38!
39! Description:
40! ------------
41!> ...
42!------------------------------------------------------------------------------!
43 MODULE diagnostic_output_quantities_mod
44 
45
46!     USE arrays_3d,                                                             &
47!         ONLY:  dzw, e, heatflux_output_conversion, nc, nr, p, pt,              &
48!                precipitation_amount, prr, q, qc, ql, ql_c, ql_v, qr, s, tend,  &
49!                u, v, vpt, w, zu, zw, waterflux_output_conversion, hyrho, d_exner
50!
51!     USE averaging
52!
53!     USE basic_constants_and_equations_mod,                                     &
54!         ONLY:  c_p, lv_d_cp, l_v
55!
56!     USE bulk_cloud_model_mod,                                                  &
57!         ONLY:  bulk_cloud_model
58!
59    USE control_parameters,                                                    &
60        ONLY:  current_timestep_number, varnamelength
61!
62!     USE cpulog,                                                                &
63!         ONLY:  cpu_log, log_point
64!
65!     USE indices,                                                               &
66!         ONLY:  nbgp, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg,       &
67!                nzb, nzt, wall_flags_0
68!
69    USE kinds
70!
71!     USE land_surface_model_mod,                                                &
72!         ONLY:  zs
73!
74!     USE module_interface,                                                      &
75!         ONLY:  module_interface_data_output_2d
76!
77! #if defined( __netcdf )
78!     USE NETCDF
79! #endif
80!
81!     USE netcdf_interface,                                                      &
82!         ONLY:  fill_value, id_set_xy, id_set_xz, id_set_yz, id_var_do2d,       &
83!                id_var_time_xy, id_var_time_xz, id_var_time_yz, nc_stat,        &
84!                netcdf_data_format, netcdf_handle_error
85!
86!     USE particle_attributes,                                                   &
87!         ONLY:  grid_particles, number_of_particles, particle_advection_start,  &
88!                particles, prt_count
89!     
90!     USE pegrid
91!
92!     USE surface_mod,                                                           &
93!         ONLY:  ind_pav_green, ind_veg_wall, ind_wat_win, surf_def_h,           &
94!                surf_lsm_h, surf_usm_h
95!
96!     USE turbulence_closure_mod,                                                &
97!         ONLY:  tcm_data_output_2d
98
99
100    IMPLICIT NONE
101
102    CHARACTER(LEN=varnamelength), DIMENSION(500) ::  do_all = ' '
103
104    INTEGER(iwp) ::  timestep_number_at_prev_calc = 0  !< ...at previous diagnostic output calculation
105
106    LOGICAL ::  initialized_diagnostic_output_quantities = .FALSE.
107    LOGICAL ::  prepared_diagnostic_output_quantities = .FALSE.
108
109    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ti     !< rotation(u,v,w) aka turbulence intensity
110    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET ::  ti_av  !< avg. rotation(u,v,w) aka turbulence intensity
111
112
113    SAVE
114
115    PRIVATE
116
117!
118!-- Public variables
119    PUBLIC do_all,                                                                                 &
120           initialized_diagnostic_output_quantities,                                               &
121           prepared_diagnostic_output_quantities,                                                  &
122           ti, ti_av,                                                                              &
123           timestep_number_at_prev_calc
124!
125!-- Public routines
126    PUBLIC diagnostic_output_quantities_calculate,                                                 &
127           diagnostic_output_quantities_init,                                                      &
128           diagnostic_output_quantities_prepare
129
130
131    INTERFACE diagnostic_output_quantities_init
132       MODULE PROCEDURE diagnostic_output_quantities_init
133    END INTERFACE diagnostic_output_quantities_init
134
135    INTERFACE diagnostic_output_quantities_calculate
136       MODULE PROCEDURE diagnostic_output_quantities_calculate
137    END INTERFACE diagnostic_output_quantities_calculate
138
139    INTERFACE diagnostic_output_quantities_prepare
140       MODULE PROCEDURE diagnostic_output_quantities_prepare
141    END INTERFACE diagnostic_output_quantities_prepare
142
143
144 CONTAINS
145
146!------------------------------------------------------------------------------!
147! Description:
148! ------------
149!> ...
150!------------------------------------------------------------------------------!
151 SUBROUTINE diagnostic_output_quantities_init
152
153
154    USE indices,                                                                                   &
155        ONLY:  nxl, nxr, nyn, nys, nzb, nzt
156
157    IMPLICIT NONE
158!
159!-- Next line is to avoid compiler warnings about unused variables
160    IF ( timestep_number_at_prev_calc == 0 )  CONTINUE
161
162    initialized_diagnostic_output_quantities = .FALSE.
163
164    IF ( .NOT. ALLOCATED( ti ) )  THEN
165       ALLOCATE( ti(nzb:nzt+1,nys:nyn,nxl:nxr) )
166       ti = 0.0_wp
167    ENDIF
168
169    initialized_diagnostic_output_quantities = .TRUE.
170
171 END SUBROUTINE diagnostic_output_quantities_init
172
173
174!--------------------------------------------------------------------------------------------------!
175! Description:
176! ------------
177!> ...
178!--------------------------------------------------------------------------------------------------!
179 SUBROUTINE diagnostic_output_quantities_calculate
180
181
182    USE arrays_3d,                                                                                 &
183        ONLY:  ddzu, u, v, w
184
185    USE grid_variables,                                                                            &
186        ONLY:  ddx, ddy
187
188    USE indices,                                                                                   &
189        ONLY:  nxl, nxr, nyn, nys, nzb, nzt, wall_flags_0
190
191    IMPLICIT NONE
192
193    INTEGER(iwp) ::  i, j, k    !< grid loop index in x-, y-, z-direction
194    INTEGER(iwp) ::  ivar_all   !< loop index over all 2d/3d/mask output quantities
195
196
197!     CALL cpu_log( log_point(41), 'calculate_quantities', 'start' )
198!
199!-- Preparatory steps and initialization of output arrays
200    IF ( .NOT.  prepared_diagnostic_output_quantities )     CALL diagnostic_output_quantities_prepare
201    IF ( .NOT.  initialized_diagnostic_output_quantities )  CALL diagnostic_output_quantities_init
202!
203!-- Save timestep number to check in time_integration if diagnostic_output_quantities_calculate
204!-- has been called already, since the CALL occurs at two locations, but the calculations need to be
205!-- done only once per timestep.
206    timestep_number_at_prev_calc = current_timestep_number
207
208
209    ivar_all = 1
210
211    DO  WHILE ( do_all(ivar_all)(1:1) /= ' ' )
212
213       SELECT CASE ( TRIM( do_all(ivar_all) ) )
214!
215!--       Calculate 'turbulence intensity' from rot[(u,v,w)] at scalar grid point
216          CASE ( 'ti' )
217             DO  i = nxl, nxr
218                DO  j = nys, nyn
219                   DO  k = nzb+1, nzt
220
221                      ti(k,j,i) = 0.25_wp * SQRT(                                                  &
222                        (   (   w(k,j+1,i) + w(k-1,j+1,i)                                          &
223                              - w(k,j-1,i) - w(k-1,j-1,i) ) * ddy                                  &
224                          - (   v(k+1,j,i) + v(k+1,j+1,i)                                          &
225                              - v(k-1,j,i) - v(k-1,j+1,i) ) * ddzu(k) )**2                         &
226                      + (   (   u(k+1,j,i) + u(k+1,j,i+1)                                          &
227                              - u(k-1,j,i) - u(k-1,j,i+1) ) * ddzu(k)                              &
228                          - (   w(k,j,i+1) + w(k-1,j,i+1)                                          &
229                              - w(k,j,i-1) - w(k-1,j,i-1) ) * ddx     )**2                         &
230                      + (   (   v(k,j,i+1) + v(k,j+1,i+1)                                          &
231                              - v(k,j,i-1) - v(k,j+1,i-1) ) * ddx                                  &
232                          - (   u(k,j+1,i) + u(k,j+1,i+1)                                          &
233                              - u(k,j-1,i) - u(k,j-1,i+1) ) * ddy     )**2  )                      &
234                                          * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 0) )
235                   ENDDO
236                ENDDO
237             ENDDO
238
239          END SELECT
240
241          ivar_all = ivar_all + 1
242    ENDDO
243
244!     CALL cpu_log( log_point(41), 'calculate_quantities', 'stop' )
245
246 END SUBROUTINE diagnostic_output_quantities_calculate
247
248
249!------------------------------------------------------------------------------!
250! Description:
251! ------------
252!> ...
253!------------------------------------------------------------------------------!
254 SUBROUTINE diagnostic_output_quantities_prepare
255
256
257    USE control_parameters,                                                                        &
258        ONLY:  do2d, do3d, domask, masks, mid
259
260    IMPLICIT NONE
261
262    CHARACTER (LEN=varnamelength), DIMENSION(0:1,500) ::  do2d_var = ' '  !<
263                                                          !< label array for 2d output quantities
264
265    INTEGER(iwp) ::  av         !< index defining type of output, av=0 instantaneous, av=1 averaged
266    INTEGER(iwp) ::  ivar       !< loop index
267    INTEGER(iwp) ::  ivar_all   !< loop index
268    INTEGER(iwp) ::  l          !< index for cutting string
269
270
271    prepared_diagnostic_output_quantities = .FALSE.
272
273    ivar     = 1
274    ivar_all = 1
275
276    DO  av = 0, 1
277!
278!--    Remove _xy, _xz, or _yz from string
279       l = MAX( 3, LEN_TRIM( do2d(av,ivar) ) )
280       do2d_var(av,ivar)(1:l-3) = do2d(av,ivar)(1:l-3)
281!
282!--    Gather 2d output quantity names, check for double occurrence of output quantity
283       DO  WHILE ( do2d_var(av,ivar)(1:1) /= ' ' )
284          IF ( .NOT.  ANY( do_all == do2d_var(av,ivar) ) )  THEN
285             do_all(ivar_all) = do2d_var(av,ivar)
286          ENDIF
287          ivar = ivar + 1
288          ivar_all = ivar_all + 1
289          l = MAX( 3, LEN_TRIM( do2d(av,ivar) ) )
290          do2d_var(av,ivar)(1:l-3) = do2d(av,ivar)(1:l-3)
291       ENDDO
292
293       ivar = 1
294!
295!--    Gather 3d output quantity names, check for double occurrence of output quantity
296       DO  WHILE ( do3d(av,ivar)(1:1) /= ' ' )
297          do_all(ivar_all) = do3d(av,ivar)
298
299          ivar = ivar + 1
300          ivar_all = ivar_all + 1
301       ENDDO
302
303       ivar = 1
304!
305!--    Gather masked output quantity names, check for double occurrence of output quantity
306       DO  mid = 1, masks
307          DO  WHILE ( domask(mid,av,ivar)(1:1) /= ' ' )
308             do_all(ivar_all) = domask(mid,av,ivar)
309
310             ivar = ivar + 1
311             ivar_all = ivar_all + 1
312          ENDDO
313          ivar = 1
314       ENDDO
315
316    ENDDO
317
318    prepared_diagnostic_output_quantities = .TRUE.
319
320 END SUBROUTINE diagnostic_output_quantities_prepare
321
322 END MODULE diagnostic_output_quantities_mod
Note: See TracBrowser for help on using the repository browser.