source: palm/trunk/SOURCE/sum_up_3d_data.f90 @ 1007

Last change on this file since 1007 was 1007, checked in by franke, 9 years ago

Bugfixes


Missing calculation of mean particle weighting factor for output added. (data_output_2d, data_output_3d, data_output_mask, sum_up_3d_data)
Calculation of mean particle radius for output now considers the weighting factor. (data_output_mask)
Calculation of sugrid-scale buoyancy flux for humidity and cloud droplets corrected. (flow_statistics)
Factor in calculation of enhancement factor for collision efficencies corrected. (lpm_collision_kernels)
Calculation of buoyancy production now considers the liquid water mixing ratio in case of cloud droplets. (production_e)

Changes


Calculation of buoyancy flux for humidity in case of WS-scheme is now using turbulent fluxes of WS-scheme. (flow_statistics)
Calculation of the collision kernels now in SI units. (lpm_collision_kernels)

  • Property svn:keywords set to Id
File size: 15.4 KB
Line 
1 SUBROUTINE sum_up_3d_data
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6! Bugfix in calculation of ql_vp
7!
8! Former revisions:
9! -----------------
10! $Id: sum_up_3d_data.f90 1007 2012-09-19 14:30:36Z franke $
11!
12! 978 2012-08-09 08:28:32Z fricke
13! +z0h*
14!
15! 790 2011-11-29 03:11:20Z raasch
16! bugfix: calculation of 'pr' must depend on the particle weighting factor
17!
18! 771 2011-10-27 10:56:21Z heinze
19! +lpt_av
20!
21! 667 2010-12-23 12:06:00Z suehring/gryschka
22! nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng
23!
24! 402 2009-10-21 11:59:41Z maronga
25! Bugfix in calculation of shf*_av, qsws*_av
26!
27! 2009-08-25 08:35:52Z maronga
28! +shf*, qsws*
29!
30! 96 2007-06-04 08:07:41Z raasch
31! +sum-up of density and salinity
32!
33! 72 2007-03-19 08:20:46Z raasch
34! +sum-up of precipitation rate and roughness length (prr*, z0*)
35!
36! RCS Log replace by Id keyword, revision history cleaned up
37!
38! Revision 1.1  2006/02/23 12:55:23  raasch
39! Initial revision
40!
41!
42! Description:
43! ------------
44! Sum-up the values of 3d-arrays. The real averaging is later done in routine
45! average_3d_data.
46!------------------------------------------------------------------------------!
47
48    USE arrays_3d
49    USE averaging
50    USE cloud_parameters
51    USE control_parameters
52    USE cpulog
53    USE indices
54    USE interfaces
55    USE particle_attributes
56
57    IMPLICIT NONE
58
59    INTEGER ::  i, ii, j, k, n, psi
60
61    REAL    ::  mean_r, s_r3, s_r4
62
63
64    CALL cpu_log (log_point(34),'sum_up_3d_data','start')
65
66!
67!-- Allocate and initialize the summation arrays if called for the very first
68!-- time or the first time after average_3d_data has been called
69!-- (some or all of the arrays may have been already allocated
70!-- in read_3d_binary)
71    IF ( average_count_3d == 0 )  THEN
72
73       DO  ii = 1, doav_n
74
75          SELECT CASE ( TRIM( doav(ii) ) )
76
77             CASE ( 'e' )
78                IF ( .NOT. ALLOCATED( e_av ) )  THEN
79                   ALLOCATE( e_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
80                ENDIF
81                e_av = 0.0
82
83             CASE ( 'lpt' )
84                IF ( .NOT. ALLOCATED( lpt_av ) )  THEN
85                   ALLOCATE( lpt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
86                ENDIF
87                lpt_av = 0.0
88
89             CASE ( 'lwp*' )
90                IF ( .NOT. ALLOCATED( lwp_av ) )  THEN
91                   ALLOCATE( lwp_av(nysg:nyng,nxlg:nxrg) )
92                ENDIF
93                lwp_av = 0.0
94
95             CASE ( 'p' )
96                IF ( .NOT. ALLOCATED( p_av ) )  THEN
97                   ALLOCATE( p_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
98                ENDIF
99                p_av = 0.0
100
101             CASE ( 'pc' )
102                IF ( .NOT. ALLOCATED( pc_av ) )  THEN
103                   ALLOCATE( pc_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
104                ENDIF
105                pc_av = 0.0
106
107             CASE ( 'pr' )
108                IF ( .NOT. ALLOCATED( pr_av ) )  THEN
109                   ALLOCATE( pr_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
110                ENDIF
111                pr_av = 0.0
112
113             CASE ( 'prr*' )
114                IF ( .NOT. ALLOCATED( precipitation_rate_av ) )  THEN
115                   ALLOCATE( precipitation_rate_av(nysg:nyng,nxlg:nxrg) )
116                ENDIF
117                precipitation_rate_av = 0.0
118
119             CASE ( 'pt' )
120                IF ( .NOT. ALLOCATED( pt_av ) )  THEN
121                   ALLOCATE( pt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
122                ENDIF
123                pt_av = 0.0
124
125             CASE ( 'q' )
126                IF ( .NOT. ALLOCATED( q_av ) )  THEN
127                   ALLOCATE( q_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
128                ENDIF
129                q_av = 0.0
130
131             CASE ( 'ql' )
132                IF ( .NOT. ALLOCATED( ql_av ) )  THEN
133                   ALLOCATE( ql_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
134                ENDIF
135                ql_av = 0.0
136
137             CASE ( 'ql_c' )
138                IF ( .NOT. ALLOCATED( ql_c_av ) )  THEN
139                   ALLOCATE( ql_c_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
140                ENDIF
141                ql_c_av = 0.0
142
143             CASE ( 'ql_v' )
144                IF ( .NOT. ALLOCATED( ql_v_av ) )  THEN
145                   ALLOCATE( ql_v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
146                ENDIF
147                ql_v_av = 0.0
148
149             CASE ( 'ql_vp' )
150                IF ( .NOT. ALLOCATED( ql_vp_av ) )  THEN
151                   ALLOCATE( ql_vp_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
152                ENDIF
153                ql_vp_av = 0.0
154
155             CASE ( 'qsws*' )
156                IF ( .NOT. ALLOCATED( qsws_av ) )  THEN
157                   ALLOCATE( qsws_av(nysg:nyng,nxlg:nxrg) )
158                ENDIF
159                qsws_av = 0.0
160
161             CASE ( 'qv' )
162                IF ( .NOT. ALLOCATED( qv_av ) )  THEN
163                   ALLOCATE( qv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
164                ENDIF
165                qv_av = 0.0
166
167             CASE ( 'rho' )
168                IF ( .NOT. ALLOCATED( rho_av ) )  THEN
169                   ALLOCATE( rho_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
170                ENDIF
171                rho_av = 0.0
172
173             CASE ( 's' )
174                IF ( .NOT. ALLOCATED( s_av ) )  THEN
175                   ALLOCATE( s_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
176                ENDIF
177                s_av = 0.0
178
179             CASE ( 'sa' )
180                IF ( .NOT. ALLOCATED( sa_av ) )  THEN
181                   ALLOCATE( sa_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
182                ENDIF
183                sa_av = 0.0
184
185             CASE ( 'shf*' )
186                IF ( .NOT. ALLOCATED( shf_av ) )  THEN
187                   ALLOCATE( shf_av(nysg:nyng,nxlg:nxrg) )
188                ENDIF
189                shf_av = 0.0
190
191             CASE ( 't*' )
192                IF ( .NOT. ALLOCATED( ts_av ) )  THEN
193                   ALLOCATE( ts_av(nysg:nyng,nxlg:nxrg) )
194                ENDIF
195                ts_av = 0.0
196
197             CASE ( 'u' )
198                IF ( .NOT. ALLOCATED( u_av ) )  THEN
199                   ALLOCATE( u_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
200                ENDIF
201                u_av = 0.0
202
203             CASE ( 'u*' )
204                IF ( .NOT. ALLOCATED( us_av ) )  THEN
205                   ALLOCATE( us_av(nysg:nyng,nxlg:nxrg) )
206                ENDIF
207                us_av = 0.0
208
209             CASE ( 'v' )
210                IF ( .NOT. ALLOCATED( v_av ) )  THEN
211                   ALLOCATE( v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
212                ENDIF
213                v_av = 0.0
214
215             CASE ( 'vpt' )
216                IF ( .NOT. ALLOCATED( vpt_av ) )  THEN
217                   ALLOCATE( vpt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
218                ENDIF
219                vpt_av = 0.0
220
221             CASE ( 'w' )
222                IF ( .NOT. ALLOCATED( w_av ) )  THEN
223                   ALLOCATE( w_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
224                ENDIF
225                w_av = 0.0
226
227             CASE ( 'z0*' )
228                IF ( .NOT. ALLOCATED( z0_av ) )  THEN
229                   ALLOCATE( z0_av(nysg:nyng,nxlg:nxrg) )
230                ENDIF
231                z0_av = 0.0
232
233             CASE ( 'z0h*' )
234                IF ( .NOT. ALLOCATED( z0h_av ) )  THEN
235                   ALLOCATE( z0h_av(nysg:nyng,nxlg:nxrg) )
236                ENDIF
237                z0h_av = 0.0
238
239             CASE DEFAULT
240!
241!--             User-defined quantity
242                CALL user_3d_data_averaging( 'allocate', doav(ii) )
243
244          END SELECT
245
246       ENDDO
247
248    ENDIF
249
250!
251!-- Loop of all variables to be averaged.
252    DO  ii = 1, doav_n
253
254!
255!--    Store the array chosen on the temporary array.
256       SELECT CASE ( TRIM( doav(ii) ) )
257
258          CASE ( 'e' )
259             DO  i = nxlg, nxrg
260                DO  j = nysg, nyng
261                   DO  k = nzb, nzt+1
262                      e_av(k,j,i) = e_av(k,j,i) + e(k,j,i)
263                   ENDDO
264                ENDDO
265             ENDDO
266
267          CASE ( 'lpt' )
268             DO  i = nxlg, nxrg
269                DO  j = nysg, nyng
270                   DO  k = nzb, nzt+1
271                      lpt_av(k,j,i) = lpt_av(k,j,i) + pt(k,j,i)
272                   ENDDO
273                ENDDO
274             ENDDO
275
276          CASE ( 'lwp*' )
277             DO  i = nxlg, nxrg
278                DO  j = nysg, nyng
279                   lwp_av(j,i) = lwp_av(j,i) + SUM( ql(nzb:nzt,j,i) * &
280                                                    dzw(1:nzt+1) )
281                ENDDO
282             ENDDO
283
284          CASE ( 'p' )
285             DO  i = nxlg, nxrg
286                DO  j = nysg, nyng
287                   DO  k = nzb, nzt+1
288                      p_av(k,j,i) = p_av(k,j,i) + p(k,j,i)
289                   ENDDO
290                ENDDO
291             ENDDO
292
293          CASE ( 'pc' )
294             DO  i = nxl, nxr
295                DO  j = nys, nyn
296                   DO  k = nzb, nzt+1
297                      pc_av(k,j,i) = pc_av(k,j,i) + prt_count(k,j,i)
298                   ENDDO
299                ENDDO
300             ENDDO
301
302          CASE ( 'pr' )
303             DO  i = nxl, nxr
304                DO  j = nys, nyn
305                   DO  k = nzb, nzt+1
306                      psi = prt_start_index(k,j,i)
307                      s_r3 = 0.0
308                      s_r4 = 0.0
309                      DO  n = psi, psi+prt_count(k,j,i)-1
310                         s_r3 = s_r3 + particles(n)%radius**3 * &
311                                       particles(n)%weight_factor
312                         s_r4 = s_r4 + particles(n)%radius**4 * &
313                                       particles(n)%weight_factor
314                      ENDDO
315                      IF ( s_r3 /= 0.0 )  THEN
316                         mean_r = s_r4 / s_r3
317                      ELSE
318                         mean_r = 0.0
319                      ENDIF
320                      pr_av(k,j,i) = pr_av(k,j,i) + mean_r
321                   ENDDO
322                ENDDO
323             ENDDO
324
325          CASE ( 'pr*' )
326             DO  i = nxlg, nxrg
327                DO  j = nysg, nyng
328                   precipitation_rate_av(j,i) = precipitation_rate_av(j,i) + &
329                                                precipitation_rate(j,i)
330                ENDDO
331             ENDDO
332
333          CASE ( 'pt' )
334             IF ( .NOT. cloud_physics ) THEN
335             DO  i = nxlg, nxrg
336                DO  j = nysg, nyng
337                   DO  k = nzb, nzt+1
338                         pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i)
339                      ENDDO
340                   ENDDO
341                ENDDO
342             ELSE
343             DO  i = nxlg, nxrg
344                DO  j = nysg, nyng
345                   DO  k = nzb, nzt+1
346                         pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i) + l_d_cp * &
347                                                       pt_d_t(k) * ql(k,j,i)
348                      ENDDO
349                   ENDDO
350                ENDDO
351             ENDIF
352
353          CASE ( 'q' )
354             DO  i = nxlg, nxrg
355                DO  j = nysg, nyng
356                   DO  k = nzb, nzt+1
357                      q_av(k,j,i) = q_av(k,j,i) + q(k,j,i)
358                   ENDDO
359                ENDDO
360             ENDDO
361
362          CASE ( 'ql' )
363             DO  i = nxlg, nxrg
364                DO  j = nysg, nyng
365                   DO  k = nzb, nzt+1
366                      ql_av(k,j,i) = ql_av(k,j,i) + ql(k,j,i)
367                   ENDDO
368                ENDDO
369             ENDDO
370
371          CASE ( 'ql_c' )
372             DO  i = nxlg, nxrg
373                DO  j = nysg, nyng
374                   DO  k = nzb, nzt+1
375                      ql_c_av(k,j,i) = ql_c_av(k,j,i) + ql_c(k,j,i)
376                   ENDDO
377                ENDDO
378             ENDDO
379
380          CASE ( 'ql_v' )
381             DO  i = nxlg, nxrg
382                DO  j = nysg, nyng
383                   DO  k = nzb, nzt+1
384                      ql_v_av(k,j,i) = ql_v_av(k,j,i) + ql_v(k,j,i)
385                   ENDDO
386                ENDDO
387             ENDDO
388
389          CASE ( 'ql_vp' )
390             DO  i = nxl, nxr
391                DO  j = nys, nyn
392                   DO  k = nzb, nzt+1
393                      psi = prt_start_index(k,j,i)
394                      DO  n = psi, psi+prt_count(k,j,i)-1
395                         ql_vp_av(k,j,i) = ql_vp_av(k,j,i) + &
396                                           particles(n)%weight_factor / &
397                                           prt_count(k,j,i)
398                      ENDDO
399                   ENDDO
400                ENDDO
401             ENDDO
402
403          CASE ( 'qsws*' )
404             DO  i = nxlg, nxrg
405                DO  j = nysg, nyng
406                   qsws_av(j,i) = qsws_av(j,i) + qsws(j,i)
407                ENDDO
408             ENDDO
409
410          CASE ( 'qv' )
411             DO  i = nxlg, nxrg
412                DO  j = nysg, nyng
413                   DO  k = nzb, nzt+1
414                      qv_av(k,j,i) = qv_av(k,j,i) + q(k,j,i) - ql(k,j,i)
415                   ENDDO
416                ENDDO
417             ENDDO
418
419          CASE ( 'rho' )
420             DO  i = nxlg, nxrg
421                DO  j = nysg, nyng
422                   DO  k = nzb, nzt+1
423                      rho_av(k,j,i) = rho_av(k,j,i) + rho(k,j,i)
424                   ENDDO
425                ENDDO
426             ENDDO
427
428          CASE ( 's' )
429             DO  i = nxlg, nxrg
430                DO  j = nysg, nyng
431                   DO  k = nzb, nzt+1
432                      s_av(k,j,i) = s_av(k,j,i) + q(k,j,i)
433                   ENDDO
434                ENDDO
435             ENDDO
436
437          CASE ( 'sa' )
438             DO  i = nxlg, nxrg
439                DO  j = nysg, nyng
440                   DO  k = nzb, nzt+1
441                      sa_av(k,j,i) = sa_av(k,j,i) + sa(k,j,i)
442                   ENDDO
443                ENDDO
444             ENDDO
445
446          CASE ( 'shf*' )
447             DO  i = nxlg, nxrg
448                DO  j = nysg, nyng
449                   shf_av(j,i) = shf_av(j,i) + shf(j,i)
450                ENDDO
451             ENDDO
452
453          CASE ( 't*' )
454             DO  i = nxlg, nxrg
455                DO  j = nysg, nyng
456                   ts_av(j,i) = ts_av(j,i) + ts(j,i)
457                ENDDO
458             ENDDO
459
460          CASE ( 'u' )
461             DO  i = nxlg, nxrg
462                DO  j = nysg, nyng
463                   DO  k = nzb, nzt+1
464                      u_av(k,j,i) = u_av(k,j,i) + u(k,j,i)
465                   ENDDO
466                ENDDO
467             ENDDO
468
469          CASE ( 'u*' )
470             DO  i = nxlg, nxrg
471                DO  j = nysg, nyng
472                   us_av(j,i) = us_av(j,i) + us(j,i)
473                ENDDO
474             ENDDO
475
476          CASE ( 'v' )
477             DO  i = nxlg, nxrg
478                DO  j = nysg, nyng
479                   DO  k = nzb, nzt+1
480                      v_av(k,j,i) = v_av(k,j,i) + v(k,j,i)
481                   ENDDO
482                ENDDO
483             ENDDO
484
485          CASE ( 'vpt' )
486             DO  i = nxlg, nxrg
487                DO  j = nysg, nyng
488                   DO  k = nzb, nzt+1
489                      vpt_av(k,j,i) = vpt_av(k,j,i) + vpt(k,j,i)
490                   ENDDO
491                ENDDO
492             ENDDO
493
494          CASE ( 'w' )
495             DO  i = nxlg, nxrg
496                DO  j = nysg, nyng
497                   DO  k = nzb, nzt+1
498                      w_av(k,j,i) = w_av(k,j,i) + w(k,j,i)
499                   ENDDO
500                ENDDO
501             ENDDO
502
503          CASE ( 'z0*' )
504             DO  i = nxlg, nxrg
505                DO  j = nysg, nyng
506                   z0_av(j,i) = z0_av(j,i) + z0(j,i)
507                ENDDO
508             ENDDO
509
510          CASE ( 'z0h*' )
511             DO  i = nxlg, nxrg
512                DO  j = nysg, nyng
513                   z0h_av(j,i) = z0h_av(j,i) + z0h(j,i)
514                ENDDO
515             ENDDO
516
517          CASE DEFAULT
518!
519!--          User-defined quantity
520             CALL user_3d_data_averaging( 'sum', doav(ii) )
521
522       END SELECT
523
524    ENDDO
525
526    CALL cpu_log (log_point(34),'sum_up_3d_data','stop','nobarrier')
527
528
529 END SUBROUTINE sum_up_3d_data
Note: See TracBrowser for help on using the repository browser.