source: palm/trunk/SOURCE/data_output_mask.f90 @ 583

Last change on this file since 583 was 565, checked in by helmke, 14 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 15.1 KB
RevLine 
[298]1 SUBROUTINE data_output_mask( av )
2
3!------------------------------------------------------------------------------!
[484]4! Current revisions:
[298]5! -----------------
[565]6!
[298]7!
8! Former revisions:
9! -----------------
10! $Id: data_output_mask.f90 565 2010-09-30 13:34:53Z heinze $
11!
[565]12! 564 2010-09-30 13:18:59Z helmke
13! start number of mask output files changed to 201, netcdf message identifiers
14! of masked output changed, palm message identifiers of masked output changed
15!
[494]16! 493 2010-03-01 08:30:24Z raasch
17! netcdf_format_mask* and format_parallel_io replaced by netcdf_data_format
18!
[482]19! 475 2010-02-04 02:26:16Z raasch
20! Bugfix in serial branch: arguments from array local_pf removed in N90_PUT_VAR
21!
[449]22! 410 2009-12-04 17:05:40Z letzel
23! Initial version
[298]24!
25! Description:
26! ------------
27! Masked data output in NetCDF format for current mask (current value of mid).
28!------------------------------------------------------------------------------!
29
30#if defined( __netcdf )
31    USE arrays_3d
32    USE averaging
33    USE cloud_parameters
34    USE control_parameters
35    USE cpulog
36    USE grid_variables
37    USE indices
38    USE interfaces
39    USE netcdf
40    USE netcdf_control
41    USE particle_attributes
42    USE pegrid
43
44    IMPLICIT NONE
45
46    INTEGER ::  av, ngp, file_id, i, if, is, j, k, l, n, psi, s, sender, &
47                ind(6)
48    LOGICAL ::  found, resorted
49    REAL    ::  mean_r, s_r3, s_r4
50    REAL, DIMENSION(:,:,:), ALLOCATABLE ::  local_pf
51#if defined( __parallel )
52    REAL, DIMENSION(:,:,:), ALLOCATABLE ::  total_pf
53#endif
54    REAL, DIMENSION(:,:,:), POINTER ::  to_be_resorted
55
56!
57!-- Return, if nothing to output
58    IF ( domask_no(mid,av) == 0 )  RETURN
59
60    CALL cpu_log (log_point(49),'data_output_mask','start')
61
[493]62!
[298]63!-- Open output file.
[493]64    IF ( netcdf_output  .AND.  ( myid == 0  .OR.  netcdf_data_format > 2 ) ) &
65    THEN
[564]66       CALL check_open( 200+mid+av*max_masks )
[409]67    ENDIF 
[298]68
69!
70!-- Allocate total and local output arrays.
71#if defined( __parallel )
72    IF ( myid == 0 )  THEN
73       ALLOCATE( total_pf(mask_size(mid,1),mask_size(mid,2),mask_size(mid,3)) )
74    ENDIF
75#endif
76    ALLOCATE( local_pf(mask_size_l(mid,1),mask_size_l(mid,2), &
77                       mask_size_l(mid,3)) )
78
79!
80!-- Update the NetCDF time axis.
81    domask_time_count(mid,av) = domask_time_count(mid,av) + 1
[493]82    IF ( netcdf_output  .AND.  ( myid == 0  .OR.  netcdf_data_format > 2 ) ) &
83    THEN
[298]84       nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_time_mask(mid,av), &
85                               (/ simulated_time /),                          &
86                               start = (/ domask_time_count(mid,av) /),       &
87                               count = (/ 1 /) )
[564]88       CALL handle_netcdf_error( 'data_output_mask', 460 )
[298]89    ENDIF
90
91!
92!-- Loop over all variables to be written.
93    if = 1
94
95    DO  WHILE ( domask(mid,av,if)(1:1) /= ' ' )
96!
97!--    Reallocate local_pf on PE 0 since its shape changes during MPI exchange
[493]98       IF ( netcdf_data_format < 3   .AND.  myid == 0  .AND.  if > 1 )  THEN
[298]99          DEALLOCATE( local_pf )
100          ALLOCATE( local_pf(mask_size_l(mid,1),mask_size_l(mid,2), &
101                             mask_size_l(mid,3)) )
102       ENDIF
103!
104!--    Store the variable chosen.
105       resorted = .FALSE.
106       SELECT CASE ( TRIM( domask(mid,av,if) ) )
107
108          CASE ( 'e' )
109             IF ( av == 0 )  THEN
110                to_be_resorted => e
111             ELSE
112                to_be_resorted => e_av
113             ENDIF
114
115          CASE ( 'p' )
116             IF ( av == 0 )  THEN
117                to_be_resorted => p
118             ELSE
119                to_be_resorted => p_av
120             ENDIF
121
122          CASE ( 'pc' )  ! particle concentration (requires ghostpoint exchange)
123             IF ( av == 0 )  THEN
124                tend = prt_count
125                CALL exchange_horiz( tend )
126                DO  i = 1, mask_size_l(mid,1)
127                   DO  j = 1, mask_size_l(mid,2)
128                      DO  k = 1, mask_size_l(mid,3)
129                         local_pf(i,j,k) =  tend(mask_k(mid,k), &
130                                   mask_j(mid,j),mask_i(mid,i))
131                      ENDDO
132                   ENDDO
133                ENDDO
134                resorted = .TRUE.
135             ELSE
136                CALL exchange_horiz( pc_av )
137                to_be_resorted => pc_av
138             ENDIF
139
140          CASE ( 'pr' )  ! mean particle radius
141             IF ( av == 0 )  THEN
142                DO  i = nxl, nxr
143                   DO  j = nys, nyn
144                      DO  k = nzb, nzt+1
145                         psi = prt_start_index(k,j,i)
146                         s_r3 = 0.0
147                         s_r4 = 0.0
148                         DO  n = psi, psi+prt_count(k,j,i)-1
149                            s_r3 = s_r3 + particles(n)%radius**3
150                            s_r4 = s_r4 + particles(n)%radius**4
151                         ENDDO
152                         IF ( s_r3 /= 0.0 )  THEN
153                            mean_r = s_r4 / s_r3
154                         ELSE
155                            mean_r = 0.0
156                         ENDIF
157                         tend(k,j,i) = mean_r
158                      ENDDO
159                   ENDDO
160                ENDDO
161                CALL exchange_horiz( tend )
162                DO  i = 1, mask_size_l(mid,1)
163                   DO  j = 1, mask_size_l(mid,2)
164                      DO  k = 1, mask_size_l(mid,3)
165                         local_pf(i,j,k) =  tend(mask_k(mid,k), &
166                                   mask_j(mid,j),mask_i(mid,i))
167                      ENDDO
168                   ENDDO
169                ENDDO
170                resorted = .TRUE.
171             ELSE
172                CALL exchange_horiz( pr_av )
173                to_be_resorted => pr_av
174             ENDIF
175
176          CASE ( 'pt' )
177             IF ( av == 0 )  THEN
178                IF ( .NOT. cloud_physics ) THEN
179                   to_be_resorted => pt
180                ELSE
181                   DO  i = 1, mask_size_l(mid,1)
182                      DO  j = 1, mask_size_l(mid,2)
183                         DO  k = 1, mask_size_l(mid,3)
184                            local_pf(i,j,k) =  &
185                                 pt(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i)) &
186                                 + l_d_cp * pt_d_t(mask_k(mid,k)) * &
187                                   ql(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i))
188                         ENDDO
189                      ENDDO
190                   ENDDO
191                   resorted = .TRUE.
192                ENDIF
193             ELSE
194                to_be_resorted => pt_av
195             ENDIF
196
197          CASE ( 'q' )
198             IF ( av == 0 )  THEN
199                to_be_resorted => q
200             ELSE
201                to_be_resorted => q_av
202             ENDIF
203             
204          CASE ( 'ql' )
205             IF ( av == 0 )  THEN
206                to_be_resorted => ql
207             ELSE
208                to_be_resorted => ql_av
209             ENDIF
210
211          CASE ( 'ql_c' )
212             IF ( av == 0 )  THEN
213                to_be_resorted => ql_c
214             ELSE
215                to_be_resorted => ql_c_av
216             ENDIF
217
218          CASE ( 'ql_v' )
219             IF ( av == 0 )  THEN
220                to_be_resorted => ql_v
221             ELSE
222                to_be_resorted => ql_v_av
223             ENDIF
224
225          CASE ( 'ql_vp' )
226             IF ( av == 0 )  THEN
227                to_be_resorted => ql_vp
228             ELSE
229                to_be_resorted => ql_vp_av
230             ENDIF
231
232          CASE ( 'qv' )
233             IF ( av == 0 )  THEN
234                DO  i = 1, mask_size_l(mid,1)
235                   DO  j = 1, mask_size_l(mid,2)
236                      DO  k = 1, mask_size_l(mid,3)
237                         local_pf(i,j,k) =  &
238                              q(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i)) -  &
239                              ql(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i))
240                      ENDDO
241                   ENDDO
242                ENDDO
243                resorted = .TRUE.
244             ELSE
245                to_be_resorted => qv_av
246             ENDIF
247
248          CASE ( 'rho' )
249             IF ( av == 0 )  THEN
250                to_be_resorted => rho
251             ELSE
252                to_be_resorted => rho_av
253             ENDIF
254             
255          CASE ( 's' )
256             IF ( av == 0 )  THEN
257                to_be_resorted => q
258             ELSE
[356]259                to_be_resorted => s_av
[298]260             ENDIF
261             
262          CASE ( 'sa' )
263             IF ( av == 0 )  THEN
264                to_be_resorted => sa
265             ELSE
266                to_be_resorted => sa_av
267             ENDIF
268             
269          CASE ( 'u' )
270             IF ( av == 0 )  THEN
271                to_be_resorted => u
272             ELSE
273                to_be_resorted => u_av
274             ENDIF
275
276          CASE ( 'v' )
277             IF ( av == 0 )  THEN
278                to_be_resorted => v
279             ELSE
280                to_be_resorted => v_av
281             ENDIF
282
283          CASE ( 'vpt' )
284             IF ( av == 0 )  THEN
285                to_be_resorted => vpt
286             ELSE
287                to_be_resorted => vpt_av
288             ENDIF
289
290          CASE ( 'w' )
291             IF ( av == 0 )  THEN
292                to_be_resorted => w
293             ELSE
294                to_be_resorted => w_av
295             ENDIF
296
297          CASE DEFAULT
298!
299!--          User defined quantity
300             CALL user_data_output_mask(av, domask(mid,av,if), found, local_pf )
301             resorted = .TRUE.
302
303             IF ( .NOT. found )  THEN
304                WRITE ( message_string, * ) 'no output available for: ', &
305                                            TRIM( domask(mid,av,if) )
[564]306                CALL message( 'data_output_mask', 'PA0327', 0, 0, 0, 6, 0 )
[298]307             ENDIF
308
309       END SELECT
310
311!
312!--    Resort the array to be output, if not done above
313       IF ( .NOT. resorted )  THEN
314          DO  i = 1, mask_size_l(mid,1)
315             DO  j = 1, mask_size_l(mid,2)
316                DO  k = 1, mask_size_l(mid,3)
317                   local_pf(i,j,k) =  to_be_resorted(mask_k(mid,k), &
318                                      mask_j(mid,j),mask_i(mid,i))
319                ENDDO
320             ENDDO
321          ENDDO
322       ENDIF
323
324!
325!--    I/O block. I/O methods are implemented
326!--    (1) for parallel execution
327!--     a. with NetCDF 4 parallel I/O-enabled library
328!--     b. with NetCDF 3 library
329!--    (2) for serial execution.
330!--    The choice of method depends on the correct setting of preprocessor
331!--    directives __parallel and __netcdf4 as well as on the parameter
[493]332!--    netcdf_data_format.
[298]333#if defined( __parallel )
334#if defined( __netcdf4 )
[493]335       IF ( netcdf_data_format > 2 )  THEN
[298]336!
337!--       (1) a. Parallel I/O using NetCDF 4 (not yet tested)
338          nc_stat = NF90_PUT_VAR( id_set_mask(mid,av),  &
[340]339               id_var_domask(mid,av,if),  &
[409]340               local_pf,  &
[340]341               start = (/ mask_start_l(mid,1), mask_start_l(mid,2),  &
[409]342                          mask_start_l(mid,3), domask_time_count(mid,av) /),  &
343               count = (/ mask_size_l(mid,1), mask_size_l(mid,2),  &
344                          mask_size_l(mid,3), 1 /) )
[564]345          CALL handle_netcdf_error( 'data_output_mask', 461 )
[298]346       ELSE
347#endif
348!
349!--       (1) b. Conventional I/O only through PE0
350!--       PE0 receives partial arrays from all processors of the respective mask
351!--       and outputs them. Here a barrier has to be set, because otherwise
352!--       "-MPI- FATAL: Remote protocol queue full" may occur.
353          CALL MPI_BARRIER( comm2d, ierr )
354
355          ngp = mask_size_l(mid,1) * mask_size_l(mid,2) * mask_size_l(mid,3)
356          IF ( myid == 0 )  THEN
357!
358!--          Local array can be relocated directly.
359             total_pf( &
360               mask_start_l(mid,1):mask_start_l(mid,1)+mask_size_l(mid,1)-1, &
361               mask_start_l(mid,2):mask_start_l(mid,2)+mask_size_l(mid,2)-1, &
362               mask_start_l(mid,3):mask_start_l(mid,3)+mask_size_l(mid,3)-1 ) &
363               = local_pf
364!
365!--          Receive data from all other PEs.
366             DO  n = 1, numprocs-1
367!
368!--             Receive index limits first, then array.
369!--             Index limits are received in arbitrary order from the PEs.
370                CALL MPI_RECV( ind(1), 6, MPI_INTEGER, MPI_ANY_SOURCE, 0,  &
371                     comm2d, status, ierr )
372!
373!--             Not all PEs have data for the mask
374                IF ( ind(1) /= -9999 )  THEN
375                   ngp = ( ind(2)-ind(1)+1 ) * (ind(4)-ind(3)+1 ) *  &
376                         ( ind(6)-ind(5)+1 )
377                   sender = status(MPI_SOURCE)
378                   DEALLOCATE( local_pf )
379                   ALLOCATE(local_pf(ind(1):ind(2),ind(3):ind(4),ind(5):ind(6)))
380                   CALL MPI_RECV( local_pf(ind(1),ind(3),ind(5)), ngp,  &
381                        MPI_REAL, sender, 1, comm2d, status, ierr )
382                   total_pf(ind(1):ind(2),ind(3):ind(4),ind(5):ind(6)) &
383                        = local_pf
384                ENDIF
385             ENDDO
386
387             nc_stat = NF90_PUT_VAR( id_set_mask(mid,av),  &
388                  id_var_domask(mid,av,if), total_pf, &
389                  start = (/ 1, 1, 1, domask_time_count(mid,av) /), &
390                  count = (/ mask_size(mid,1), mask_size(mid,2), &
391                             mask_size(mid,3), 1 /) )
[564]392             CALL handle_netcdf_error( 'data_output_mask', 462 )
[298]393
394          ELSE
395!
396!--          If at least part of the mask resides on the PE, send the index
397!--          limits for the target array, otherwise send -9999 to PE0.
398             IF ( mask_size_l(mid,1) > 0 .AND.  mask_size_l(mid,2) > 0 .AND. &
399                  mask_size_l(mid,3) > 0  ) &
400                  THEN
401                ind(1) = mask_start_l(mid,1)
402                ind(2) = mask_start_l(mid,1) + mask_size_l(mid,1) - 1
403                ind(3) = mask_start_l(mid,2)
404                ind(4) = mask_start_l(mid,2) + mask_size_l(mid,2) - 1
405                ind(5) = mask_start_l(mid,3)
406                ind(6) = mask_start_l(mid,3) + mask_size_l(mid,3) - 1
407             ELSE
408                ind(1) = -9999; ind(2) = -9999
409                ind(3) = -9999; ind(4) = -9999
410                ind(5) = -9999; ind(6) = -9999
411             ENDIF
412             CALL MPI_SEND( ind(1), 6, MPI_INTEGER, 0, 0, comm2d, ierr )
413!
414!--          If applicable, send data to PE0.
415             IF ( ind(1) /= -9999 )  THEN
416                CALL MPI_SEND( local_pf(1,1,1), ngp, MPI_REAL, 0, 1, comm2d, &
417                     ierr )
418             ENDIF
419          ENDIF
420!
421!--       A barrier has to be set, because otherwise some PEs may proceed too
422!--       fast so that PE0 may receive wrong data on tag 0.
423          CALL MPI_BARRIER( comm2d, ierr )
424#if defined( __netcdf4 )
425       ENDIF
426#endif
427#else
428!
429!--    (2) For serial execution of PALM, the single processor (PE0) holds all
430!--    data and writes them directly to file.
431       nc_stat = NF90_PUT_VAR( id_set_mask(mid,av),  &
432            id_var_domask(mid,av,if),       &
[475]433            local_pf, &
[298]434            start = (/ 1, 1, 1, domask_time_count(mid,av) /), &
435            count = (/ mask_size_l(mid,1), mask_size_l(mid,2), &
436                       mask_size_l(mid,3), 1 /) )
[564]437       CALL handle_netcdf_error( 'data_output_mask', 463 )
[298]438#endif
439
440       if = if + 1
441    ENDDO
442
443!
444!-- Deallocate temporary arrays.
445    DEALLOCATE( local_pf )
446#if defined( __parallel )
447    IF ( myid == 0 )  THEN
448       DEALLOCATE( total_pf )
449    ENDIF
450#endif
451
452
453    CALL cpu_log (log_point(49),'data_output_mask','stop','nobarrier')
454#endif
455
456 END SUBROUTINE data_output_mask
Note: See TracBrowser for help on using the repository browser.