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

Last change on this file since 1692 was 1692, checked in by maronga, 8 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 23.8 KB
Line 
1!> @file data_output_mask.f90
2!--------------------------------------------------------------------------------!
3! This file is part of PALM.
4!
5! PALM is free software: you can redistribute it and/or modify it under the terms
6! of the GNU General Public License as published by the Free Software Foundation,
7! either version 3 of the License, or (at your option) any later version.
8!
9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
10! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
11! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
12!
13! You should have received a copy of the GNU General Public License along with
14! PALM. If not, see <http://www.gnu.org/licenses/>.
15!
16! Copyright 1997-2015 Leibniz Universitaet Hannover
17!--------------------------------------------------------------------------------!
18!
19! Current revisions:
20! -----------------
21!
22!
23! Former revisions:
24! -----------------
25! $Id: data_output_mask.f90 1692 2015-10-26 16:29:17Z maronga $
26!
27! 1691 2015-10-26 16:17:44Z maronga
28! Added output of radiative heating rates for RRTMG
29!
30! 1682 2015-10-07 23:56:08Z knoop
31! Code annotations made doxygen readable
32!
33! 1585 2015-04-30 07:05:52Z maronga
34! Added support for RRTMG
35!
36! 1438 2014-07-22 14:14:06Z heinze
37! +nr, qc, qr
38!
39! 1359 2014-04-11 17:15:14Z hoffmann
40! New particle structure integrated.
41!
42! 1353 2014-04-08 15:21:23Z heinze
43! REAL constants provided with KIND-attribute
44!
45! 1327 2014-03-21 11:00:16Z raasch
46!
47!
48! 1320 2014-03-20 08:40:49Z raasch
49! ONLY-attribute added to USE-statements,
50! kind-parameters added to all INTEGER and REAL declaration statements,
51! kinds are defined in new module kinds,
52! revision history before 2012 removed,
53! comment fields (!:) to be used for variable explanations added to
54! all variable declaration statements
55!
56! 1318 2014-03-17 13:35:16Z raasch
57! barrier argument removed from cpu_log,
58! module interfaces removed
59!
60! 1092 2013-02-02 11:24:22Z raasch
61! unused variables removed
62!
63! 1036 2012-10-22 13:43:42Z raasch
64! code put under GPL (PALM 3.9)
65!
66! 1031 2012-10-19 14:35:30Z raasch
67! netCDF4 without parallel file support implemented
68!
69! 1007 2012-09-19 14:30:36Z franke
70! Bugfix: calculation of pr must depend on the particle weighting factor,
71! missing calculation of ql_vp added
72!
73! 410 2009-12-04 17:05:40Z letzel
74! Initial version
75!
76! Description:
77! ------------
78!> Masked data output in netCDF format for current mask (current value of mid).
79!------------------------------------------------------------------------------!
80 SUBROUTINE data_output_mask( av )
81
82 
83
84#if defined( __netcdf )
85    USE arrays_3d,                                                             &
86        ONLY:  e, nr, p, pt, q, qc, ql, ql_c, ql_v, qr, rho, sa, tend, u,      &
87               v, vpt, w
88   
89    USE averaging,                                                             &
90        ONLY:  e_av, lpt_av, nr_av, p_av, pc_av, pr_av, pt_av, q_av, qc_av,    &
91               ql_av, ql_c_av, ql_v_av, ql_vp_av, qv_av, qr_av, rho_av, s_av,  &
92               sa_av, u_av, v_av, vpt_av, w_av 
93   
94    USE cloud_parameters,                                                      &
95        ONLY:  l_d_cp, pt_d_t
96   
97    USE control_parameters,                                                    &
98        ONLY:  cloud_physics, domask, domask_no, domask_time_count,            &
99               mask_i, mask_j, mask_k, mask_size, mask_size_l,                 &
100               mask_start_l, max_masks, message_string, mid,                   &
101               netcdf_data_format, nz_do3d, simulated_time
102    USE cpulog,                                                                &
103        ONLY:  cpu_log, log_point
104
105
106   
107    USE indices,                                                               &
108        ONLY:  nbgp, nxl, nxr, nyn, nys, nzb, nzt
109       
110    USE kinds
111   
112    USE netcdf
113   
114    USE netcdf_control
115   
116    USE particle_attributes,                                                   &
117        ONLY:  grid_particles, number_of_particles, particles,                 &
118               particle_advection_start, prt_count
119   
120    USE pegrid
121
122    USE radiation_model_mod,                                                   &
123        ONLY:  rad_lw_in, rad_lw_in_av, rad_lw_out, rad_lw_out_av,             &
124               rad_lw_cs_hr, rad_lw_cs_hr_av, rad_lw_hr, rad_lw_hr_av,         &
125               rad_sw_in, rad_sw_in_av, rad_sw_out, rad_sw_out_av,             &
126               rad_sw_cs_hr, rad_sw_cs_hr_av, rad_sw_hr, rad_sw_hr_av
127
128
129    IMPLICIT NONE
130
131    INTEGER(iwp) ::  av       !<
132    INTEGER(iwp) ::  ngp      !<
133    INTEGER(iwp) ::  i        !<
134    INTEGER(iwp) ::  if       !<
135    INTEGER(iwp) ::  j        !<
136    INTEGER(iwp) ::  k        !<
137    INTEGER(iwp) ::  n        !<
138    INTEGER(iwp) ::  psi      !<
139    INTEGER(iwp) ::  sender   !<
140    INTEGER(iwp) ::  ind(6)   !<
141   
142    LOGICAL ::  found         !<
143    LOGICAL ::  resorted      !<
144   
145    REAL(wp) ::  mean_r       !<
146    REAL(wp) ::  s_r2         !<
147    REAL(wp) ::  s_r3         !<
148   
149    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  local_pf    !<
150#if defined( __parallel )
151    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  total_pf    !<
152#endif
153    REAL(wp), DIMENSION(:,:,:), POINTER ::  to_be_resorted  !<
154
155!
156!-- Return, if nothing to output
157    IF ( domask_no(mid,av) == 0 )  RETURN
158
159    CALL cpu_log (log_point(49),'data_output_mask','start')
160
161!
162!-- Open output file.
163    IF ( myid == 0  .OR.  netcdf_data_format > 4 )  THEN
164       CALL check_open( 200+mid+av*max_masks )
165    ENDIF 
166
167!
168!-- Allocate total and local output arrays.
169#if defined( __parallel )
170    IF ( myid == 0 )  THEN
171       ALLOCATE( total_pf(mask_size(mid,1),mask_size(mid,2),mask_size(mid,3)) )
172    ENDIF
173#endif
174    ALLOCATE( local_pf(mask_size_l(mid,1),mask_size_l(mid,2), &
175                       mask_size_l(mid,3)) )
176
177!
178!-- Update the netCDF time axis.
179    domask_time_count(mid,av) = domask_time_count(mid,av) + 1
180    IF ( myid == 0  .OR.  netcdf_data_format > 4 )  THEN
181       nc_stat = NF90_PUT_VAR( id_set_mask(mid,av), id_var_time_mask(mid,av), &
182                               (/ simulated_time /),                          &
183                               start = (/ domask_time_count(mid,av) /),       &
184                               count = (/ 1 /) )
185       CALL handle_netcdf_error( 'data_output_mask', 460 )
186    ENDIF
187
188!
189!-- Loop over all variables to be written.
190    if = 1
191
192    DO  WHILE ( domask(mid,av,if)(1:1) /= ' ' )
193!
194!--    Reallocate local_pf on PE 0 since its shape changes during MPI exchange
195       IF ( netcdf_data_format < 5   .AND.  myid == 0  .AND.  if > 1 )  THEN
196          DEALLOCATE( local_pf )
197          ALLOCATE( local_pf(mask_size_l(mid,1),mask_size_l(mid,2), &
198                             mask_size_l(mid,3)) )
199       ENDIF
200!
201!--    Store the variable chosen.
202       resorted = .FALSE.
203       SELECT CASE ( TRIM( domask(mid,av,if) ) )
204
205          CASE ( 'e' )
206             IF ( av == 0 )  THEN
207                to_be_resorted => e
208             ELSE
209                to_be_resorted => e_av
210             ENDIF
211
212          CASE ( 'lpt' )
213             IF ( av == 0 )  THEN
214                to_be_resorted => pt
215             ELSE
216                to_be_resorted => lpt_av
217             ENDIF
218
219          CASE ( 'nr' )
220             IF ( av == 0 )  THEN
221                to_be_resorted => nr
222             ELSE
223                to_be_resorted => nr_av
224             ENDIF
225
226          CASE ( 'p' )
227             IF ( av == 0 )  THEN
228                to_be_resorted => p
229             ELSE
230                to_be_resorted => p_av
231             ENDIF
232
233          CASE ( 'pc' )  ! particle concentration (requires ghostpoint exchange)
234             IF ( av == 0 )  THEN
235                tend = prt_count
236                CALL exchange_horiz( tend, nbgp )
237                DO  i = 1, mask_size_l(mid,1)
238                   DO  j = 1, mask_size_l(mid,2)
239                      DO  k = 1, mask_size_l(mid,3)
240                         local_pf(i,j,k) =  tend(mask_k(mid,k), &
241                                   mask_j(mid,j),mask_i(mid,i))
242                      ENDDO
243                   ENDDO
244                ENDDO
245                resorted = .TRUE.
246             ELSE
247                CALL exchange_horiz( pc_av, nbgp )
248                to_be_resorted => pc_av
249             ENDIF
250
251          CASE ( 'pr' )  ! mean particle radius (effective radius)
252             IF ( av == 0 )  THEN
253                IF ( simulated_time >= particle_advection_start )  THEN
254                   DO  i = nxl, nxr
255                      DO  j = nys, nyn
256                         DO  k = nzb, nz_do3d
257                            number_of_particles = prt_count(k,j,i)
258                            IF (number_of_particles <= 0)  CYCLE
259                            particles => grid_particles(k,j,i)%particles(1:number_of_particles)
260                            s_r2 = 0.0_wp
261                            s_r3 = 0.0_wp
262                            DO  n = 1, number_of_particles
263                               IF ( particles(n)%particle_mask )  THEN
264                                  s_r2 = s_r2 + grid_particles(k,j,i)%particles(n)%radius**2 * &
265                                         grid_particles(k,j,i)%particles(n)%weight_factor
266                                  s_r3 = s_r3 + grid_particles(k,j,i)%particles(n)%radius**3 * &
267                                         grid_particles(k,j,i)%particles(n)%weight_factor
268                               ENDIF
269                            ENDDO
270                            IF ( s_r2 > 0.0_wp )  THEN
271                               mean_r = s_r3 / s_r2
272                            ELSE
273                               mean_r = 0.0_wp
274                            ENDIF
275                            tend(k,j,i) = mean_r
276                         ENDDO
277                      ENDDO
278                   ENDDO
279                   CALL exchange_horiz( tend, nbgp )
280                ELSE
281                   tend = 0.0_wp
282                ENDIF
283                DO  i = 1, mask_size_l(mid,1)
284                   DO  j = 1, mask_size_l(mid,2)
285                      DO  k = 1, mask_size_l(mid,3)
286                         local_pf(i,j,k) =  tend(mask_k(mid,k), &
287                                   mask_j(mid,j),mask_i(mid,i))
288                      ENDDO
289                   ENDDO
290                ENDDO
291                resorted = .TRUE.
292             ELSE
293                CALL exchange_horiz( pr_av, nbgp )
294                to_be_resorted => pr_av
295             ENDIF
296
297          CASE ( 'pt' )
298             IF ( av == 0 )  THEN
299                IF ( .NOT. cloud_physics ) THEN
300                   to_be_resorted => pt
301                ELSE
302                   DO  i = 1, mask_size_l(mid,1)
303                      DO  j = 1, mask_size_l(mid,2)
304                         DO  k = 1, mask_size_l(mid,3)
305                            local_pf(i,j,k) =  &
306                                 pt(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i)) &
307                                 + l_d_cp * pt_d_t(mask_k(mid,k)) * &
308                                   ql(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i))
309                         ENDDO
310                      ENDDO
311                   ENDDO
312                   resorted = .TRUE.
313                ENDIF
314             ELSE
315                to_be_resorted => pt_av
316             ENDIF
317
318          CASE ( 'q' )
319             IF ( av == 0 )  THEN
320                to_be_resorted => q
321             ELSE
322                to_be_resorted => q_av
323             ENDIF
324
325          CASE ( 'qc' )
326             IF ( av == 0 )  THEN
327                to_be_resorted => qc
328             ELSE
329                to_be_resorted => qc_av
330             ENDIF
331
332          CASE ( 'ql' )
333             IF ( av == 0 )  THEN
334                to_be_resorted => ql
335             ELSE
336                to_be_resorted => ql_av
337             ENDIF
338
339          CASE ( 'ql_c' )
340             IF ( av == 0 )  THEN
341                to_be_resorted => ql_c
342             ELSE
343                to_be_resorted => ql_c_av
344             ENDIF
345
346          CASE ( 'ql_v' )
347             IF ( av == 0 )  THEN
348                to_be_resorted => ql_v
349             ELSE
350                to_be_resorted => ql_v_av
351             ENDIF
352
353          CASE ( 'ql_vp' )
354             IF ( av == 0 )  THEN
355                IF ( simulated_time >= particle_advection_start )  THEN
356                   DO  i = nxl, nxr
357                      DO  j = nys, nyn
358                         DO  k = nzb, nz_do3d
359                            number_of_particles = prt_count(k,j,i)
360                            IF (number_of_particles <= 0)  CYCLE
361                            particles => grid_particles(k,j,i)%particles(1:number_of_particles)
362                            DO  n = 1, number_of_particles
363                               IF ( particles(n)%particle_mask )  THEN
364                                  tend(k,j,i) = tend(k,j,i) + &
365                                          particles(n)%weight_factor / &
366                                          prt_count(k,j,i)
367                               ENDIF
368                            ENDDO
369                         ENDDO
370                      ENDDO
371                   ENDDO
372                   CALL exchange_horiz( tend, nbgp )
373                ELSE
374                   tend = 0.0_wp
375                ENDIF
376                DO  i = 1, mask_size_l(mid,1)
377                   DO  j = 1, mask_size_l(mid,2)
378                      DO  k = 1, mask_size_l(mid,3)
379                         local_pf(i,j,k) =  tend(mask_k(mid,k), &
380                                   mask_j(mid,j),mask_i(mid,i))
381                      ENDDO
382                   ENDDO
383                ENDDO
384                resorted = .TRUE.
385             ELSE
386                CALL exchange_horiz( ql_vp_av, nbgp )
387                to_be_resorted => ql_vp_av
388             ENDIF
389
390          CASE ( 'qv' )
391             IF ( av == 0 )  THEN
392                DO  i = 1, mask_size_l(mid,1)
393                   DO  j = 1, mask_size_l(mid,2)
394                      DO  k = 1, mask_size_l(mid,3)
395                         local_pf(i,j,k) =  &
396                              q(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i)) -  &
397                              ql(mask_k(mid,k),mask_j(mid,j),mask_i(mid,i))
398                      ENDDO
399                   ENDDO
400                ENDDO
401                resorted = .TRUE.
402             ELSE
403                to_be_resorted => qv_av
404             ENDIF
405
406          CASE ( 'qr' )
407             IF ( av == 0 )  THEN
408                to_be_resorted => qr
409             ELSE
410                to_be_resorted => qr_av
411             ENDIF
412
413          CASE ( 'rad_lw_in' )
414             IF ( av == 0 )  THEN
415                to_be_resorted => rad_lw_in
416             ELSE
417                to_be_resorted => rad_lw_in_av
418             ENDIF
419
420          CASE ( 'rad_lw_out' )
421             IF ( av == 0 )  THEN
422                to_be_resorted => rad_lw_out
423             ELSE
424                to_be_resorted => rad_lw_out_av
425             ENDIF
426
427          CASE ( 'rad_lw_cs_hr' )
428             IF ( av == 0 )  THEN
429                to_be_resorted => rad_lw_cs_hr
430             ELSE
431                to_be_resorted => rad_lw_cs_hr_av
432             ENDIF
433
434          CASE ( 'rad_lw_hr' )
435             IF ( av == 0 )  THEN
436                to_be_resorted => rad_lw_hr
437             ELSE
438                to_be_resorted => rad_lw_hr_av
439             ENDIF
440
441          CASE ( 'rad_sw_in' )
442             IF ( av == 0 )  THEN
443                to_be_resorted => rad_sw_in
444             ELSE
445                to_be_resorted => rad_sw_in_av
446             ENDIF
447
448          CASE ( 'rad_sw_out' )
449             IF ( av == 0 )  THEN
450                to_be_resorted => rad_sw_out
451             ELSE
452                to_be_resorted => rad_sw_out_av
453             ENDIF
454
455          CASE ( 'rad_sw_cs_hr' )
456             IF ( av == 0 )  THEN
457                to_be_resorted => rad_sw_cs_hr
458             ELSE
459                to_be_resorted => rad_sw_cs_hr_av
460             ENDIF
461
462          CASE ( 'rad_sw_hr' )
463             IF ( av == 0 )  THEN
464                to_be_resorted => rad_sw_hr
465             ELSE
466                to_be_resorted => rad_sw_hr_av
467             ENDIF
468
469          CASE ( 'rho' )
470             IF ( av == 0 )  THEN
471                to_be_resorted => rho
472             ELSE
473                to_be_resorted => rho_av
474             ENDIF
475
476          CASE ( 's' )
477             IF ( av == 0 )  THEN
478                to_be_resorted => q
479             ELSE
480                to_be_resorted => s_av
481             ENDIF
482
483          CASE ( 'sa' )
484             IF ( av == 0 )  THEN
485                to_be_resorted => sa
486             ELSE
487                to_be_resorted => sa_av
488             ENDIF
489
490          CASE ( 'u' )
491             IF ( av == 0 )  THEN
492                to_be_resorted => u
493             ELSE
494                to_be_resorted => u_av
495             ENDIF
496
497          CASE ( 'v' )
498             IF ( av == 0 )  THEN
499                to_be_resorted => v
500             ELSE
501                to_be_resorted => v_av
502             ENDIF
503
504          CASE ( 'vpt' )
505             IF ( av == 0 )  THEN
506                to_be_resorted => vpt
507             ELSE
508                to_be_resorted => vpt_av
509             ENDIF
510
511          CASE ( 'w' )
512             IF ( av == 0 )  THEN
513                to_be_resorted => w
514             ELSE
515                to_be_resorted => w_av
516             ENDIF
517
518          CASE DEFAULT
519!
520!--          User defined quantity
521             CALL user_data_output_mask(av, domask(mid,av,if), found, local_pf )
522             resorted = .TRUE.
523
524             IF ( .NOT. found )  THEN
525                WRITE ( message_string, * ) 'no output available for: ', &
526                                            TRIM( domask(mid,av,if) )
527                CALL message( 'data_output_mask', 'PA0327', 0, 0, 0, 6, 0 )
528             ENDIF
529
530       END SELECT
531
532!
533!--    Resort the array to be output, if not done above
534       IF ( .NOT. resorted )  THEN
535          DO  i = 1, mask_size_l(mid,1)
536             DO  j = 1, mask_size_l(mid,2)
537                DO  k = 1, mask_size_l(mid,3)
538                   local_pf(i,j,k) =  to_be_resorted(mask_k(mid,k), &
539                                      mask_j(mid,j),mask_i(mid,i))
540                ENDDO
541             ENDDO
542          ENDDO
543       ENDIF
544
545!
546!--    I/O block. I/O methods are implemented
547!--    (1) for parallel execution
548!--     a. with netCDF 4 parallel I/O-enabled library
549!--     b. with netCDF 3 library
550!--    (2) for serial execution.
551!--    The choice of method depends on the correct setting of preprocessor
552!--    directives __parallel and __netcdf4_parallel as well as on the parameter
553!--    netcdf_data_format.
554#if defined( __parallel )
555#if defined( __netcdf4_parallel )
556       IF ( netcdf_data_format > 4 )  THEN
557!
558!--       (1) a. Parallel I/O using netCDF 4 (not yet tested)
559          nc_stat = NF90_PUT_VAR( id_set_mask(mid,av),  &
560               id_var_domask(mid,av,if),  &
561               local_pf,  &
562               start = (/ mask_start_l(mid,1), mask_start_l(mid,2),  &
563                          mask_start_l(mid,3), domask_time_count(mid,av) /),  &
564               count = (/ mask_size_l(mid,1), mask_size_l(mid,2),  &
565                          mask_size_l(mid,3), 1 /) )
566          CALL handle_netcdf_error( 'data_output_mask', 461 )
567       ELSE
568#endif
569!
570!--       (1) b. Conventional I/O only through PE0
571!--       PE0 receives partial arrays from all processors of the respective mask
572!--       and outputs them. Here a barrier has to be set, because otherwise
573!--       "-MPI- FATAL: Remote protocol queue full" may occur.
574          CALL MPI_BARRIER( comm2d, ierr )
575
576          ngp = mask_size_l(mid,1) * mask_size_l(mid,2) * mask_size_l(mid,3)
577          IF ( myid == 0 )  THEN
578!
579!--          Local array can be relocated directly.
580             total_pf( &
581               mask_start_l(mid,1):mask_start_l(mid,1)+mask_size_l(mid,1)-1, &
582               mask_start_l(mid,2):mask_start_l(mid,2)+mask_size_l(mid,2)-1, &
583               mask_start_l(mid,3):mask_start_l(mid,3)+mask_size_l(mid,3)-1 ) &
584               = local_pf
585!
586!--          Receive data from all other PEs.
587             DO  n = 1, numprocs-1
588!
589!--             Receive index limits first, then array.
590!--             Index limits are received in arbitrary order from the PEs.
591                CALL MPI_RECV( ind(1), 6, MPI_INTEGER, MPI_ANY_SOURCE, 0,  &
592                     comm2d, status, ierr )
593!
594!--             Not all PEs have data for the mask
595                IF ( ind(1) /= -9999 )  THEN
596                   ngp = ( ind(2)-ind(1)+1 ) * (ind(4)-ind(3)+1 ) *  &
597                         ( ind(6)-ind(5)+1 )
598                   sender = status(MPI_SOURCE)
599                   DEALLOCATE( local_pf )
600                   ALLOCATE(local_pf(ind(1):ind(2),ind(3):ind(4),ind(5):ind(6)))
601                   CALL MPI_RECV( local_pf(ind(1),ind(3),ind(5)), ngp,  &
602                        MPI_REAL, sender, 1, comm2d, status, ierr )
603                   total_pf(ind(1):ind(2),ind(3):ind(4),ind(5):ind(6)) &
604                        = local_pf
605                ENDIF
606             ENDDO
607
608             nc_stat = NF90_PUT_VAR( id_set_mask(mid,av),  &
609                  id_var_domask(mid,av,if), total_pf, &
610                  start = (/ 1, 1, 1, domask_time_count(mid,av) /), &
611                  count = (/ mask_size(mid,1), mask_size(mid,2), &
612                             mask_size(mid,3), 1 /) )
613             CALL handle_netcdf_error( 'data_output_mask', 462 )
614
615          ELSE
616!
617!--          If at least part of the mask resides on the PE, send the index
618!--          limits for the target array, otherwise send -9999 to PE0.
619             IF ( mask_size_l(mid,1) > 0 .AND.  mask_size_l(mid,2) > 0 .AND. &
620                  mask_size_l(mid,3) > 0  ) &
621                  THEN
622                ind(1) = mask_start_l(mid,1)
623                ind(2) = mask_start_l(mid,1) + mask_size_l(mid,1) - 1
624                ind(3) = mask_start_l(mid,2)
625                ind(4) = mask_start_l(mid,2) + mask_size_l(mid,2) - 1
626                ind(5) = mask_start_l(mid,3)
627                ind(6) = mask_start_l(mid,3) + mask_size_l(mid,3) - 1
628             ELSE
629                ind(1) = -9999; ind(2) = -9999
630                ind(3) = -9999; ind(4) = -9999
631                ind(5) = -9999; ind(6) = -9999
632             ENDIF
633             CALL MPI_SEND( ind(1), 6, MPI_INTEGER, 0, 0, comm2d, ierr )
634!
635!--          If applicable, send data to PE0.
636             IF ( ind(1) /= -9999 )  THEN
637                CALL MPI_SEND( local_pf(1,1,1), ngp, MPI_REAL, 0, 1, comm2d, &
638                     ierr )
639             ENDIF
640          ENDIF
641!
642!--       A barrier has to be set, because otherwise some PEs may proceed too
643!--       fast so that PE0 may receive wrong data on tag 0.
644          CALL MPI_BARRIER( comm2d, ierr )
645#if defined( __netcdf4_parallel )
646       ENDIF
647#endif
648#else
649!
650!--    (2) For serial execution of PALM, the single processor (PE0) holds all
651!--    data and writes them directly to file.
652       nc_stat = NF90_PUT_VAR( id_set_mask(mid,av),  &
653            id_var_domask(mid,av,if),       &
654            local_pf, &
655            start = (/ 1, 1, 1, domask_time_count(mid,av) /), &
656            count = (/ mask_size_l(mid,1), mask_size_l(mid,2), &
657                       mask_size_l(mid,3), 1 /) )
658       CALL handle_netcdf_error( 'data_output_mask', 463 )
659#endif
660
661       if = if + 1
662
663    ENDDO
664
665!
666!-- Deallocate temporary arrays.
667    DEALLOCATE( local_pf )
668#if defined( __parallel )
669    IF ( myid == 0 )  THEN
670       DEALLOCATE( total_pf )
671    ENDIF
672#endif
673
674
675    CALL cpu_log( log_point(49), 'data_output_mask', 'stop' )
676#endif
677
678 END SUBROUTINE data_output_mask
Note: See TracBrowser for help on using the repository browser.