source: palm/tags/release-3.3/SOURCE/data_output_ptseries.f90 @ 4418

Last change on this file since 4418 was 77, checked in by raasch, 17 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: 10.8 KB
Line 
1 SUBROUTINE data_output_ptseries
2
3!------------------------------------------------------------------------------!
4! Actual revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: data_output_ptseries.f90 77 2007-03-29 04:26:56Z raasch $
11!
12! 60 2007-03-11 11:50:04Z raasch
13! Particles-package is now part of the default code.
14!
15! RCS Log replace by Id keyword, revision history cleaned up
16!
17! Revision 1.2  2006/08/22 13:51:13  raasch
18! Seperate output for particle groups
19!
20! Revision 1.1  2006/08/04 14:24:18  raasch
21! Initial revision
22!
23!
24! Description:
25! ------------
26! Output of particle data timeseries in NetCDF format.
27!------------------------------------------------------------------------------!
28
29    USE control_parameters
30    USE cpulog
31    USE indices
32    USE interfaces
33    USE netcdf_control
34    USE particle_attributes
35    USE pegrid
36
37    IMPLICIT NONE
38
39
40    INTEGER ::  i, inum, j, n
41
42    REAL, DIMENSION(0:number_of_particle_groups,30) ::  pts_value, pts_value_l
43
44
45
46    CALL cpu_log( log_point(36), 'data_output_ptseries', 'start' )
47
48    IF ( myid == 0  .AND.  netcdf_output )  THEN
49!
50!--    Open file for time series output in NetCDF format
51       dopts_time_count = dopts_time_count + 1
52       CALL check_open( 109 )
53#if defined( __netcdf )
54!
55!--    Update the particle time series time axis
56       nc_stat = NF90_PUT_VAR( id_set_pts, id_var_time_pts,    &
57                               (/ simulated_time /),           &
58                               start = (/ dopts_time_count /), count = (/ 1 /) )
59       IF (nc_stat /= NF90_NOERR)  CALL handle_netcdf_error( 391 )
60#endif
61
62    ENDIF
63
64    pts_value_l = 0.0
65
66!
67!-- Calculate or collect the particle time series quantities for all particles
68!-- and seperately for each particle group (if there is more than one group)
69    DO  n = 1, number_of_particles
70
71       pts_value_l(0,1)  = number_of_particles            ! total # of particles
72       pts_value_l(0,2)  = pts_value_l(0,2) + &
73                           ( particles(n)%x - particles(n)%origin_x )  ! mean x
74       pts_value_l(0,3)  = pts_value_l(0,3) + &
75                           ( particles(n)%y - particles(n)%origin_y )  ! mean y
76       pts_value_l(0,4)  = pts_value_l(0,4) + &
77                           ( particles(n)%z - particles(n)%origin_z )  ! mean z
78       pts_value_l(0,5)  = pts_value_l(0,5) + particles(n)%z ! mean z (absolute)
79       pts_value_l(0,6)  = pts_value_l(0,6) + particles(n)%speed_x     ! mean u
80       pts_value_l(0,7)  = pts_value_l(0,7) + particles(n)%speed_y     ! mean v
81       pts_value_l(0,8)  = pts_value_l(0,8) + particles(n)%speed_z     ! mean w
82       pts_value_l(0,9)  = pts_value_l(0,9) + &
83                                            particles(n)%speed_x_sgs ! mean sgsu
84       pts_value_l(0,10) = pts_value_l(0,10) + &
85                                            particles(n)%speed_y_sgs ! mean sgsv
86       pts_value_l(0,11) = pts_value_l(0,11) + &
87                                            particles(n)%speed_z_sgs ! mean sgsw
88       IF ( particles(n)%speed_z > 0.0 )  THEN
89          pts_value_l(0,12) = pts_value_l(0,12) + 1.0  ! # of upward moving prts
90          pts_value_l(0,13) = pts_value_l(0,13) + &
91                                            particles(n)%speed_z ! mean w upw.
92       ELSE
93          pts_value_l(0,14) = pts_value_l(0,14) + &
94                                            particles(n)%speed_z ! mean w down
95       ENDIF
96       pts_value_l(0,15) = number_of_particles
97       pts_value_l(0,16) = number_of_particles
98
99!
100!--    Repeat the same for the respective particle group
101       IF ( number_of_particle_groups > 1 )  THEN
102          j = particles(n)%group
103
104          pts_value_l(j,1)  = pts_value_l(j,1)  + 1
105          pts_value_l(j,2)  = pts_value_l(j,2)  + &
106                              ( particles(n)%x - particles(n)%origin_x )
107          pts_value_l(j,3)  = pts_value_l(j,3)  + &
108                              ( particles(n)%y - particles(n)%origin_y )
109          pts_value_l(j,4)  = pts_value_l(j,4)  + &
110                              ( particles(n)%z - particles(n)%origin_z )
111          pts_value_l(j,5)  = pts_value_l(j,5)  + particles(n)%z
112          pts_value_l(j,6)  = pts_value_l(j,6)  + particles(n)%speed_x
113          pts_value_l(j,7)  = pts_value_l(j,7)  + particles(n)%speed_y
114          pts_value_l(j,8)  = pts_value_l(j,8)  + particles(n)%speed_z
115          pts_value_l(j,9)  = pts_value_l(j,9)  + particles(n)%speed_x_sgs
116          pts_value_l(j,10) = pts_value_l(j,10) + particles(n)%speed_y_sgs
117          pts_value_l(j,11) = pts_value_l(j,11) + particles(n)%speed_z_sgs
118          IF ( particles(n)%speed_z > 0.0 )  THEN
119             pts_value_l(j,12) = pts_value_l(j,12) + 1.0
120             pts_value_l(j,13) = pts_value_l(j,13) + particles(n)%speed_z
121          ELSE
122             pts_value_l(j,14) = pts_value_l(j,14) + particles(n)%speed_z
123          ENDIF
124          pts_value_l(j,15) = pts_value_l(j,15) + 1.0
125          pts_value_l(j,16) = pts_value_l(j,16) + 1.0
126
127       ENDIF
128
129    ENDDO
130
131#if defined( __parallel )
132!
133!-- Sum values of the subdomains
134    inum = number_of_particle_groups + 1
135
136    CALL MPI_ALLREDUCE( pts_value_l(0,1), pts_value(0,1), 14*inum, MPI_REAL, &
137                        MPI_SUM, comm2d, ierr )
138    CALL MPI_ALLREDUCE( pts_value_l(0,15), pts_value(0,15), inum, MPI_REAL, &
139                        MPI_MAX, comm2d, ierr )
140    CALL MPI_ALLREDUCE( pts_value_l(0,16), pts_value(0,16), inum, MPI_REAL, &
141                        MPI_MIN, comm2d, ierr )
142#else
143    pts_value(:,1:16) = pts_value_l(:,1:16)
144#endif
145
146!
147!-- Normalize the above calculated quantities with the total number of
148!-- particles
149    IF ( number_of_particle_groups > 1 )  THEN
150       inum = number_of_particle_groups
151    ELSE
152       inum = 0
153    ENDIF
154
155    DO  j = 0, inum
156
157       IF ( pts_value(j,1) > 0.0 )  THEN
158
159          pts_value(j,2:14) = pts_value(j,2:14) / pts_value(j,1)
160          IF ( pts_value(j,12) > 0.0  .AND.  pts_value(j,12) < 1.0 )  THEN
161             pts_value(j,13) = pts_value(j,13) / pts_value(j,12)
162             pts_value(j,14) = pts_value(j,14) / ( 1.0 - pts_value(j,12) )
163          ELSEIF ( pts_value(j,12) == 0.0 )  THEN
164             pts_value(j,13) = -1.0
165          ELSE
166             pts_value(j,14) = -1.0
167          ENDIF
168
169       ENDIF
170
171    ENDDO
172
173!
174!-- Calculate higher order moments of particle time series quantities,
175!-- seperately for each particle group (if there is more than one group)
176    DO  n = 1, number_of_particles
177
178       pts_value_l(0,17) = pts_value_l(0,17) + ( particles(n)%x - &
179                           particles(n)%origin_x - pts_value(0,2) )**2 ! x*2
180       pts_value_l(0,18) = pts_value_l(0,18) + ( particles(n)%y - &
181                           particles(n)%origin_y - pts_value(0,3) )**2 ! y*2
182       pts_value_l(0,19) = pts_value_l(0,19) + ( particles(n)%z - &
183                           particles(n)%origin_z - pts_value(0,4) )**2 ! z*2
184       pts_value_l(0,20) = pts_value_l(0,20) + ( particles(n)%speed_x - &
185                                               pts_value(0,6) )**2     ! u*2
186       pts_value_l(0,21) = pts_value_l(0,21) + ( particles(n)%speed_y - &
187                                               pts_value(0,7) )**2     ! v*2
188       pts_value_l(0,22) = pts_value_l(0,22) + ( particles(n)%speed_z - &
189                                               pts_value(0,8) )**2     ! w*2
190       pts_value_l(0,23) = pts_value_l(0,23) + ( particles(n)%speed_x_sgs - &
191                                               pts_value(0,9) )**2     ! u"2
192       pts_value_l(0,24) = pts_value_l(0,24) + ( particles(n)%speed_y_sgs - &
193                                               pts_value(0,10) )**2    ! v"2
194       pts_value_l(0,25) = pts_value_l(0,25) + ( particles(n)%speed_z_sgs - &
195                                               pts_value(0,11) )**2    ! w"2
196!
197!--    Repeat the same for the respective particle group
198       IF ( number_of_particle_groups > 1 )  THEN
199          j = particles(n)%group
200
201          pts_value_l(j,17) = pts_value_l(j,17) + ( particles(n)%x - &
202                              particles(n)%origin_x - pts_value(j,2) )**2
203          pts_value_l(j,18) = pts_value_l(j,18) + ( particles(n)%y - &
204                              particles(n)%origin_y - pts_value(j,3) )**2
205          pts_value_l(j,19) = pts_value_l(j,19) + ( particles(n)%z - &
206                              particles(n)%origin_z - pts_value(j,4) )**2
207          pts_value_l(j,20) = pts_value_l(j,20) + ( particles(n)%speed_x - &
208                                                  pts_value(j,6) )**2
209          pts_value_l(j,21) = pts_value_l(j,21) + ( particles(n)%speed_y - &
210                                                  pts_value(j,7) )**2
211          pts_value_l(j,22) = pts_value_l(j,22) + ( particles(n)%speed_z - &
212                                                  pts_value(j,8) )**2
213          pts_value_l(j,23) = pts_value_l(j,23) + ( particles(n)%speed_x_sgs - &
214                                                  pts_value(j,9) )**2
215          pts_value_l(j,24) = pts_value_l(j,24) + ( particles(n)%speed_y_sgs - &
216                                                  pts_value(j,10) )**2
217          pts_value_l(j,25) = pts_value_l(j,25) + ( particles(n)%speed_z_sgs - &
218                                                  pts_value(j,11) )**2
219       ENDIF
220
221    ENDDO
222
223    pts_value_l(0,26) = ( number_of_particles - pts_value(0,1) / numprocs )**2
224                                                 ! variance of particle numbers
225    IF ( number_of_particle_groups > 1 )  THEN
226       DO  j = 1, number_of_particle_groups
227          pts_value_l(j,26) = ( pts_value_l(j,1) - &
228                                pts_value(j,1) / numprocs )**2
229       ENDDO
230    ENDIF
231
232#if defined( __parallel )
233!
234!-- Sum values of the subdomains
235    inum = number_of_particle_groups + 1
236
237    CALL MPI_ALLREDUCE( pts_value_l(0,17), pts_value(0,17), inum*10, MPI_REAL, &
238                        MPI_SUM, comm2d, ierr )
239#else
240    pts_value(:,17:26) = pts_value_l(:,17:26)
241#endif
242
243!
244!-- Normalize the above calculated quantities with the total number of
245!-- particles
246    IF ( number_of_particle_groups > 1 )  THEN
247       inum = number_of_particle_groups
248    ELSE
249       inum = 0
250    ENDIF
251
252    DO  j = 0, inum
253
254       IF ( pts_value(j,1) > 0.0 )  THEN
255          pts_value(j,17:25) = pts_value(j,17:25) / pts_value(j,1)
256       ENDIF
257       pts_value(j,26) = pts_value(j,26) / numprocs
258
259    ENDDO
260
261#if defined( __netcdf )
262!
263!-- Output particle time series quantities in NetCDF format
264    IF ( myid == 0  .AND.  netcdf_output )  THEN
265       DO  j = 0, inum
266          DO  i = 1, dopts_num
267             nc_stat = NF90_PUT_VAR( id_set_pts, id_var_dopts(i,j),  &
268                                     (/ pts_value(j,i) /),           &
269                                     start = (/ dopts_time_count /), &
270                                     count = (/ 1 /) )
271             IF ( nc_stat /= NF90_NOERR )  CALL handle_netcdf_error( 392 )
272          ENDDO
273       ENDDO
274    ENDIF
275#endif
276
277    CALL cpu_log( log_point(36), 'data_output_ptseries','stop', 'nobarrier' )
278
279 END SUBROUTINE data_output_ptseries
Note: See TracBrowser for help on using the repository browser.