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

Last change on this file since 985 was 979, checked in by fricke, 12 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 15.2 KB
Line 
1 SUBROUTINE sum_up_3d_data
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: sum_up_3d_data.f90 979 2012-08-09 08:50:11Z maronga $
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 = nxlg, nxrg
391                DO  j = nysg, nyng
392                   DO  k = nzb, nzt+1
393                      ql_vp_av(k,j,i) = ql_vp_av(k,j,i) + ql_vp(k,j,i)
394                   ENDDO
395                ENDDO
396             ENDDO
397
398          CASE ( 'qsws*' )
399             DO  i = nxlg, nxrg
400                DO  j = nysg, nyng
401                   qsws_av(j,i) = qsws_av(j,i) + qsws(j,i)
402                ENDDO
403             ENDDO
404
405          CASE ( 'qv' )
406             DO  i = nxlg, nxrg
407                DO  j = nysg, nyng
408                   DO  k = nzb, nzt+1
409                      qv_av(k,j,i) = qv_av(k,j,i) + q(k,j,i) - ql(k,j,i)
410                   ENDDO
411                ENDDO
412             ENDDO
413
414          CASE ( 'rho' )
415             DO  i = nxlg, nxrg
416                DO  j = nysg, nyng
417                   DO  k = nzb, nzt+1
418                      rho_av(k,j,i) = rho_av(k,j,i) + rho(k,j,i)
419                   ENDDO
420                ENDDO
421             ENDDO
422
423          CASE ( 's' )
424             DO  i = nxlg, nxrg
425                DO  j = nysg, nyng
426                   DO  k = nzb, nzt+1
427                      s_av(k,j,i) = s_av(k,j,i) + q(k,j,i)
428                   ENDDO
429                ENDDO
430             ENDDO
431
432          CASE ( 'sa' )
433             DO  i = nxlg, nxrg
434                DO  j = nysg, nyng
435                   DO  k = nzb, nzt+1
436                      sa_av(k,j,i) = sa_av(k,j,i) + sa(k,j,i)
437                   ENDDO
438                ENDDO
439             ENDDO
440
441          CASE ( 'shf*' )
442             DO  i = nxlg, nxrg
443                DO  j = nysg, nyng
444                   shf_av(j,i) = shf_av(j,i) + shf(j,i)
445                ENDDO
446             ENDDO
447
448          CASE ( 't*' )
449             DO  i = nxlg, nxrg
450                DO  j = nysg, nyng
451                   ts_av(j,i) = ts_av(j,i) + ts(j,i)
452                ENDDO
453             ENDDO
454
455          CASE ( 'u' )
456             DO  i = nxlg, nxrg
457                DO  j = nysg, nyng
458                   DO  k = nzb, nzt+1
459                      u_av(k,j,i) = u_av(k,j,i) + u(k,j,i)
460                   ENDDO
461                ENDDO
462             ENDDO
463
464          CASE ( 'u*' )
465             DO  i = nxlg, nxrg
466                DO  j = nysg, nyng
467                   us_av(j,i) = us_av(j,i) + us(j,i)
468                ENDDO
469             ENDDO
470
471          CASE ( 'v' )
472             DO  i = nxlg, nxrg
473                DO  j = nysg, nyng
474                   DO  k = nzb, nzt+1
475                      v_av(k,j,i) = v_av(k,j,i) + v(k,j,i)
476                   ENDDO
477                ENDDO
478             ENDDO
479
480          CASE ( 'vpt' )
481             DO  i = nxlg, nxrg
482                DO  j = nysg, nyng
483                   DO  k = nzb, nzt+1
484                      vpt_av(k,j,i) = vpt_av(k,j,i) + vpt(k,j,i)
485                   ENDDO
486                ENDDO
487             ENDDO
488
489          CASE ( 'w' )
490             DO  i = nxlg, nxrg
491                DO  j = nysg, nyng
492                   DO  k = nzb, nzt+1
493                      w_av(k,j,i) = w_av(k,j,i) + w(k,j,i)
494                   ENDDO
495                ENDDO
496             ENDDO
497
498          CASE ( 'z0*' )
499             DO  i = nxlg, nxrg
500                DO  j = nysg, nyng
501                   z0_av(j,i) = z0_av(j,i) + z0(j,i)
502                ENDDO
503             ENDDO
504
505          CASE ( 'z0h*' )
506             DO  i = nxlg, nxrg
507                DO  j = nysg, nyng
508                   z0h_av(j,i) = z0h_av(j,i) + z0h(j,i)
509                ENDDO
510             ENDDO
511
512          CASE DEFAULT
513!
514!--          User-defined quantity
515             CALL user_3d_data_averaging( 'sum', doav(ii) )
516
517       END SELECT
518
519    ENDDO
520
521    CALL cpu_log (log_point(34),'sum_up_3d_data','stop','nobarrier')
522
523
524 END SUBROUTINE sum_up_3d_data
Note: See TracBrowser for help on using the repository browser.