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

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

last commit documented

  • Property svn:keywords set to Id
File size: 15.2 KB
RevLine 
[1]1 SUBROUTINE sum_up_3d_data
2
3!------------------------------------------------------------------------------!
[484]4! Current revisions:
[1]5! -----------------
6!
[979]7!
[1]8! Former revisions:
9! -----------------
[3]10! $Id: sum_up_3d_data.f90 979 2012-08-09 08:50:11Z raasch $
[77]11!
[979]12! 978 2012-08-09 08:28:32Z fricke
13! +z0h*
14!
[791]15! 790 2011-11-29 03:11:20Z raasch
16! bugfix: calculation of 'pr' must depend on the particle weighting factor
17!
[772]18! 771 2011-10-27 10:56:21Z heinze
19! +lpt_av
20!
[668]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!
[482]24! 402 2009-10-21 11:59:41Z maronga
25! Bugfix in calculation of shf*_av, qsws*_av
26!
[392]27! 2009-08-25 08:35:52Z maronga
28! +shf*, qsws*
29!
[98]30! 96 2007-06-04 08:07:41Z raasch
31! +sum-up of density and salinity
32!
[77]33! 72 2007-03-19 08:20:46Z raasch
34! +sum-up of precipitation rate and roughness length (prr*, z0*)
35!
[3]36! RCS Log replace by Id keyword, revision history cleaned up
37!
[1]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
[667]79                   ALLOCATE( e_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]80                ENDIF
81                e_av = 0.0
82
[771]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
[1]89             CASE ( 'lwp*' )
90                IF ( .NOT. ALLOCATED( lwp_av ) )  THEN
[667]91                   ALLOCATE( lwp_av(nysg:nyng,nxlg:nxrg) )
[1]92                ENDIF
93                lwp_av = 0.0
94
95             CASE ( 'p' )
96                IF ( .NOT. ALLOCATED( p_av ) )  THEN
[667]97                   ALLOCATE( p_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]98                ENDIF
99                p_av = 0.0
100
101             CASE ( 'pc' )
102                IF ( .NOT. ALLOCATED( pc_av ) )  THEN
[667]103                   ALLOCATE( pc_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]104                ENDIF
105                pc_av = 0.0
106
107             CASE ( 'pr' )
108                IF ( .NOT. ALLOCATED( pr_av ) )  THEN
[667]109                   ALLOCATE( pr_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]110                ENDIF
111                pr_av = 0.0
112
[72]113             CASE ( 'prr*' )
114                IF ( .NOT. ALLOCATED( precipitation_rate_av ) )  THEN
[667]115                   ALLOCATE( precipitation_rate_av(nysg:nyng,nxlg:nxrg) )
[72]116                ENDIF
117                precipitation_rate_av = 0.0
118
[1]119             CASE ( 'pt' )
120                IF ( .NOT. ALLOCATED( pt_av ) )  THEN
[667]121                   ALLOCATE( pt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]122                ENDIF
123                pt_av = 0.0
124
125             CASE ( 'q' )
126                IF ( .NOT. ALLOCATED( q_av ) )  THEN
[667]127                   ALLOCATE( q_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]128                ENDIF
129                q_av = 0.0
130
131             CASE ( 'ql' )
132                IF ( .NOT. ALLOCATED( ql_av ) )  THEN
[667]133                   ALLOCATE( ql_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]134                ENDIF
135                ql_av = 0.0
136
137             CASE ( 'ql_c' )
138                IF ( .NOT. ALLOCATED( ql_c_av ) )  THEN
[667]139                   ALLOCATE( ql_c_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]140                ENDIF
141                ql_c_av = 0.0
142
143             CASE ( 'ql_v' )
144                IF ( .NOT. ALLOCATED( ql_v_av ) )  THEN
[667]145                   ALLOCATE( ql_v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]146                ENDIF
147                ql_v_av = 0.0
148
149             CASE ( 'ql_vp' )
150                IF ( .NOT. ALLOCATED( ql_vp_av ) )  THEN
[667]151                   ALLOCATE( ql_vp_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]152                ENDIF
153                ql_vp_av = 0.0
154
[354]155             CASE ( 'qsws*' )
156                IF ( .NOT. ALLOCATED( qsws_av ) )  THEN
[667]157                   ALLOCATE( qsws_av(nysg:nyng,nxlg:nxrg) )
[354]158                ENDIF
159                qsws_av = 0.0
160
[1]161             CASE ( 'qv' )
162                IF ( .NOT. ALLOCATED( qv_av ) )  THEN
[667]163                   ALLOCATE( qv_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]164                ENDIF
165                qv_av = 0.0
166
[96]167             CASE ( 'rho' )
168                IF ( .NOT. ALLOCATED( rho_av ) )  THEN
[667]169                   ALLOCATE( rho_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[96]170                ENDIF
171                rho_av = 0.0
172
[1]173             CASE ( 's' )
174                IF ( .NOT. ALLOCATED( s_av ) )  THEN
[667]175                   ALLOCATE( s_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]176                ENDIF
177                s_av = 0.0
178
[96]179             CASE ( 'sa' )
180                IF ( .NOT. ALLOCATED( sa_av ) )  THEN
[667]181                   ALLOCATE( sa_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[96]182                ENDIF
183                sa_av = 0.0
184
[354]185             CASE ( 'shf*' )
186                IF ( .NOT. ALLOCATED( shf_av ) )  THEN
[667]187                   ALLOCATE( shf_av(nysg:nyng,nxlg:nxrg) )
[354]188                ENDIF
189                shf_av = 0.0
190
[1]191             CASE ( 't*' )
192                IF ( .NOT. ALLOCATED( ts_av ) )  THEN
[667]193                   ALLOCATE( ts_av(nysg:nyng,nxlg:nxrg) )
[1]194                ENDIF
195                ts_av = 0.0
196
197             CASE ( 'u' )
198                IF ( .NOT. ALLOCATED( u_av ) )  THEN
[667]199                   ALLOCATE( u_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]200                ENDIF
201                u_av = 0.0
202
203             CASE ( 'u*' )
204                IF ( .NOT. ALLOCATED( us_av ) )  THEN
[667]205                   ALLOCATE( us_av(nysg:nyng,nxlg:nxrg) )
[1]206                ENDIF
207                us_av = 0.0
208
209             CASE ( 'v' )
210                IF ( .NOT. ALLOCATED( v_av ) )  THEN
[667]211                   ALLOCATE( v_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]212                ENDIF
213                v_av = 0.0
214
215             CASE ( 'vpt' )
216                IF ( .NOT. ALLOCATED( vpt_av ) )  THEN
[667]217                   ALLOCATE( vpt_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]218                ENDIF
219                vpt_av = 0.0
220
221             CASE ( 'w' )
222                IF ( .NOT. ALLOCATED( w_av ) )  THEN
[667]223                   ALLOCATE( w_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1]224                ENDIF
225                w_av = 0.0
226
[72]227             CASE ( 'z0*' )
228                IF ( .NOT. ALLOCATED( z0_av ) )  THEN
[667]229                   ALLOCATE( z0_av(nysg:nyng,nxlg:nxrg) )
[72]230                ENDIF
231                z0_av = 0.0
232
[978]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
[1]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' )
[667]259             DO  i = nxlg, nxrg
260                DO  j = nysg, nyng
[1]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
[771]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
[1]276          CASE ( 'lwp*' )
[667]277             DO  i = nxlg, nxrg
278                DO  j = nysg, nyng
[1]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' )
[667]285             DO  i = nxlg, nxrg
286                DO  j = nysg, nyng
[1]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
[790]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
[1]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
[72]325          CASE ( 'pr*' )
[667]326             DO  i = nxlg, nxrg
327                DO  j = nysg, nyng
[72]328                   precipitation_rate_av(j,i) = precipitation_rate_av(j,i) + &
329                                                precipitation_rate(j,i)
330                ENDDO
331             ENDDO
332
[1]333          CASE ( 'pt' )
334             IF ( .NOT. cloud_physics ) THEN
[667]335             DO  i = nxlg, nxrg
336                DO  j = nysg, nyng
337                   DO  k = nzb, nzt+1
[1]338                         pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i)
339                      ENDDO
340                   ENDDO
341                ENDDO
342             ELSE
[667]343             DO  i = nxlg, nxrg
344                DO  j = nysg, nyng
345                   DO  k = nzb, nzt+1
[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' )
[667]354             DO  i = nxlg, nxrg
355                DO  j = nysg, nyng
[1]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
[402]361
[1]362          CASE ( 'ql' )
[667]363             DO  i = nxlg, nxrg
364                DO  j = nysg, nyng
[1]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' )
[667]372             DO  i = nxlg, nxrg
373                DO  j = nysg, nyng
[1]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' )
[667]381             DO  i = nxlg, nxrg
382                DO  j = nysg, nyng
[1]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' )
[667]390             DO  i = nxlg, nxrg
391                DO  j = nysg, nyng
[1]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
[402]398          CASE ( 'qsws*' )
[667]399             DO  i = nxlg, nxrg
400                DO  j = nysg, nyng
[402]401                   qsws_av(j,i) = qsws_av(j,i) + qsws(j,i)
402                ENDDO
403             ENDDO
404
[1]405          CASE ( 'qv' )
[667]406             DO  i = nxlg, nxrg
407                DO  j = nysg, nyng
[1]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
[96]414          CASE ( 'rho' )
[667]415             DO  i = nxlg, nxrg
416                DO  j = nysg, nyng
[96]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
[402]422
[1]423          CASE ( 's' )
[667]424             DO  i = nxlg, nxrg
425                DO  j = nysg, nyng
[1]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
[402]431
[96]432          CASE ( 'sa' )
[667]433             DO  i = nxlg, nxrg
434                DO  j = nysg, nyng
[96]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
[402]440
441          CASE ( 'shf*' )
[667]442             DO  i = nxlg, nxrg
443                DO  j = nysg, nyng
[402]444                   shf_av(j,i) = shf_av(j,i) + shf(j,i)
445                ENDDO
446             ENDDO
447
[1]448          CASE ( 't*' )
[667]449             DO  i = nxlg, nxrg
450                DO  j = nysg, nyng
[1]451                   ts_av(j,i) = ts_av(j,i) + ts(j,i)
452                ENDDO
453             ENDDO
454
455          CASE ( 'u' )
[667]456             DO  i = nxlg, nxrg
457                DO  j = nysg, nyng
[1]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*' )
[667]465             DO  i = nxlg, nxrg
466                DO  j = nysg, nyng
[1]467                   us_av(j,i) = us_av(j,i) + us(j,i)
468                ENDDO
469             ENDDO
470
471          CASE ( 'v' )
[667]472             DO  i = nxlg, nxrg
473                DO  j = nysg, nyng
[1]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' )
[667]481             DO  i = nxlg, nxrg
482                DO  j = nysg, nyng
[1]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' )
[667]490             DO  i = nxlg, nxrg
491                DO  j = nysg, nyng
[1]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
[72]498          CASE ( 'z0*' )
[667]499             DO  i = nxlg, nxrg
500                DO  j = nysg, nyng
[72]501                   z0_av(j,i) = z0_av(j,i) + z0(j,i)
502                ENDDO
503             ENDDO
504
[978]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
[1]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.