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

Last change on this file since 77 was 77, checked in by raasch, 14 years ago

New:
---

particle reflection from vertical walls implemented, particle SGS model adjusted to walls

Wall functions for vertical walls now include diabatic conditions. New subroutines wall_fluxes, wall_fluxes_e. New 4D-array rif_wall.

new d3par-parameter netcdf_64bit_3d to switch on 64bit offset only for 3D files

new d3par-parameter dt_max to define the maximum value for the allowed timestep

new inipar-parameter loop_optimization to control the loop optimization method

new inipar-parameter pt_refrence. If given, this value is used as the reference that in buoyancy terms (otherwise, the instantaneous horizontally averaged temperature is used).

new user interface user_advec_particles

new initializing action "by_user" calls user_init_3d_model and allows the initial setting of all 3d arrays

topography height informations are stored on arrays zu_s_inner and zw_w_inner and output to the 2d/3d NetCDF files

samples added to the user interface which show how to add user-define time series quantities.

calculation/output of precipitation amount, precipitation rate and z0 (by setting "pra*", "prr*", "z0*" with data_output). The time interval on which the precipitation amount is defined is set by new d3par-parameter precipitation_amount_interval

unit 9 opened for debug output (file DEBUG_<pe#>)

Makefile, advec_particles, average_3d_data, buoyancy, calc_precipitation, check_open, check_parameters, data_output_2d, diffusion_e, diffusion_u, diffusion_v, diffusion_w, diffusivities, header, impact_of_latent_heat, init_particles, init_3d_model, modules, netcdf, parin, production_e, read_var_list, read_3d_binary, sum_up_3d_data, user_interface, write_var_list, write_3d_binary

New: wall_fluxes

Changed:


General revision of non-cyclic horizontal boundary conditions: radiation boundary conditions are now used instead of Neumann conditions at the outflow (calculation needs velocity values for t-dt, which are stored on new arrays u_m_l, u_m_r, etc.), calculation of mean outflow is not needed any more, volume flow control is added for the outflow boundary (currently only for the north boundary!!), additional gridpoints along x and y (uxrp, vynp) are not needed any more, routine "boundary_conds" now operates on timelevel t+dt and is not split in two parts (main, uvw_outflow) any more, Neumann boundary conditions at inflow/outflow in case of non-cyclic boundary conditions for all 2d-arrays that are handled by exchange_horiz_2d

The FFT-method for solving the Poisson-equation is now working with Neumann boundary conditions both at the bottom and the top. This requires adjustments of the tridiagonal coefficients and subtracting the horizontally averaged mean from the vertical velocity field.

+age_m in particle_type

Particles-package is now part of the default code ("-p particles" is not needed any more).

Move call of user_actions( 'after_integration' ) below increment of times
and counters. user_actions is now called for each statistic region and has as an argument the number of the respective region (sr)

d3par-parameter data_output_ts removed. Timeseries output for "profil" removed. Timeseries are now switched on by dt_dots. Timeseries data is collected in flow_statistics.

Initial velocities at nzb+1 are regarded for volume flow control in case they have been set zero before (to avoid small timesteps); see new internal parameters u/v_nzb_p1_for_vfc.

q is not allowed to become negative (prognostic_equations).

poisfft_init is only called if fft-solver is switched on (init_pegrid).

d3par-parameter moisture renamed to humidity.

Subversion global revision number is read from mrun and added to the run description header and to the run control (_rc) file.

vtk directives removed from main program.

The uitility routine interpret_config reads PALM environment variables from NAMELIST instead using the system call GETENV.

advec_u_pw, advec_u_up, advec_v_pw, advec_v_up, asselin_filter, check_parameters, coriolis, data_output_dvrp, data_output_ptseries, data_output_ts, data_output_2d, data_output_3d, diffusion_u, diffusion_v, exchange_horiz, exchange_horiz_2d, flow_statistics, header, init_grid, init_particles, init_pegrid, init_rankine, init_pt_anomaly, init_1d_model, init_3d_model, modules, palm, package_parin, parin, poisfft, poismg, prandtl_fluxes, pres, production_e, prognostic_equations, read_var_list, read_3d_binary, sor, swap_timelevel, time_integration, write_var_list, write_3d_binary

Errors:


Bugfix: preset of tendencies te_em, te_um, te_vm in init_1d_model

Bugfix in sample for reading user defined data from restart file (user_init)

Bugfix in setting diffusivities for cases with the outflow damping layer extending over more than one subdomain (init_3d_model)

Check for possible negative humidities in the initial humidity profile.

in Makefile, default suffixes removed from the suffix list to avoid calling of m2c in
# case of .mod files

Makefile
check_parameters, init_1d_model, init_3d_model, user_interface

  • Property svn:keywords set to Id
File size: 12.2 KB
Line 
1 SUBROUTINE sum_up_3d_data
2
3!------------------------------------------------------------------------------!
4! Actual revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: sum_up_3d_data.f90 77 2007-03-29 04:26:56Z raasch $
11!
12! 72 2007-03-19 08:20:46Z raasch
13! +sum-up of precipitation rate and roughness length (prr*, z0*)
14!
15! RCS Log replace by Id keyword, revision history cleaned up
16!
17! Revision 1.1  2006/02/23 12:55:23  raasch
18! Initial revision
19!
20!
21! Description:
22! ------------
23! Sum-up the values of 3d-arrays. The real averaging is later done in routine
24! average_3d_data.
25!------------------------------------------------------------------------------!
26
27    USE arrays_3d
28    USE averaging
29    USE cloud_parameters
30    USE control_parameters
31    USE cpulog
32    USE indices
33    USE interfaces
34    USE particle_attributes
35
36    IMPLICIT NONE
37
38    INTEGER ::  i, ii, j, k, n, psi
39
40    REAL    ::  mean_r, s_r3, s_r4
41
42
43    CALL cpu_log (log_point(34),'sum_up_3d_data','start')
44
45!
46!-- Allocate and initialize the summation arrays if called for the very first
47!-- time or the first time after average_3d_data has been called
48!-- (some or all of the arrays may have been already allocated
49!-- in read_3d_binary)
50    IF ( average_count_3d == 0 )  THEN
51
52       DO  ii = 1, doav_n
53
54          SELECT CASE ( TRIM( doav(ii) ) )
55
56             CASE ( 'e' )
57                IF ( .NOT. ALLOCATED( e_av ) )  THEN
58                   ALLOCATE( e_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
59                ENDIF
60                e_av = 0.0
61
62             CASE ( 'lwp*' )
63                IF ( .NOT. ALLOCATED( lwp_av ) )  THEN
64                   ALLOCATE( lwp_av(nys-1:nyn+1,nxl-1:nxr+1) )
65                ENDIF
66                lwp_av = 0.0
67
68             CASE ( 'p' )
69                IF ( .NOT. ALLOCATED( p_av ) )  THEN
70                   ALLOCATE( p_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
71                ENDIF
72                p_av = 0.0
73
74             CASE ( 'pc' )
75                IF ( .NOT. ALLOCATED( pc_av ) )  THEN
76                   ALLOCATE( pc_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
77                ENDIF
78                pc_av = 0.0
79
80             CASE ( 'pr' )
81                IF ( .NOT. ALLOCATED( pr_av ) )  THEN
82                   ALLOCATE( pr_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
83                ENDIF
84                pr_av = 0.0
85
86             CASE ( 'prr*' )
87                IF ( .NOT. ALLOCATED( precipitation_rate_av ) )  THEN
88                   ALLOCATE( precipitation_rate_av(nys-1:nyn+1,nxl-1:nxr+1) )
89                ENDIF
90                precipitation_rate_av = 0.0
91
92             CASE ( 'pt' )
93                IF ( .NOT. ALLOCATED( pt_av ) )  THEN
94                   ALLOCATE( pt_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
95                ENDIF
96                pt_av = 0.0
97
98             CASE ( 'q' )
99                IF ( .NOT. ALLOCATED( q_av ) )  THEN
100                   ALLOCATE( q_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
101                ENDIF
102                q_av = 0.0
103
104             CASE ( 'ql' )
105                IF ( .NOT. ALLOCATED( ql_av ) )  THEN
106                   ALLOCATE( ql_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
107                ENDIF
108                ql_av = 0.0
109
110             CASE ( 'ql_c' )
111                IF ( .NOT. ALLOCATED( ql_c_av ) )  THEN
112                   ALLOCATE( ql_c_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
113                ENDIF
114                ql_c_av = 0.0
115
116             CASE ( 'ql_v' )
117                IF ( .NOT. ALLOCATED( ql_v_av ) )  THEN
118                   ALLOCATE( ql_v_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
119                ENDIF
120                ql_v_av = 0.0
121
122             CASE ( 'ql_vp' )
123                IF ( .NOT. ALLOCATED( ql_vp_av ) )  THEN
124                   ALLOCATE( ql_vp_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
125                ENDIF
126                ql_vp_av = 0.0
127
128             CASE ( 'qv' )
129                IF ( .NOT. ALLOCATED( qv_av ) )  THEN
130                   ALLOCATE( qv_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
131                ENDIF
132                qv_av = 0.0
133
134             CASE ( 's' )
135                IF ( .NOT. ALLOCATED( s_av ) )  THEN
136                   ALLOCATE( s_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
137                ENDIF
138                s_av = 0.0
139
140             CASE ( 't*' )
141                IF ( .NOT. ALLOCATED( ts_av ) )  THEN
142                   ALLOCATE( ts_av(nys-1:nyn+1,nxl-1:nxr+1) )
143                ENDIF
144                ts_av = 0.0
145
146             CASE ( 'u' )
147                IF ( .NOT. ALLOCATED( u_av ) )  THEN
148                   ALLOCATE( u_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
149                ENDIF
150                u_av = 0.0
151
152             CASE ( 'u*' )
153                IF ( .NOT. ALLOCATED( us_av ) )  THEN
154                   ALLOCATE( us_av(nys-1:nyn+1,nxl-1:nxr+1) )
155                ENDIF
156                us_av = 0.0
157
158             CASE ( 'v' )
159                IF ( .NOT. ALLOCATED( v_av ) )  THEN
160                   ALLOCATE( v_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
161                ENDIF
162                v_av = 0.0
163
164             CASE ( 'vpt' )
165                IF ( .NOT. ALLOCATED( vpt_av ) )  THEN
166                   ALLOCATE( vpt_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
167                ENDIF
168                vpt_av = 0.0
169
170             CASE ( 'w' )
171                IF ( .NOT. ALLOCATED( w_av ) )  THEN
172                   ALLOCATE( w_av(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
173                ENDIF
174                w_av = 0.0
175
176             CASE ( 'z0*' )
177                IF ( .NOT. ALLOCATED( z0_av ) )  THEN
178                   ALLOCATE( z0_av(nys-1:nyn+1,nxl-1:nxr+1) )
179                ENDIF
180                z0_av = 0.0
181
182             CASE DEFAULT
183!
184!--             User-defined quantity
185                CALL user_3d_data_averaging( 'allocate', doav(ii) )
186
187          END SELECT
188
189       ENDDO
190
191    ENDIF
192
193!
194!-- Loop of all variables to be averaged.
195    DO  ii = 1, doav_n
196
197!
198!--    Store the array chosen on the temporary array.
199       SELECT CASE ( TRIM( doav(ii) ) )
200
201          CASE ( 'e' )
202             DO  i = nxl-1, nxr+1
203                DO  j = nys-1, nyn+1
204                   DO  k = nzb, nzt+1
205                      e_av(k,j,i) = e_av(k,j,i) + e(k,j,i)
206                   ENDDO
207                ENDDO
208             ENDDO
209
210          CASE ( 'lwp*' )
211             DO  i = nxl-1, nxr+1
212                DO  j = nys-1, nyn+1
213                   lwp_av(j,i) = lwp_av(j,i) + SUM( ql(nzb:nzt,j,i) * &
214                                                    dzw(1:nzt+1) )
215                ENDDO
216             ENDDO
217
218          CASE ( 'p' )
219             DO  i = nxl-1, nxr+1
220                DO  j = nys-1, nyn+1
221                   DO  k = nzb, nzt+1
222                      p_av(k,j,i) = p_av(k,j,i) + p(k,j,i)
223                   ENDDO
224                ENDDO
225             ENDDO
226
227          CASE ( 'pc' )
228             DO  i = nxl, nxr
229                DO  j = nys, nyn
230                   DO  k = nzb, nzt+1
231                      pc_av(k,j,i) = pc_av(k,j,i) + prt_count(k,j,i)
232                   ENDDO
233                ENDDO
234             ENDDO
235
236          CASE ( 'pr' )
237             DO  i = nxl, nxr
238                DO  j = nys, nyn
239                   DO  k = nzb, nzt+1
240                      psi = prt_start_index(k,j,i)
241                      s_r3 = 0.0
242                      s_r4 = 0.0
243                      DO  n = psi, psi+prt_count(k,j,i)-1
244                         s_r3 = s_r3 + particles(n)%radius**3
245                         s_r4 = s_r4 + particles(n)%radius**4
246                      ENDDO
247                      IF ( s_r3 /= 0.0 )  THEN
248                         mean_r = s_r4 / s_r3
249                      ELSE
250                         mean_r = 0.0
251                      ENDIF
252                      pr_av(k,j,i) = pr_av(k,j,i) + mean_r
253                   ENDDO
254                ENDDO
255             ENDDO
256
257          CASE ( 'pr*' )
258             DO  i = nxl-1, nxr+1
259                DO  j = nys-1, nyn+1
260                   precipitation_rate_av(j,i) = precipitation_rate_av(j,i) + &
261                                                precipitation_rate(j,i)
262                ENDDO
263             ENDDO
264
265          CASE ( 'pt' )
266             IF ( .NOT. cloud_physics ) THEN
267                DO  i = nxl-1, nxr+1
268                   DO  j = nys-1, nyn+1
269                      DO  k = nzb, nzt+1
270                         pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i)
271                      ENDDO
272                   ENDDO
273                ENDDO
274             ELSE
275                DO  i = nxl-1, nxr+1
276                   DO  j = nys-1, nyn+1
277                      DO  k = nzb, nzt+1
278                         pt_av(k,j,i) = pt_av(k,j,i) + pt(k,j,i) + l_d_cp * &
279                                                       pt_d_t(k) * ql(k,j,i)
280                      ENDDO
281                   ENDDO
282                ENDDO
283             ENDIF
284
285          CASE ( 'q' )
286             DO  i = nxl-1, nxr+1
287                DO  j = nys-1, nyn+1
288                   DO  k = nzb, nzt+1
289                      q_av(k,j,i) = q_av(k,j,i) + q(k,j,i)
290                   ENDDO
291                ENDDO
292             ENDDO
293             
294          CASE ( 'ql' )
295             DO  i = nxl-1, nxr+1
296                DO  j = nys-1, nyn+1
297                   DO  k = nzb, nzt+1
298                      ql_av(k,j,i) = ql_av(k,j,i) + ql(k,j,i)
299                   ENDDO
300                ENDDO
301             ENDDO
302
303          CASE ( 'ql_c' )
304             DO  i = nxl-1, nxr+1
305                DO  j = nys-1, nyn+1
306                   DO  k = nzb, nzt+1
307                      ql_c_av(k,j,i) = ql_c_av(k,j,i) + ql_c(k,j,i)
308                   ENDDO
309                ENDDO
310             ENDDO
311
312          CASE ( 'ql_v' )
313             DO  i = nxl-1, nxr+1
314                DO  j = nys-1, nyn+1
315                   DO  k = nzb, nzt+1
316                      ql_v_av(k,j,i) = ql_v_av(k,j,i) + ql_v(k,j,i)
317                   ENDDO
318                ENDDO
319             ENDDO
320
321          CASE ( 'ql_vp' )
322             DO  i = nxl-1, nxr+1
323                DO  j = nys-1, nyn+1
324                   DO  k = nzb, nzt+1
325                      ql_vp_av(k,j,i) = ql_vp_av(k,j,i) + ql_vp(k,j,i)
326                   ENDDO
327                ENDDO
328             ENDDO
329
330          CASE ( 'qv' )
331             DO  i = nxl-1, nxr+1
332                DO  j = nys-1, nyn+1
333                   DO  k = nzb, nzt+1
334                      qv_av(k,j,i) = qv_av(k,j,i) + q(k,j,i) - ql(k,j,i)
335                   ENDDO
336                ENDDO
337             ENDDO
338
339          CASE ( 's' )
340             DO  i = nxl-1, nxr+1
341                DO  j = nys-1, nyn+1
342                   DO  k = nzb, nzt+1
343                      s_av(k,j,i) = s_av(k,j,i) + q(k,j,i)
344                   ENDDO
345                ENDDO
346             ENDDO
347             
348          CASE ( 't*' )
349             DO  i = nxl-1, nxr+1
350                DO  j = nys-1, nyn+1
351                   ts_av(j,i) = ts_av(j,i) + ts(j,i)
352                ENDDO
353             ENDDO
354
355          CASE ( 'u' )
356             DO  i = nxl-1, nxr+1
357                DO  j = nys-1, nyn+1
358                   DO  k = nzb, nzt+1
359                      u_av(k,j,i) = u_av(k,j,i) + u(k,j,i)
360                   ENDDO
361                ENDDO
362             ENDDO
363
364          CASE ( 'u*' )
365             DO  i = nxl-1, nxr+1
366                DO  j = nys-1, nyn+1
367                   us_av(j,i) = us_av(j,i) + us(j,i)
368                ENDDO
369             ENDDO
370
371          CASE ( 'v' )
372             DO  i = nxl-1, nxr+1
373                DO  j = nys-1, nyn+1
374                   DO  k = nzb, nzt+1
375                      v_av(k,j,i) = v_av(k,j,i) + v(k,j,i)
376                   ENDDO
377                ENDDO
378             ENDDO
379
380          CASE ( 'vpt' )
381             DO  i = nxl-1, nxr+1
382                DO  j = nys-1, nyn+1
383                   DO  k = nzb, nzt+1
384                      vpt_av(k,j,i) = vpt_av(k,j,i) + vpt(k,j,i)
385                   ENDDO
386                ENDDO
387             ENDDO
388
389          CASE ( 'w' )
390             DO  i = nxl-1, nxr+1
391                DO  j = nys-1, nyn+1
392                   DO  k = nzb, nzt+1
393                      w_av(k,j,i) = w_av(k,j,i) + w(k,j,i)
394                   ENDDO
395                ENDDO
396             ENDDO
397
398          CASE ( 'z0*' )
399             DO  i = nxl-1, nxr+1
400                DO  j = nys-1, nyn+1
401                   z0_av(j,i) = z0_av(j,i) + z0(j,i)
402                ENDDO
403             ENDDO
404
405          CASE DEFAULT
406!
407!--          User-defined quantity
408             CALL user_3d_data_averaging( 'sum', doav(ii) )
409
410       END SELECT
411
412    ENDDO
413
414    CALL cpu_log (log_point(34),'sum_up_3d_data','stop','nobarrier')
415
416
417 END SUBROUTINE sum_up_3d_data
Note: See TracBrowser for help on using the repository browser.