source: palm/trunk/SOURCE/data_output_2d.f90 @ 622

Last change on this file since 622 was 622, checked in by raasch, 11 years ago

New:
---

Optional barriers included in order to speed up collective operations
MPI_ALLTOALL and MPI_ALLREDUCE. This feature is controlled with new initial
parameter collective_wait. Default is .FALSE, but .TRUE. on SGI-type
systems. (advec_particles, advec_s_bc, buoyancy, check_for_restart,
cpu_statistics, data_output_2d, data_output_ptseries, flow_statistics,
global_min_max, inflow_turbulence, init_3d_model, init_particles, init_pegrid,
init_slope, parin, pres, poismg, set_particle_attributes, timestep,
read_var_list, user_statistics, write_compressed, write_var_list)

Adjustments for Kyushu Univ. (lcrte, ibmku). Concerning hybrid
(MPI/openMP) runs, the number of openMP threads per MPI tasks can now
be given as an argument to mrun-option -O. (mbuild, mrun, subjob)

Changed:


Initialization of the module command changed for SGI-ICE/lcsgi (mbuild, subjob)

Errors:


  • Property svn:keywords set to Id
File size: 63.3 KB
Line 
1 SUBROUTINE data_output_2d( mode, av )
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6! optional barriers included in order to speed up collective operations
7!
8! Former revisions:
9! -----------------
10! $Id: data_output_2d.f90 622 2010-12-10 08:08:13Z raasch $
11!
12! 493 2010-03-01 08:30:24Z raasch
13! NetCDF4 support (parallel output)
14!
15! 367 2009-08-25 08:35:52Z maronga
16! simulated_time in NetCDF output replaced by time_since_reference_point.
17! Output of NetCDF messages with aid of message handling routine.
18! Bugfix: averaging along z is not allowed for 2d quantities (e.g. u* and z0)
19! Output of messages replaced by message handling routine.
20! Output of user defined 2D (XY) arrays at z=nzb+1 is now possible
21! Bugfix: to_be_resorted => s_av for time-averaged scalars
22! Calculation of shf* and qsws* added.
23!
24! 215 2008-11-18 09:54:31Z raasch
25! Bugfix: no output of particle concentration and radius unless particles
26! have been started
27!
28! 96 2007-06-04 08:07:41Z raasch
29! Output of density and salinity
30!
31! 75 2007-03-22 09:54:05Z raasch
32! Output of precipitation amount/rate and roughness length,
33! 2nd+3rd argument removed from exchange horiz
34!
35! RCS Log replace by Id keyword, revision history cleaned up
36!
37! Revision 1.5  2006/08/22 13:50:29  raasch
38! xz and yz cross sections now up to nzt+1
39!
40! Revision 1.2  2006/02/23 10:19:22  raasch
41! Output of time-averaged data, output of averages along x, y, or z,
42! output of user-defined quantities,
43! section data are copied from local_pf to local_2d before they are output,
44! output of particle concentration and mean radius,
45! Former subroutine plot_2d renamed data_output_2d, pl2d.. renamed do2d..,
46! anz renamed ngp, ebene renamed section, pl2d_.._anz renamed do2d_.._n
47!
48! Revision 1.1  1997/08/11 06:24:09  raasch
49! Initial revision
50!
51!
52! Description:
53! ------------
54! Data output of horizontal cross-sections in NetCDF format or binary format
55! compatible to old graphic software iso2d.
56! Attention: The position of the sectional planes is still not always computed
57! ---------  correctly. (zu is used always)!
58!------------------------------------------------------------------------------!
59
60    USE arrays_3d
61    USE averaging
62    USE cloud_parameters
63    USE control_parameters
64    USE cpulog
65    USE grid_variables
66    USE indices
67    USE interfaces
68    USE netcdf_control
69    USE particle_attributes
70    USE pegrid
71
72    IMPLICIT NONE
73
74    CHARACTER (LEN=2)  ::  do2d_mode, mode
75    CHARACTER (LEN=4)  ::  grid
76    CHARACTER (LEN=25) ::  section_chr
77    CHARACTER (LEN=50) ::  rtext
78    INTEGER ::  av, ngp, file_id, i, if, is, iis, j, k, l, layer_xy, n, psi, &
79                s, sender, &
80                ind(4)
81    LOGICAL ::  found, resorted, two_d
82    REAL    ::  mean_r, s_r3, s_r4
83    REAL, DIMENSION(:), ALLOCATABLE ::      level_z
84    REAL, DIMENSION(:,:), ALLOCATABLE ::    local_2d, local_2d_l
85    REAL, DIMENSION(:,:,:), ALLOCATABLE ::  local_pf
86#if defined( __parallel )
87    REAL, DIMENSION(:,:),   ALLOCATABLE ::  total_2d
88#endif
89    REAL, DIMENSION(:,:,:), POINTER ::  to_be_resorted
90
91    NAMELIST /LOCAL/  rtext
92
93    CALL cpu_log (log_point(3),'data_output_2d','start')
94
95!
96!-- Immediate return, if no output is requested (no respective sections
97!-- found in parameter data_output)
98    IF ( mode == 'xy'  .AND.  .NOT. data_output_xy(av) )  RETURN
99    IF ( mode == 'xz'  .AND.  .NOT. data_output_xz(av) )  RETURN
100    IF ( mode == 'yz'  .AND.  .NOT. data_output_yz(av) )  RETURN
101
102    two_d = .FALSE.    ! local variable to distinguish between output of pure 2D
103                       ! arrays and cross-sections of 3D arrays.
104
105!
106!-- Depending on the orientation of the cross-section, the respective output
107!-- files have to be opened.
108    SELECT CASE ( mode )
109
110       CASE ( 'xy' )
111
112          s = 1
113          ALLOCATE( level_z(0:nzt+1), local_2d(nxl-1:nxr+1,nys-1:nyn+1) )
114
115!
116!--       Classic and 64bit offset NetCDF output is done only on PE0.
117!--       netCDF4/HDF5 output is done in parallel on all PEs.
118          IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
119          THEN
120             CALL check_open( 101+av*10 )
121          ENDIF
122
123          IF ( data_output_2d_on_each_pe )  THEN
124             CALL check_open( 21 )
125          ELSE
126             IF ( myid == 0 )  THEN
127                IF ( iso2d_output )  CALL check_open( 21 )
128#if defined( __parallel )
129                ALLOCATE( total_2d(-1:nx+1,-1:ny+1) )
130#endif
131             ENDIF
132          ENDIF
133
134       CASE ( 'xz' )
135
136          s = 2
137          ALLOCATE( local_2d(nxl-1:nxr+1,nzb:nzt+1) )
138
139!
140!--       Classic and 64bit offset NetCDF output is done only on PE0.
141!--       netCDF4/HDF5 output may be done in parallel on all PEs.
142          IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
143          THEN
144             CALL check_open( 102+av*10 )
145          ENDIF
146
147          IF ( data_output_2d_on_each_pe )  THEN
148             CALL check_open( 22 )
149          ELSE
150             IF ( myid == 0 )  THEN
151                IF ( iso2d_output )  CALL check_open( 22 )
152#if defined( __parallel )
153                ALLOCATE( total_2d(-1:nx+1,nzb:nzt+1) )
154#endif
155             ENDIF
156          ENDIF
157
158       CASE ( 'yz' )
159
160          s = 3
161          ALLOCATE( local_2d(nys-1:nyn+1,nzb:nzt+1) )
162
163!
164!--       Classic and 64bit offset NetCDF output is done only on PE0.
165!--       netCDF4/HDF5 output may be done in parallel on all PEs.
166          IF ( netcdf_output .AND. ( myid == 0 .OR. netcdf_data_format > 2 ) ) &
167          THEN
168             CALL check_open( 103+av*10 )
169          ENDIF
170
171          IF ( data_output_2d_on_each_pe )  THEN
172             CALL check_open( 23 )
173          ELSE
174             IF ( myid == 0 )  THEN
175                IF ( iso2d_output )  CALL check_open( 23 )
176#if defined( __parallel )
177                ALLOCATE( total_2d(-1:ny+1,nzb:nzt+1) )
178#endif
179             ENDIF
180          ENDIF
181
182       CASE DEFAULT
183
184          message_string = 'unknown cross-section: ' // TRIM( mode )
185          CALL message( 'data_output_2d', 'PA0180', 1, 2, 0, 6, 0 )
186
187    END SELECT
188
189!
190!-- Allocate a temporary array for resorting (kji -> ijk).
191    ALLOCATE( local_pf(nxl-1:nxr+1,nys-1:nyn+1,nzb:nzt+1) )
192
193!
194!-- Loop of all variables to be written.
195!-- Output dimensions chosen
196    if = 1
197    l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
198    do2d_mode = do2d(av,if)(l-1:l)
199
200    DO  WHILE ( do2d(av,if)(1:1) /= ' ' )
201
202       IF ( do2d_mode == mode )  THEN
203!
204!--       Store the array chosen on the temporary array.
205          resorted = .FALSE.
206          SELECT CASE ( TRIM( do2d(av,if) ) )
207
208             CASE ( 'e_xy', 'e_xz', 'e_yz' )
209                IF ( av == 0 )  THEN
210                   to_be_resorted => e
211                ELSE
212                   to_be_resorted => e_av
213                ENDIF
214                IF ( mode == 'xy' )  level_z = zu
215
216             CASE ( 'lwp*_xy' )        ! 2d-array
217                IF ( av == 0 )  THEN
218                   DO  i = nxl-1, nxr+1
219                      DO  j = nys-1, nyn+1
220                         local_pf(i,j,nzb+1) = SUM( ql(nzb:nzt,j,i) * &
221                                                    dzw(1:nzt+1) )
222                      ENDDO
223                   ENDDO
224                ELSE
225                   DO  i = nxl-1, nxr+1
226                      DO  j = nys-1, nyn+1
227                         local_pf(i,j,nzb+1) = lwp_av(j,i)
228                      ENDDO
229                   ENDDO
230                ENDIF
231                resorted = .TRUE.
232                two_d = .TRUE.
233                level_z(nzb+1) = zu(nzb+1)
234
235             CASE ( 'p_xy', 'p_xz', 'p_yz' )
236                IF ( av == 0 )  THEN
237                   to_be_resorted => p
238                ELSE
239                   to_be_resorted => p_av
240                ENDIF
241                IF ( mode == 'xy' )  level_z = zu
242
243             CASE ( 'pc_xy', 'pc_xz', 'pc_yz' )  ! particle concentration
244                IF ( av == 0 )  THEN
245                   IF ( simulated_time >= particle_advection_start )  THEN
246                      tend = prt_count
247                      CALL exchange_horiz( tend )
248                   ELSE
249                      tend = 0.0
250                   ENDIF
251                   DO  i = nxl-1, nxr+1
252                      DO  j = nys-1, nyn+1
253                         DO  k = nzb, nzt+1
254                            local_pf(i,j,k) = tend(k,j,i)
255                         ENDDO
256                      ENDDO
257                   ENDDO
258                   resorted = .TRUE.
259                ELSE
260                   CALL exchange_horiz( pc_av )
261                   to_be_resorted => pc_av
262                ENDIF
263
264             CASE ( 'pr_xy', 'pr_xz', 'pr_yz' )  ! mean particle radius
265                IF ( av == 0 )  THEN
266                   IF ( simulated_time >= particle_advection_start )  THEN
267                      DO  i = nxl, nxr
268                         DO  j = nys, nyn
269                            DO  k = nzb, nzt+1
270                               psi = prt_start_index(k,j,i)
271                               s_r3 = 0.0
272                               s_r4 = 0.0
273                               DO  n = psi, psi+prt_count(k,j,i)-1
274                                  s_r3 = s_r3 + particles(n)%radius**3
275                                  s_r4 = s_r4 + particles(n)%radius**4
276                               ENDDO
277                               IF ( s_r3 /= 0.0 )  THEN
278                                  mean_r = s_r4 / s_r3
279                               ELSE
280                                  mean_r = 0.0
281                               ENDIF
282                               tend(k,j,i) = mean_r
283                            ENDDO
284                         ENDDO
285                      ENDDO
286                      CALL exchange_horiz( tend )
287                   ELSE
288                      tend = 0.0
289                   ENDIF
290                   DO  i = nxl-1, nxr+1
291                      DO  j = nys-1, nyn+1
292                         DO  k = nzb, nzt+1
293                            local_pf(i,j,k) = tend(k,j,i)
294                         ENDDO
295                      ENDDO
296                   ENDDO
297                   resorted = .TRUE.
298                ELSE
299                   CALL exchange_horiz( pr_av )
300                   to_be_resorted => pr_av
301                ENDIF
302
303             CASE ( 'pra*_xy' )        ! 2d-array / integral quantity => no av
304                CALL exchange_horiz_2d( precipitation_amount )
305                DO  i = nxl-1, nxr+1
306                   DO  j = nys-1, nyn+1
307                      local_pf(i,j,nzb+1) =  precipitation_amount(j,i)
308                   ENDDO
309                ENDDO
310                precipitation_amount = 0.0   ! reset for next integ. interval
311                resorted = .TRUE.
312                two_d = .TRUE.
313                level_z(nzb+1) = zu(nzb+1)
314
315             CASE ( 'prr*_xy' )        ! 2d-array
316                IF ( av == 0 )  THEN
317                   CALL exchange_horiz_2d( precipitation_rate )
318                   DO  i = nxl-1, nxr+1
319                      DO  j = nys-1, nyn+1 
320                         local_pf(i,j,nzb+1) =  precipitation_rate(j,i)
321                      ENDDO
322                   ENDDO
323                ELSE
324                   CALL exchange_horiz_2d( precipitation_rate_av )
325                   DO  i = nxl-1, nxr+1
326                      DO  j = nys-1, nyn+1 
327                         local_pf(i,j,nzb+1) =  precipitation_rate_av(j,i)
328                      ENDDO
329                   ENDDO
330                ENDIF
331                resorted = .TRUE.
332                two_d = .TRUE.
333                level_z(nzb+1) = zu(nzb+1)
334
335             CASE ( 'pt_xy', 'pt_xz', 'pt_yz' )
336                IF ( av == 0 )  THEN
337                   IF ( .NOT. cloud_physics ) THEN
338                      to_be_resorted => pt
339                   ELSE
340                      DO  i = nxl-1, nxr+1
341                         DO  j = nys-1, nyn+1
342                            DO  k = nzb, nzt+1
343                               local_pf(i,j,k) = pt(k,j,i) + l_d_cp *    &
344                                                             pt_d_t(k) * &
345                                                             ql(k,j,i)
346                            ENDDO
347                         ENDDO
348                      ENDDO
349                      resorted = .TRUE.
350                   ENDIF
351                ELSE
352                   to_be_resorted => pt_av
353                ENDIF
354                IF ( mode == 'xy' )  level_z = zu
355
356             CASE ( 'q_xy', 'q_xz', 'q_yz' )
357                IF ( av == 0 )  THEN
358                   to_be_resorted => q
359                ELSE
360                   to_be_resorted => q_av
361                ENDIF
362                IF ( mode == 'xy' )  level_z = zu
363
364             CASE ( 'ql_xy', 'ql_xz', 'ql_yz' )
365                IF ( av == 0 )  THEN
366                   to_be_resorted => ql
367                ELSE
368                   to_be_resorted => ql_av
369                ENDIF
370                IF ( mode == 'xy' )  level_z = zu
371
372             CASE ( 'ql_c_xy', 'ql_c_xz', 'ql_c_yz' )
373                IF ( av == 0 )  THEN
374                   to_be_resorted => ql_c
375                ELSE
376                   to_be_resorted => ql_c_av
377                ENDIF
378                IF ( mode == 'xy' )  level_z = zu
379
380             CASE ( 'ql_v_xy', 'ql_v_xz', 'ql_v_yz' )
381                IF ( av == 0 )  THEN
382                   to_be_resorted => ql_v
383                ELSE
384                   to_be_resorted => ql_v_av
385                ENDIF
386                IF ( mode == 'xy' )  level_z = zu
387
388             CASE ( 'ql_vp_xy', 'ql_vp_xz', 'ql_vp_yz' )
389                IF ( av == 0 )  THEN
390                   to_be_resorted => ql_vp
391                ELSE
392                   to_be_resorted => ql_vp_av
393                ENDIF
394                IF ( mode == 'xy' )  level_z = zu
395
396             CASE ( 'qsws*_xy' )        ! 2d-array
397                IF ( av == 0 ) THEN
398                   DO  i = nxl-1, nxr+1
399                      DO  j = nys-1, nyn+1 
400                         local_pf(i,j,nzb+1) =  qsws(j,i)
401                      ENDDO
402                   ENDDO
403                ELSE
404                   DO  i = nxl-1, nxr+1
405                      DO  j = nys-1, nyn+1 
406                         local_pf(i,j,nzb+1) =  qsws_av(j,i)
407                      ENDDO
408                   ENDDO
409                ENDIF
410                resorted = .TRUE.
411                two_d = .TRUE.
412                level_z(nzb+1) = zu(nzb+1)
413
414             CASE ( 'qv_xy', 'qv_xz', 'qv_yz' )
415                IF ( av == 0 )  THEN
416                   DO  i = nxl-1, nxr+1
417                      DO  j = nys-1, nyn+1
418                         DO  k = nzb, nzt+1
419                            local_pf(i,j,k) = q(k,j,i) - ql(k,j,i)
420                         ENDDO
421                      ENDDO
422                   ENDDO
423                   resorted = .TRUE.
424                ELSE
425                   to_be_resorted => qv_av
426                ENDIF
427                IF ( mode == 'xy' )  level_z = zu
428
429             CASE ( 'rho_xy', 'rho_xz', 'rho_yz' )
430                IF ( av == 0 )  THEN
431                   to_be_resorted => rho
432                ELSE
433                   to_be_resorted => rho_av
434                ENDIF
435
436             CASE ( 's_xy', 's_xz', 's_yz' )
437                IF ( av == 0 )  THEN
438                   to_be_resorted => q
439                ELSE
440                   to_be_resorted => s_av
441                ENDIF
442
443             CASE ( 'sa_xy', 'sa_xz', 'sa_yz' )
444                IF ( av == 0 )  THEN
445                   to_be_resorted => sa
446                ELSE
447                   to_be_resorted => sa_av
448                ENDIF
449
450             CASE ( 'shf*_xy' )        ! 2d-array
451                IF ( av == 0 ) THEN
452                   DO  i = nxl-1, nxr+1
453                      DO  j = nys-1, nyn+1 
454                         local_pf(i,j,nzb+1) =  shf(j,i)
455                      ENDDO
456                   ENDDO
457                ELSE
458                   DO  i = nxl-1, nxr+1
459                      DO  j = nys-1, nyn+1 
460                         local_pf(i,j,nzb+1) =  shf_av(j,i)
461                      ENDDO
462                   ENDDO
463                ENDIF
464                resorted = .TRUE.
465                two_d = .TRUE.
466                level_z(nzb+1) = zu(nzb+1)
467
468             CASE ( 't*_xy' )        ! 2d-array
469                IF ( av == 0 )  THEN
470                   DO  i = nxl-1, nxr+1
471                      DO  j = nys-1, nyn+1
472                         local_pf(i,j,nzb+1) = ts(j,i)
473                      ENDDO
474                   ENDDO
475                ELSE
476                   DO  i = nxl-1, nxr+1
477                      DO  j = nys-1, nyn+1
478                         local_pf(i,j,nzb+1) = ts_av(j,i)
479                      ENDDO
480                   ENDDO
481                ENDIF
482                resorted = .TRUE.
483                two_d = .TRUE.
484                level_z(nzb+1) = zu(nzb+1)
485
486             CASE ( 'u_xy', 'u_xz', 'u_yz' )
487                IF ( av == 0 )  THEN
488                   to_be_resorted => u
489                ELSE
490                   to_be_resorted => u_av
491                ENDIF
492                IF ( mode == 'xy' )  level_z = zu
493!
494!--             Substitute the values generated by "mirror" boundary condition
495!--             at the bottom boundary by the real surface values.
496                IF ( do2d(av,if) == 'u_xz'  .OR.  do2d(av,if) == 'u_yz' )  THEN
497                   IF ( ibc_uv_b == 0 )  local_pf(:,:,nzb) = 0.0
498                ENDIF
499
500             CASE ( 'u*_xy' )        ! 2d-array
501                IF ( av == 0 )  THEN
502                   DO  i = nxl-1, nxr+1
503                      DO  j = nys-1, nyn+1
504                         local_pf(i,j,nzb+1) = us(j,i)
505                      ENDDO
506                   ENDDO
507                ELSE
508                   DO  i = nxl-1, nxr+1
509                      DO  j = nys-1, nyn+1
510                         local_pf(i,j,nzb+1) = us_av(j,i)
511                      ENDDO
512                   ENDDO
513                ENDIF
514                resorted = .TRUE.
515                two_d = .TRUE.
516                level_z(nzb+1) = zu(nzb+1)
517
518             CASE ( 'v_xy', 'v_xz', 'v_yz' )
519                IF ( av == 0 )  THEN
520                   to_be_resorted => v
521                ELSE
522                   to_be_resorted => v_av
523                ENDIF
524                IF ( mode == 'xy' )  level_z = zu
525!
526!--             Substitute the values generated by "mirror" boundary condition
527!--             at the bottom boundary by the real surface values.
528                IF ( do2d(av,if) == 'v_xz'  .OR.  do2d(av,if) == 'v_yz' )  THEN
529                   IF ( ibc_uv_b == 0 )  local_pf(:,:,nzb) = 0.0
530                ENDIF
531
532             CASE ( 'vpt_xy', 'vpt_xz', 'vpt_yz' )
533                IF ( av == 0 )  THEN
534                   to_be_resorted => vpt
535                ELSE
536                   to_be_resorted => vpt_av
537                ENDIF
538                IF ( mode == 'xy' )  level_z = zu
539
540             CASE ( 'w_xy', 'w_xz', 'w_yz' )
541                IF ( av == 0 )  THEN
542                   to_be_resorted => w
543                ELSE
544                   to_be_resorted => w_av
545                ENDIF
546                IF ( mode == 'xy' )  level_z = zw
547
548             CASE ( 'z0*_xy' )        ! 2d-array
549                IF ( av == 0 ) THEN
550                   DO  i = nxl-1, nxr+1
551                      DO  j = nys-1, nyn+1 
552                         local_pf(i,j,nzb+1) =  z0(j,i)
553                      ENDDO
554                   ENDDO
555                ELSE
556                   DO  i = nxl-1, nxr+1
557                      DO  j = nys-1, nyn+1 
558                         local_pf(i,j,nzb+1) =  z0_av(j,i)
559                      ENDDO
560                   ENDDO
561                ENDIF
562                resorted = .TRUE.
563                two_d = .TRUE.
564                level_z(nzb+1) = zu(nzb+1)
565
566             CASE DEFAULT
567!
568!--             User defined quantity
569                CALL user_data_output_2d( av, do2d(av,if), found, grid, &
570                                          local_pf, two_d )
571                resorted = .TRUE.
572
573                IF ( grid == 'zu' )  THEN
574                   IF ( mode == 'xy' )  level_z = zu
575                ELSEIF ( grid == 'zw' )  THEN
576                   IF ( mode == 'xy' )  level_z = zw
577                ELSEIF ( grid == 'zu1' ) THEN
578                   IF ( mode == 'xy' )  level_z(nzb+1) = zu(nzb+1)
579                ENDIF
580
581                IF ( .NOT. found )  THEN
582                   message_string = 'no output provided for: ' //    &
583                                    TRIM( do2d(av,if) )
584                   CALL message( 'data_output_2d', 'PA0181', 0, 0, 0, 6, 0 )
585                ENDIF
586
587          END SELECT
588
589!
590!--       Resort the array to be output, if not done above
591          IF ( .NOT. resorted )  THEN
592             DO  i = nxl-1, nxr+1
593                DO  j = nys-1, nyn+1
594                   DO  k = nzb, nzt+1
595                      local_pf(i,j,k) = to_be_resorted(k,j,i)
596                   ENDDO
597                ENDDO
598             ENDDO
599          ENDIF
600
601!
602!--       Output of the individual cross-sections, depending on the cross-
603!--       section mode chosen.
604          is = 1
605   loop1: DO  WHILE ( section(is,s) /= -9999  .OR.  two_d )
606
607             SELECT CASE ( mode )
608
609                CASE ( 'xy' )
610!
611!--                Determine the cross section index
612                   IF ( two_d )  THEN
613                      layer_xy = nzb+1
614                   ELSE
615                      layer_xy = section(is,s)
616                   ENDIF
617
618!
619!--                Update the NetCDF xy cross section time axis
620                   IF ( myid == 0  .OR.  netcdf_data_format > 2 )  THEN
621                      IF ( simulated_time /= do2d_xy_last_time(av) )  THEN
622                         do2d_xy_time_count(av) = do2d_xy_time_count(av) + 1
623                         do2d_xy_last_time(av)  = simulated_time
624                         IF ( ( .NOT. data_output_2d_on_each_pe  .AND. &
625                              netcdf_output )  .OR.  netcdf_data_format > 2 ) &
626                         THEN
627#if defined( __netcdf )
628                            nc_stat = NF90_PUT_VAR( id_set_xy(av),             &
629                                                    id_var_time_xy(av),        &
630                                             (/ time_since_reference_point /), &
631                                         start = (/ do2d_xy_time_count(av) /), &
632                                                    count = (/ 1 /) )
633                            CALL handle_netcdf_error( 'data_output_2d', 53 )
634#endif
635                         ENDIF
636                      ENDIF
637                   ENDIF
638!
639!--                If required, carry out averaging along z
640                   IF ( section(is,s) == -1  .AND.  .NOT. two_d )  THEN
641
642                      local_2d = 0.0
643!
644!--                   Carry out the averaging (all data are on the PE)
645                      DO  k = nzb, nzt+1
646                         DO  j = nys-1, nyn+1
647                            DO  i = nxl-1, nxr+1
648                               local_2d(i,j) = local_2d(i,j) + local_pf(i,j,k)
649                            ENDDO
650                         ENDDO
651                      ENDDO
652
653                      local_2d = local_2d / ( nzt -nzb + 2.0 )
654
655                   ELSE
656!
657!--                   Just store the respective section on the local array
658                      local_2d = local_pf(:,:,layer_xy)
659
660                   ENDIF
661
662#if defined( __parallel )
663                   IF ( netcdf_output  .AND.  netcdf_data_format > 2 )  THEN
664!
665!--                   Output in NetCDF4/HDF5 format.
666!--                   Do not output redundant ghost point data except for the
667!--                   boundaries of the total domain.
668                      IF ( two_d ) THEN
669                         iis = 1
670                      ELSE
671                         iis = is
672                      ENDIF
673
674#if defined( __netcdf )
675                      IF ( nxr == nx  .AND.  nyn /= ny )  THEN
676                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
677                                                 id_var_do2d(av,if),           &
678                                                 local_2d(nxl:nxr+1,nys:nyn),  &
679                                                 start = (/ nxl+1, nys+1, iis, &
680                                                    do2d_xy_time_count(av) /), &
681                                                 count = (/ nxr-nxl+2,         &
682                                                            nyn-nys+1, 1, 1 /) )
683                      ELSEIF ( nxr /= nx  .AND.  nyn == ny )  THEN
684                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
685                                                 id_var_do2d(av,if),           &
686                                                 local_2d(nxl:nxr,nys:nyn+1),  &
687                                                 start = (/ nxl+1, nys+1, iis, &
688                                                    do2d_xy_time_count(av) /), &
689                                                 count = (/ nxr-nxl+1,         &
690                                                            nyn-nys+2, 1, 1 /) )
691                      ELSEIF ( nxr == nx  .AND.  nyn == ny )  THEN
692                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
693                                                 id_var_do2d(av,if),           &
694                                                 local_2d(nxl:nxr+1,nys:nyn+1),&
695                                                 start = (/ nxl+1, nys+1, iis, &
696                                                    do2d_xy_time_count(av) /), &
697                                                 count = (/ nxr-nxl+2,          &
698                                                            nyn-nys+2, 1, 1 /) )
699                      ELSE
700                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
701                                                 id_var_do2d(av,if),           &
702                                                 local_2d(nxl:nxr,nys:nyn),    &
703                                                 start = (/ nxl+1, nys+1, iis, &
704                                                    do2d_xy_time_count(av) /), &
705                                                 count = (/ nxr-nxl+1,         &
706                                                            nyn-nys+1, 1, 1 /) )
707                      ENDIF
708
709                      CALL handle_netcdf_error( 'data_output_2d', 55 )
710#endif
711                   ELSE
712
713                      IF ( data_output_2d_on_each_pe )  THEN
714!
715!--                      Output of partial arrays on each PE
716#if defined( __netcdf )
717                         IF ( netcdf_output  .AND.  myid == 0 )  THEN
718                            WRITE ( 21 )  simulated_time, &
719                                          do2d_xy_time_count(av), av
720                         ENDIF
721#endif
722                         WRITE ( 21 )  nxl-1, nxr+1, nys-1, nyn+1
723                         WRITE ( 21 )  local_2d
724
725                      ELSE
726!
727!--                      PE0 receives partial arrays from all processors and
728!--                      then outputs them. Here a barrier has to be set,
729!--                      because otherwise "-MPI- FATAL: Remote protocol queue
730!--                      full" may occur.
731                         CALL MPI_BARRIER( comm2d, ierr )
732
733                         ngp = ( nxr-nxl+3 ) * ( nyn-nys+3 )
734                         IF ( myid == 0 )  THEN
735!
736!--                         Local array can be relocated directly.
737                            total_2d(nxl-1:nxr+1,nys-1:nyn+1) = local_2d
738!
739!--                         Receive data from all other PEs.
740                            DO  n = 1, numprocs-1
741!
742!--                            Receive index limits first, then array.
743!--                            Index limits are received in arbitrary order from
744!--                            the PEs.
745                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,    &
746                                              MPI_ANY_SOURCE, 0, comm2d, &
747                                              status, ierr )
748                               sender = status(MPI_SOURCE)
749                               DEALLOCATE( local_2d )
750                               ALLOCATE( local_2d(ind(1):ind(2),ind(3):ind(4)) )
751                               CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp,  &
752                                              MPI_REAL, sender, 1, comm2d,   &
753                                              status, ierr )
754                               total_2d(ind(1):ind(2),ind(3):ind(4)) = local_2d
755                            ENDDO
756!
757!--                         Output of the total cross-section.
758                            IF ( iso2d_output )  THEN
759                               WRITE (21)  total_2d(0:nx+1,0:ny+1)
760                            ENDIF
761!
762!--                         Relocate the local array for the next loop increment
763                            DEALLOCATE( local_2d )
764                            ALLOCATE( local_2d(nxl-1:nxr+1,nys-1:nyn+1) )
765
766#if defined( __netcdf )
767                            IF ( netcdf_output )  THEN
768                               IF ( two_d ) THEN
769                                  nc_stat = NF90_PUT_VAR( id_set_xy(av),       &
770                                                          id_var_do2d(av,if),  &
771                                                      total_2d(0:nx+1,0:ny+1), &
772                                start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
773                                                count = (/ nx+2, ny+2, 1, 1 /) )
774                               ELSE
775                                  nc_stat = NF90_PUT_VAR( id_set_xy(av),       &
776                                                          id_var_do2d(av,if),  &
777                                                      total_2d(0:nx+1,0:ny+1), &
778                               start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
779                                                count = (/ nx+2, ny+2, 1, 1 /) )
780                               ENDIF
781                               CALL handle_netcdf_error( 'data_output_2d', 54 )
782                            ENDIF
783#endif
784
785                         ELSE
786!
787!--                         First send the local index limits to PE0
788                            ind(1) = nxl-1; ind(2) = nxr+1
789                            ind(3) = nys-1; ind(4) = nyn+1
790                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
791                                           comm2d, ierr )
792!
793!--                         Send data to PE0
794                            CALL MPI_SEND( local_2d(nxl-1,nys-1), ngp, &
795                                           MPI_REAL, 0, 1, comm2d, ierr )
796                         ENDIF
797!
798!--                      A barrier has to be set, because otherwise some PEs may
799!--                      proceed too fast so that PE0 may receive wrong data on
800!--                      tag 0
801                         CALL MPI_BARRIER( comm2d, ierr )
802                      ENDIF
803
804                   ENDIF
805#else
806                   IF ( iso2d_output )  THEN
807                      WRITE (21)  local_2d(nxl:nxr+1,nys:nyn+1)
808                   ENDIF
809#if defined( __netcdf )
810                   IF ( netcdf_output )  THEN
811                      IF ( two_d ) THEN
812                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
813                                                 id_var_do2d(av,if),           &
814                                                local_2d(nxl:nxr+1,nys:nyn+1), &
815                                start = (/ 1, 1, 1, do2d_xy_time_count(av) /), &
816                                              count = (/ nx+2, ny+2, 1, 1 /) )
817                      ELSE
818                         nc_stat = NF90_PUT_VAR( id_set_xy(av),                &
819                                                 id_var_do2d(av,if),           &
820                                                local_2d(nxl:nxr+1,nys:nyn+1), &
821                               start = (/ 1, 1, is, do2d_xy_time_count(av) /), &
822                                              count = (/ nx+2, ny+2, 1, 1 /) )
823                      ENDIF
824                      CALL handle_netcdf_error( 'data_output_2d', 447 )
825                   ENDIF
826#endif
827#endif
828                   do2d_xy_n = do2d_xy_n + 1
829!
830!--                Write LOCAL parameter set for ISO2D.
831                   IF ( myid == 0  .AND.  iso2d_output )  THEN
832                      IF ( section(is,s) /= -1 )  THEN
833                         WRITE ( section_chr, '(''z = '',F7.2,'' m  (GP '',I3, &
834                                               &'')'')'                        &
835                               )  level_z(layer_xy), layer_xy
836                      ELSE
837                         section_chr = 'averaged along z'
838                      ENDIF
839                      IF ( av == 0 )  THEN
840                         rtext = TRIM( do2d(av,if) ) // '  t = ' //    &
841                                 TRIM( simulated_time_chr ) // '  ' // &
842                                 TRIM( section_chr )
843                      ELSE
844                         rtext = TRIM( do2d(av,if) ) // '  averaged t = ' // &
845                                 TRIM( simulated_time_chr ) // '  ' //       &
846                                 TRIM( section_chr )
847                      ENDIF
848                      WRITE (27,LOCAL)
849                   ENDIF
850!
851!--                For 2D-arrays (e.g. u*) only one cross-section is available.
852!--                Hence exit loop of output levels.
853                   IF ( two_d )  THEN
854                      two_d = .FALSE.
855                      EXIT loop1
856                   ENDIF
857
858                CASE ( 'xz' )
859!
860!--                Update the NetCDF xz cross section time axis
861                   IF ( myid == 0  .OR.  netcdf_data_format > 2 )  THEN
862
863                      IF ( simulated_time /= do2d_xz_last_time(av) )  THEN
864                         do2d_xz_time_count(av) = do2d_xz_time_count(av) + 1
865                         do2d_xz_last_time(av)  = simulated_time
866                         IF ( ( .NOT. data_output_2d_on_each_pe  .AND.        &
867                              netcdf_output )  .OR.  netcdf_data_format > 2 ) &
868                         THEN
869#if defined( __netcdf )
870                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
871                                                    id_var_time_xz(av),        &
872                                             (/ time_since_reference_point /), &
873                                         start = (/ do2d_xz_time_count(av) /), &
874                                                    count = (/ 1 /) )
875                            CALL handle_netcdf_error( 'data_output_2d', 56 )
876#endif
877                         ENDIF
878                      ENDIF
879
880                   ENDIF
881!
882!--                If required, carry out averaging along y
883                   IF ( section(is,s) == -1 )  THEN
884
885                      ALLOCATE( local_2d_l(nxl-1:nxr+1,nzb:nzt+1) )
886                      local_2d_l = 0.0
887                      ngp = ( nxr-nxl+3 ) * ( nzt-nzb+2 )
888!
889!--                   First local averaging on the PE
890                      DO  k = nzb, nzt+1
891                         DO  j = nys, nyn
892                            DO  i = nxl-1, nxr+1
893                               local_2d_l(i,k) = local_2d_l(i,k) + &
894                                                 local_pf(i,j,k)
895                            ENDDO
896                         ENDDO
897                      ENDDO
898#if defined( __parallel )
899!
900!--                   Now do the averaging over all PEs along y
901                      IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
902                      CALL MPI_ALLREDUCE( local_2d_l(nxl-1,nzb),              &
903                                          local_2d(nxl-1,nzb), ngp, MPI_REAL, &
904                                          MPI_SUM, comm1dy, ierr )
905#else
906                      local_2d = local_2d_l
907#endif
908                      local_2d = local_2d / ( ny + 1.0 )
909
910                      DEALLOCATE( local_2d_l )
911
912                   ELSE
913!
914!--                   Just store the respective section on the local array
915!--                   (but only if it is available on this PE!)
916                      IF ( section(is,s) >= nys  .AND.  section(is,s) <= nyn ) &
917                      THEN
918                         local_2d = local_pf(:,section(is,s),nzb:nzt+1)
919                      ENDIF
920
921                   ENDIF
922
923#if defined( __parallel )
924                   IF ( netcdf_output  .AND.  netcdf_data_format > 2 )  THEN
925!
926!--                   ATTENTION: The following lines are a workaround, because
927!--                              independet output does not work with the
928!--                              current NetCDF4 installation. Therefore, data
929!--                              are transferred from PEs having the cross
930!--                              sections to other PEs along y having no cross
931!--                              section data. Some of these data are the
932!--                              output.
933!--                   BEGIN WORKAROUND---------------------------------------
934                      IF ( npey /= 1  .AND.  section(is,s) /= -1)  THEN
935                         ALLOCATE( local_2d_l(nxl-1:nxr+1,nzb:nzt+1) )
936                         local_2d_l = 0.0
937                         IF ( section(is,s) >= nys .AND. section(is,s) <= nyn )&
938                         THEN
939                            local_2d_l = local_2d
940                         ENDIF
941#if defined( __parallel )
942!
943!--                      Distribute data over all PEs along y
944                         ngp = ( nxr-nxl+3 ) * ( nzt-nzb+2 )
945                         IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
946                         CALL MPI_ALLREDUCE( local_2d_l(nxl-1,nzb),            &
947                                             local_2d(nxl-1,nzb), ngp,         &
948                                             MPI_REAL, MPI_SUM, comm1dy, ierr )
949#else
950                         local_2d = local_2d_l
951#endif
952                         DEALLOCATE( local_2d_l )
953                      ENDIF
954!--                   END WORKAROUND-----------------------------------------
955
956!
957!--                   Output in NetCDF4/HDF5 format.
958!--                   Output only on those PEs where the respective cross
959!--                   sections reside. Cross sections averaged along y are
960!--                   output on the respective first PE along y (myidy=0).
961                      IF ( ( section(is,s) >= nys  .AND.  &
962                             section(is,s) <= nyn )  .OR.  &
963                           ( section(is,s) == -1  .AND.  myidy == 0 ) )  THEN
964!
965!--                      Do not output redundant ghost point data except for the
966!--                      boundaries of the total domain.
967#if defined( __netcdf )
968                         IF ( nxr == nx )  THEN
969                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
970                                                id_var_do2d(av,if),            &
971                                                local_2d(nxl:nxr+1,nzb:nzt+1), &
972                                                start = (/ nxl+1, is, 1,       &
973                                                    do2d_xz_time_count(av) /), &
974                                                count = (/ nxr-nxl+2, 1,       &
975                                                           nzt+2, 1 /) )
976                         ELSE
977                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
978                                                id_var_do2d(av,if),            &
979                                                local_2d(nxl:nxr,nzb:nzt+1),   &
980                                                start = (/ nxl+1, is, 1,       &
981                                                    do2d_xz_time_count(av) /), &
982                                                count = (/ nxr-nxl+1, 1,       &
983                                                           nzt+2, 1 /) )
984                         ENDIF
985
986                         CALL handle_netcdf_error( 'data_output_2d', 57 )
987
988                      ELSE
989!
990!--                      Output on other PEs. Only one point is output!!
991!--                      ATTENTION: This is a workaround (see above)!!
992                         IF ( npey /= 1 )  THEN
993                            nc_stat = NF90_PUT_VAR( id_set_xz(av),             &
994                                                    id_var_do2d(av,if),        &
995                                                    local_2d(nxl:nxl,nzb:nzb), &
996                                                    start = (/ nxl+1, is, 1,   &
997                                                    do2d_xz_time_count(av) /), &
998                                                    count = (/ 1, 1, 1, 1 /) )
999                            CALL handle_netcdf_error( 'data_output_2d', 451 )
1000                         ENDIF
1001#endif
1002                      ENDIF
1003
1004                   ELSE
1005
1006                      IF ( data_output_2d_on_each_pe )  THEN
1007!
1008!--                      Output of partial arrays on each PE. If the cross
1009!--                      section does not reside on the PE, output special
1010!--                      index values.
1011#if defined( __netcdf )
1012                         IF ( netcdf_output  .AND.  myid == 0 )  THEN
1013                            WRITE ( 22 )  simulated_time, &
1014                                          do2d_xz_time_count(av), av
1015                         ENDIF
1016#endif
1017                         IF ( ( section(is,s) >= nys  .AND.                  &
1018                                section(is,s) <= nyn )  .OR.                 &
1019                              ( section(is,s) == -1  .AND.  nys-1 == -1 ) )  &
1020                         THEN
1021                            WRITE (22)  nxl-1, nxr+1, nzb, nzt+1
1022                            WRITE (22)  local_2d
1023                         ELSE
1024                            WRITE (22)  -1, -1, -1, -1
1025                         ENDIF
1026
1027                      ELSE
1028!
1029!--                      PE0 receives partial arrays from all processors of the
1030!--                      respective cross section and outputs them. Here a
1031!--                      barrier has to be set, because otherwise
1032!--                      "-MPI- FATAL: Remote protocol queue full" may occur.
1033                         CALL MPI_BARRIER( comm2d, ierr )
1034
1035                         ngp = ( nxr-nxl+3 ) * ( nzt-nzb+2 )
1036                         IF ( myid == 0 )  THEN
1037!
1038!--                         Local array can be relocated directly.
1039                            IF ( ( section(is,s) >= nys  .AND.                 &
1040                                   section(is,s) <= nyn )  .OR.                &
1041                                 ( section(is,s) == -1  .AND.  nys-1 == -1 ) ) &
1042                            THEN
1043                               total_2d(nxl-1:nxr+1,nzb:nzt+1) = local_2d
1044                            ENDIF
1045!
1046!--                         Receive data from all other PEs.
1047                            DO  n = 1, numprocs-1
1048!
1049!--                            Receive index limits first, then array.
1050!--                            Index limits are received in arbitrary order from
1051!--                            the PEs.
1052                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,     &
1053                                              MPI_ANY_SOURCE, 0, comm2d,  &
1054                                              status, ierr )
1055!
1056!--                            Not all PEs have data for XZ-cross-section.
1057                               IF ( ind(1) /= -9999 )  THEN
1058                                  sender = status(MPI_SOURCE)
1059                                  DEALLOCATE( local_2d )
1060                                  ALLOCATE( local_2d(ind(1):ind(2), &
1061                                                     ind(3):ind(4)) )
1062                                  CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1063                                                 MPI_REAL, sender, 1, comm2d,  &
1064                                                 status, ierr )
1065                                  total_2d(ind(1):ind(2),ind(3):ind(4)) = &
1066                                                                        local_2d
1067                               ENDIF
1068                            ENDDO
1069!
1070!--                         Output of the total cross-section.
1071                            IF ( iso2d_output )  THEN
1072                               WRITE (22)  total_2d(0:nx+1,nzb:nzt+1)
1073                            ENDIF
1074!
1075!--                         Relocate the local array for the next loop increment
1076                            DEALLOCATE( local_2d )
1077                            ALLOCATE( local_2d(nxl-1:nxr+1,nzb:nzt+1) )
1078
1079#if defined( __netcdf )
1080                            IF ( netcdf_output )  THEN
1081                               nc_stat = NF90_PUT_VAR( id_set_xz(av),          &
1082                                                    id_var_do2d(av,if),        &
1083                                                    total_2d(0:nx+1,nzb:nzt+1),&
1084                               start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1085                                                count = (/ nx+2, 1, nz+2, 1 /) )
1086                               CALL handle_netcdf_error( 'data_output_2d', 58 )
1087                            ENDIF
1088#endif
1089
1090                         ELSE
1091!
1092!--                         If the cross section resides on the PE, send the
1093!--                         local index limits, otherwise send -9999 to PE0.
1094                            IF ( ( section(is,s) >= nys  .AND.                 &
1095                                   section(is,s) <= nyn )  .OR.                &
1096                                 ( section(is,s) == -1  .AND.  nys-1 == -1 ) ) &
1097                            THEN
1098                               ind(1) = nxl-1; ind(2) = nxr+1
1099                               ind(3) = nzb;   ind(4) = nzt+1
1100                            ELSE
1101                               ind(1) = -9999; ind(2) = -9999
1102                               ind(3) = -9999; ind(4) = -9999
1103                            ENDIF
1104                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1105                                           comm2d, ierr )
1106!
1107!--                         If applicable, send data to PE0.
1108                            IF ( ind(1) /= -9999 )  THEN
1109                               CALL MPI_SEND( local_2d(nxl-1,nzb), ngp, &
1110                                              MPI_REAL, 0, 1, comm2d, ierr )
1111                            ENDIF
1112                         ENDIF
1113!
1114!--                      A barrier has to be set, because otherwise some PEs may
1115!--                      proceed too fast so that PE0 may receive wrong data on
1116!--                      tag 0
1117                         CALL MPI_BARRIER( comm2d, ierr )
1118                      ENDIF
1119
1120                   ENDIF
1121#else
1122                   IF ( iso2d_output )  THEN
1123                      WRITE (22)  local_2d(nxl:nxr+1,nzb:nzt+1)
1124                   ENDIF
1125#if defined( __netcdf )
1126                   IF ( netcdf_output )  THEN
1127                      nc_stat = NF90_PUT_VAR( id_set_xz(av),                   &
1128                                              id_var_do2d(av,if),              &
1129                                              local_2d(nxl:nxr+1,nzb:nzt+1),   &
1130                               start = (/ 1, is, 1, do2d_xz_time_count(av) /), &
1131                                              count = (/ nx+2, 1, nz+2, 1 /) )
1132                      CALL handle_netcdf_error( 'data_output_2d', 451 )
1133                   ENDIF
1134#endif
1135#endif
1136                   do2d_xz_n = do2d_xz_n + 1
1137!
1138!--                Write LOCAL-parameter set for ISO2D.
1139                   IF ( myid == 0  .AND.  iso2d_output )  THEN
1140                      IF ( section(is,s) /= -1 )  THEN
1141                         WRITE ( section_chr, '(''y = '',F8.2,'' m  (GP '',I3, &
1142                                               &'')'')'                        &
1143                               )  section(is,s)*dy, section(is,s)
1144                      ELSE
1145                         section_chr = 'averaged along y'
1146                      ENDIF
1147                      IF ( av == 0 )  THEN
1148                         rtext = TRIM( do2d(av,if) ) // '  t = ' //    &
1149                                 TRIM( simulated_time_chr ) // '  ' // &
1150                                 TRIM( section_chr )
1151                      ELSE
1152                         rtext = TRIM( do2d(av,if) ) // '  averaged t = ' // &
1153                                 TRIM( simulated_time_chr ) // '  ' //       &
1154                                 TRIM( section_chr )
1155                      ENDIF
1156                      WRITE (28,LOCAL)
1157                   ENDIF
1158
1159                CASE ( 'yz' )
1160!
1161!--                Update the NetCDF yz cross section time axis
1162                   IF ( myid == 0  .OR.  netcdf_data_format > 2 )  THEN
1163
1164                      IF ( simulated_time /= do2d_yz_last_time(av) )  THEN
1165                         do2d_yz_time_count(av) = do2d_yz_time_count(av) + 1
1166                         do2d_yz_last_time(av)  = simulated_time
1167                         IF ( ( .NOT. data_output_2d_on_each_pe  .AND.        &
1168                              netcdf_output )  .OR.  netcdf_data_format > 2 ) &
1169                         THEN
1170#if defined( __netcdf )
1171                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1172                                                    id_var_time_yz(av),        &
1173                                             (/ time_since_reference_point /), &
1174                                         start = (/ do2d_yz_time_count(av) /), &
1175                                                    count = (/ 1 /) )
1176                            CALL handle_netcdf_error( 'data_output_2d', 59 )
1177#endif
1178                         ENDIF
1179                      ENDIF
1180
1181                   ENDIF
1182!
1183!--                If required, carry out averaging along x
1184                   IF ( section(is,s) == -1 )  THEN
1185
1186                      ALLOCATE( local_2d_l(nys-1:nyn+1,nzb:nzt+1) )
1187                      local_2d_l = 0.0
1188                      ngp = ( nyn-nys+3 ) * ( nzt-nzb+2 )
1189!
1190!--                   First local averaging on the PE
1191                      DO  k = nzb, nzt+1
1192                         DO  j = nys-1, nyn+1
1193                            DO  i = nxl, nxr
1194                               local_2d_l(j,k) = local_2d_l(j,k) + &
1195                                                 local_pf(i,j,k)
1196                            ENDDO
1197                         ENDDO
1198                      ENDDO
1199#if defined( __parallel )
1200!
1201!--                   Now do the averaging over all PEs along x
1202                      IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
1203                      CALL MPI_ALLREDUCE( local_2d_l(nys-1,nzb),              &
1204                                          local_2d(nys-1,nzb), ngp, MPI_REAL, &
1205                                          MPI_SUM, comm1dx, ierr )
1206#else
1207                      local_2d = local_2d_l
1208#endif
1209                      local_2d = local_2d / ( nx + 1.0 )
1210
1211                      DEALLOCATE( local_2d_l )
1212
1213                   ELSE
1214!
1215!--                   Just store the respective section on the local array
1216!--                   (but only if it is available on this PE!)
1217                      IF ( section(is,s) >= nxl  .AND.  section(is,s) <= nxr ) &
1218                      THEN
1219                         local_2d = local_pf(section(is,s),:,nzb:nzt+1)
1220                      ENDIF
1221
1222                   ENDIF
1223
1224#if defined( __parallel )
1225                   IF ( netcdf_output  .AND.  netcdf_data_format > 2 )  THEN
1226!
1227!--                   ATTENTION: The following lines are a workaround, because
1228!--                              independet output does not work with the
1229!--                              current NetCDF4 installation. Therefore, data
1230!--                              are transferred from PEs having the cross
1231!--                              sections to other PEs along y having no cross
1232!--                              section data. Some of these data are the
1233!--                              output.
1234!--                   BEGIN WORKAROUND---------------------------------------
1235                      IF ( npex /= 1  .AND.  section(is,s) /= -1)  THEN
1236                         ALLOCATE( local_2d_l(nys-1:nyn+1,nzb:nzt+1) )
1237                         local_2d_l = 0.0
1238                         IF ( section(is,s) >= nxl .AND. section(is,s) <= nxr )&
1239                         THEN
1240                            local_2d_l = local_2d
1241                         ENDIF
1242#if defined( __parallel )
1243!
1244!--                      Distribute data over all PEs along x
1245                         ngp = ( nyn-nys+3 ) * ( nzt-nzb+2 )
1246                         IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
1247                         CALL MPI_ALLREDUCE( local_2d_l(nys-1,nzb),            &
1248                                             local_2d(nys-1,nzb), ngp,         &
1249                                             MPI_REAL, MPI_SUM, comm1dx, ierr )
1250#else
1251                         local_2d = local_2d_l
1252#endif
1253                         DEALLOCATE( local_2d_l )
1254                      ENDIF
1255!--                   END WORKAROUND-----------------------------------------
1256
1257!
1258!--                   Output in NetCDF4/HDF5 format.
1259!--                   Output only on those PEs where the respective cross
1260!--                   sections reside. Cross sections averaged along x are
1261!--                   output on the respective first PE along x (myidx=0).
1262                      IF ( ( section(is,s) >= nxl  .AND.  &
1263                             section(is,s) <= nxr )  .OR.  &
1264                           ( section(is,s) == -1  .AND.  myidx == 0 ) )  THEN
1265!
1266!--                      Do not output redundant ghost point data except for the
1267!--                      boundaries of the total domain.
1268#if defined( __netcdf )
1269                         IF ( nyn == ny )  THEN
1270                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1271                                                id_var_do2d(av,if),            &
1272                                                local_2d(nys:nyn+1,nzb:nzt+1), &
1273                                                start = (/ is, nys+1, 1,       &
1274                                                    do2d_yz_time_count(av) /), &
1275                                                count = (/ 1, nyn-nys+2,       &
1276                                                           nzt+2, 1 /) )
1277                         ELSE
1278                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1279                                                id_var_do2d(av,if),            &
1280                                                local_2d(nys:nyn,nzb:nzt+1),   &
1281                                                start = (/ is, nys+1, 1,       &
1282                                                    do2d_yz_time_count(av) /), &
1283                                                count = (/ 1, nyn-nys+1,       &
1284                                                           nzt+2, 1 /) )
1285                         ENDIF
1286
1287                         CALL handle_netcdf_error( 'data_output_2d', 60 )
1288
1289                      ELSE
1290!
1291!--                      Output on other PEs. Only one point is output!!
1292!--                      ATTENTION: This is a workaround (see above)!!
1293                         IF ( npex /= 1 )  THEN
1294                            nc_stat = NF90_PUT_VAR( id_set_yz(av),             &
1295                                                    id_var_do2d(av,if),        &
1296                                                    local_2d(nys:nys,nzb:nzb), &
1297                                                    start = (/ is, nys+1, 1,   &
1298                                                    do2d_yz_time_count(av) /), &
1299                                                    count = (/ 1, 1, 1, 1 /) )
1300                            CALL handle_netcdf_error( 'data_output_2d', 452 )
1301                         ENDIF
1302#endif
1303                      ENDIF
1304
1305                   ELSE
1306
1307                      IF ( data_output_2d_on_each_pe )  THEN
1308!
1309!--                      Output of partial arrays on each PE. If the cross
1310!--                      section does not reside on the PE, output special
1311!--                      index values.
1312#if defined( __netcdf )
1313                         IF ( netcdf_output  .AND.  myid == 0 )  THEN
1314                            WRITE ( 23 )  simulated_time, &
1315                                          do2d_yz_time_count(av), av
1316                         ENDIF
1317#endif
1318                         IF ( ( section(is,s) >= nxl  .AND.                  &
1319                                section(is,s) <= nxr )  .OR.                 &
1320                              ( section(is,s) == -1  .AND.  nxl-1 == -1 ) )  &
1321                         THEN
1322                            WRITE (23)  nys-1, nyn+1, nzb, nzt+1
1323                            WRITE (23)  local_2d
1324                         ELSE
1325                            WRITE (23)  -1, -1, -1, -1
1326                         ENDIF
1327
1328                      ELSE
1329!
1330!--                      PE0 receives partial arrays from all processors of the
1331!--                      respective cross section and outputs them. Here a
1332!--                      barrier has to be set, because otherwise
1333!--                      "-MPI- FATAL: Remote protocol queue full" may occur.
1334                         CALL MPI_BARRIER( comm2d, ierr )
1335
1336                         ngp = ( nyn-nys+3 ) * ( nzt-nzb+2 )
1337                         IF ( myid == 0 )  THEN
1338!
1339!--                         Local array can be relocated directly.
1340                            IF ( ( section(is,s) >= nxl  .AND.                 &
1341                                   section(is,s) <= nxr )   .OR.               &
1342                                 ( section(is,s) == -1  .AND.  nxl-1 == -1 ) ) &
1343                            THEN
1344                               total_2d(nys-1:nyn+1,nzb:nzt+1) = local_2d
1345                            ENDIF
1346!
1347!--                         Receive data from all other PEs.
1348                            DO  n = 1, numprocs-1
1349!
1350!--                            Receive index limits first, then array.
1351!--                            Index limits are received in arbitrary order from
1352!--                            the PEs.
1353                               CALL MPI_RECV( ind(1), 4, MPI_INTEGER,     &
1354                                              MPI_ANY_SOURCE, 0, comm2d,  &
1355                                              status, ierr )
1356!
1357!--                            Not all PEs have data for YZ-cross-section.
1358                               IF ( ind(1) /= -9999 )  THEN
1359                                  sender = status(MPI_SOURCE)
1360                                  DEALLOCATE( local_2d )
1361                                  ALLOCATE( local_2d(ind(1):ind(2), &
1362                                                     ind(3):ind(4)) )
1363                                  CALL MPI_RECV( local_2d(ind(1),ind(3)), ngp, &
1364                                                 MPI_REAL, sender, 1, comm2d,  &
1365                                                 status, ierr )
1366                                  total_2d(ind(1):ind(2),ind(3):ind(4)) = &
1367                                                                        local_2d
1368                               ENDIF
1369                            ENDDO
1370!
1371!--                         Output of the total cross-section.
1372                            IF ( iso2d_output )  THEN
1373                               WRITE (23)  total_2d(0:ny+1,nzb:nzt+1)
1374                            ENDIF
1375!
1376!--                         Relocate the local array for the next loop increment
1377                            DEALLOCATE( local_2d )
1378                            ALLOCATE( local_2d(nys-1:nyn+1,nzb:nzt+1) )
1379
1380#if defined( __netcdf )
1381                            IF ( netcdf_output )  THEN
1382                               nc_stat = NF90_PUT_VAR( id_set_yz(av),          &
1383                                                    id_var_do2d(av,if),        &
1384                                                    total_2d(0:ny+1,nzb:nzt+1),&
1385                               start = (/ is, 1, 1, do2d_yz_time_count(av) /), &
1386                                                count = (/ 1, ny+2, nz+2, 1 /) )
1387                               CALL handle_netcdf_error( 'data_output_2d', 61 )
1388                            ENDIF
1389#endif
1390
1391                         ELSE
1392!
1393!--                         If the cross section resides on the PE, send the
1394!--                         local index limits, otherwise send -9999 to PE0.
1395                            IF ( ( section(is,s) >= nxl  .AND.                 &
1396                                   section(is,s) <= nxr )  .OR.                &
1397                                 ( section(is,s) == -1  .AND.  nxl-1 == -1 ) ) &
1398                            THEN
1399                               ind(1) = nys-1; ind(2) = nyn+1
1400                               ind(3) = nzb;   ind(4) = nzt+1
1401                            ELSE
1402                               ind(1) = -9999; ind(2) = -9999
1403                               ind(3) = -9999; ind(4) = -9999
1404                            ENDIF
1405                            CALL MPI_SEND( ind(1), 4, MPI_INTEGER, 0, 0, &
1406                                           comm2d, ierr )
1407!
1408!--                         If applicable, send data to PE0.
1409                            IF ( ind(1) /= -9999 )  THEN
1410                               CALL MPI_SEND( local_2d(nys-1,nzb), ngp, &
1411                                              MPI_REAL, 0, 1, comm2d, ierr )
1412                            ENDIF
1413                         ENDIF
1414!
1415!--                      A barrier has to be set, because otherwise some PEs may
1416!--                      proceed too fast so that PE0 may receive wrong data on
1417!--                      tag 0
1418                         CALL MPI_BARRIER( comm2d, ierr )
1419                      ENDIF
1420
1421                   ENDIF
1422#else
1423                   IF ( iso2d_output )  THEN
1424                      WRITE (23)  local_2d(nys:nyn+1,nzb:nzt+1)
1425                   ENDIF
1426#if defined( __netcdf )
1427                   IF ( netcdf_output )  THEN
1428                      nc_stat = NF90_PUT_VAR( id_set_yz(av),                   &
1429                                              id_var_do2d(av,if),              &
1430                                              local_2d(nys:nyn+1,nzb:nzt+1),   &
1431                               start = (/ is, 1, 1, do2d_xz_time_count(av) /), &
1432                                              count = (/ 1, ny+2, nz+2, 1 /) )
1433                      CALL handle_netcdf_error( 'data_output_2d', 452 )
1434                   ENDIF
1435#endif
1436#endif
1437                   do2d_yz_n = do2d_yz_n + 1
1438!
1439!--                Write LOCAL-parameter set for ISO2D.
1440                   IF ( myid == 0  .AND.  iso2d_output )  THEN
1441                      IF ( section(is,s) /= -1 )  THEN
1442                         WRITE ( section_chr, '(''x = '',F8.2,'' m  (GP '',I3, &
1443                                               &'')'')'                        &
1444                               )  section(is,s)*dx, section(is,s)
1445                      ELSE
1446                         section_chr = 'averaged along x'
1447                      ENDIF
1448                      IF ( av == 0 )  THEN
1449                         rtext = TRIM( do2d(av,if) ) // '  t = ' //    &
1450                                 TRIM( simulated_time_chr ) // '  ' // &
1451                                 TRIM( section_chr )
1452                      ELSE
1453                         rtext = TRIM( do2d(av,if) ) // '  averaged t = ' // &
1454                                 TRIM( simulated_time_chr ) // '  ' //       &
1455                                 TRIM( section_chr )
1456                      ENDIF
1457                      WRITE (29,LOCAL)
1458                   ENDIF
1459
1460             END SELECT
1461
1462             is = is + 1
1463          ENDDO loop1
1464
1465       ENDIF
1466
1467       if = if + 1
1468       l = MAX( 2, LEN_TRIM( do2d(av,if) ) )
1469       do2d_mode = do2d(av,if)(l-1:l)
1470
1471    ENDDO
1472
1473!
1474!-- Deallocate temporary arrays.
1475    IF ( ALLOCATED( level_z ) )  DEALLOCATE( level_z )
1476    DEALLOCATE( local_pf, local_2d )
1477#if defined( __parallel )
1478    IF ( .NOT.  data_output_2d_on_each_pe  .AND.  myid == 0 )  THEN
1479       DEALLOCATE( total_2d )
1480    ENDIF
1481#endif
1482
1483!
1484!-- Close plot output file.
1485    file_id = 20 + s
1486
1487    IF ( data_output_2d_on_each_pe )  THEN
1488       CALL close_file( file_id )
1489    ELSE
1490       IF ( myid == 0 )  CALL close_file( file_id )
1491    ENDIF
1492
1493
1494    CALL cpu_log (log_point(3),'data_output_2d','stop','nobarrier')
1495
1496 END SUBROUTINE data_output_2d
Note: See TracBrowser for help on using the repository browser.