source: palm/trunk/SOURCE/restart_data_mpi_io_mod.f90 @ 4526

Last change on this file since 4526 was 4500, checked in by suehring, 5 years ago

Surface output: correct output of ground/wall-heat flux at USM surfaces; add conversion factor to heat- and momentum-flux outputs; data_output_2d: Unify output conversion of sensible and latent heat flux; data-output module: avoid uninitialized variables; restart_data_mpi_io: fix overlong lines

  • Property svn:keywords set to Id
File size: 81.8 KB
RevLine 
[4495]1!> @file restart_data_mpi_io_mod.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
6! Public License as published by the Free Software Foundation, either version 3 of the License, or
7! (at your option) any later version.
8!
9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
10! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
11! Public License for more details.
12!
13! You should have received a copy of the GNU General Public License along with PALM. If not, see
14! <http://www.gnu.org/licenses/>.
15!
16! Copyright 1997-2020 Leibniz Universitaet Hannover
17! -------------------------------------------------------------------------------------------------!
18!
19! Current revisions:
20! -----------------
21!
22!
23! Former revisions:
24! -----------------
25! $Id: restart_data_mpi_io_mod.f90 4500 2020-04-17 10:12:45Z oliver.maas $
[4500]26! Fix too long lines
27!
28! 4498 2020-04-15 14:26:31Z raasch
[4498]29! bugfix for creation of filetypes, argument removed from rd_mpi_io_open
30!
31! 4497 2020-04-15 10:20:51Z raasch
[4497]32! last bugfix deactivated because of compile problems
33!
34! 4496 2020-04-15 08:37:26Z raasch
[4496]35! problem with posix read arguments for surface data fixed
36!
37! 4495 2020-04-13 20:11:20Z raasch
38! Initial version (K. Ketelsen), adjusted to PALM formatting standards (s. Raasch)
[4495]39!
40!
41!
42! Description:
43! ------------
44!> Routines for restart data handling using MPI-IO.
45!--------------------------------------------------------------------------------------------------!
46 MODULE restart_data_mpi_io_mod
47
48#if defined( __parallel )
49#if defined( __mpifh )
50    INCLUDE "mpif.h"
51#else
52    USE MPI
53#endif
54#else
55    USE posix_interface,                                                                           &
56        ONLY:  posix_close, posix_lseek, posix_open, posix_read, posix_write
57#endif
58
59    USE, INTRINSIC ::  ISO_C_BINDING
60
61    USE control_parameters,                                                                        &
62        ONLY:  include_total_domain_boundaries
63
64    USE exchange_horiz_mod,                                                                        &
65        ONLY:  exchange_horiz, exchange_horiz_2d
66
67    USE indices,                                                                                   &
68        ONLY:  nbgp, nnx, nny, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg, nz, nzb, nzt
69
70    USE kinds
71
72    USE pegrid,                                                                                    &
73        ONLY:  comm1dx, comm1dy, comm2d, myid, myidx, myidy, npex, npey, numprocs, pdims
74
75
76    IMPLICIT NONE
77
78    LOGICAL ::  all_pes_write                     !< all PEs have data to write
79    LOGICAL ::  filetypes_created                 !<
80    LOGICAL ::  print_header_now = .TRUE.         !<
81    LOGICAL ::  rd_flag                           !< file is opened for read
82    LOGICAL ::  wr_flag                           !< file is opened for write
83
84#if defined( __parallel )
85    INTEGER(iwp)                     :: ierr                              !< error status of MPI-calls
86    INTEGER(iwp), PARAMETER          :: rd_offset_kind = MPI_OFFSET_KIND  !< Adress or Offset kind
87    INTEGER(iwp), PARAMETER          :: rd_status_size = MPI_STATUS_SIZE  !<
88#else
89    INTEGER(iwp), PARAMETER          :: rd_offset_kind = C_SIZE_T         !<
90    INTEGER(iwp), PARAMETER          :: rd_status_size = 1       !< Not required in sequential mode
91#endif
92
93    INTEGER(iwp)                     :: debug_level = 1 !< TODO: replace with standard debug output steering
94
95    INTEGER(iwp)                     :: fh            !< MPI-IO file handle
96    INTEGER(iwp)                     :: ft_surf = -1  !< MPI filetype surface data
97#if defined( __parallel )
98    INTEGER(iwp)                     :: ft_2di_nb     !< MPI filetype 2D array INTEGER no outer boundary
99    INTEGER(iwp)                     :: ft_2d         !< MPI filetype 2D array REAL with outer boundaries
100    INTEGER(iwp)                     :: ft_3d         !< MPI filetype 3D array REAL with outer boundaries
101    INTEGER(iwp)                     :: ft_3dsoil     !< MPI filetype for 3d-soil array
102#endif
103    INTEGER(iwp)                     :: glo_start     !< global start index on this PE
104    INTEGER(iwp)                     :: nr_val        !< local number of values in x and y direction
105    INTEGER(iwp)                     :: total_number_of_surface_values    !< total number of values for one variable
106
107    INTEGER(KIND=rd_offset_kind) ::  array_position   !<
108    INTEGER(KIND=rd_offset_kind) ::  header_position  !<
109
110    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_end_index     !<
111    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_start_index   !<
112    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_global_start  !<
113
114    REAL(KIND=wp) ::  mb_processed  !<
115
116!
117!-- Handling of outer boundaries
118    TYPE local_boundaries
119       INTEGER(iwp) ::  nnx
120       INTEGER(iwp) ::  nny
121       INTEGER(iwp) ::  nx
122       INTEGER(iwp) ::  nxl
123       INTEGER(iwp) ::  nxr
124       INTEGER(iwp) ::  ny
125       INTEGER(iwp) ::  nyn
126       INTEGER(iwp) ::  nys
127    END TYPE local_boundaries
128
129    TYPE(local_boundaries) ::  lb  !<
130
131!
132!-- General Header (first 32 byte in restart file)
133    TYPE general_header
134       INTEGER(iwp) :: nr_int         !< number of INTEGER entries in header
135       INTEGER(iwp) :: nr_char        !< number of Text strings entries in header
136       INTEGER(iwp) :: nr_real        !< number of REAL entries in header
137       INTEGER(iwp) :: nr_arrays      !< number of arrays in restart files
138       INTEGER(iwp) :: total_nx       !< total number of points in x-direction
139       INTEGER(iwp) :: total_ny       !< total number of points in y-direction
140       INTEGER(iwp) :: i_outer_bound  !< if 1, outer boundaries are stored in restart file
141       INTEGER(iwp) :: endian         !< little endian (1) or big endian (2) internal format
142    END TYPE general_header
143
144    TYPE(general_header), TARGET ::  tgh
145
146!
147!-- Declaration of varibales for file header section
148    INTEGER(KIND=rd_offset_kind)                ::  header_int_index
149    INTEGER, PARAMETER                          ::  max_nr_int=256
150    CHARACTER(LEN=32), DIMENSION(max_nr_int)    ::  int_names
151    INTEGER(KIND=iwp), DIMENSION(max_nr_int)    ::  int_values
152
153    INTEGER(KIND=rd_offset_kind)                ::  header_char_index
154    INTEGER, PARAMETER                          ::  max_nr_text=128
155    CHARACTER(LEN=128), DIMENSION(max_nr_text)  ::  text_lines
156
157    INTEGER(KIND=rd_offset_kind)                ::  header_real_index
158    INTEGER, PARAMETER                          ::  max_nr_real=256
159    CHARACTER(LEN=32), DIMENSION(max_nr_real)   ::  real_names
160    REAL(KIND=wp), DIMENSION(max_nr_real)       ::  real_values
161
162    INTEGER(KIND=rd_offset_kind)                ::  header_arr_index
163    INTEGER, PARAMETER                          ::  max_nr_arrays=600
164    CHARACTER(LEN=32), DIMENSION(max_nr_arrays) ::  array_names
165    INTEGER(KIND=rd_offset_kind), DIMENSION(max_nr_arrays) :: array_offset
166
167    SAVE
168
169    PRIVATE
170
171    PUBLIC  mb_processed, total_number_of_surface_values
172
173!
174!-- PALM interfaces
175    INTERFACE rd_mpi_io_check_array
176       MODULE PROCEDURE rd_mpi_io_check_array
177    END INTERFACE rd_mpi_io_check_array
178
179    INTERFACE rd_mpi_io_close
180       MODULE PROCEDURE rd_mpi_io_close
181    END INTERFACE rd_mpi_io_close
182
183    INTERFACE rd_mpi_io_open
184       MODULE PROCEDURE rd_mpi_io_open
185    END INTERFACE rd_mpi_io_open
186
187    INTERFACE rrd_mpi_io
188       MODULE PROCEDURE rrd_mpi_io_char
189       MODULE PROCEDURE rrd_mpi_io_int
190       MODULE PROCEDURE rrd_mpi_io_int_2d
191       MODULE PROCEDURE rrd_mpi_io_logical
192       MODULE PROCEDURE rrd_mpi_io_real
193       MODULE PROCEDURE rrd_mpi_io_real_2d
194       MODULE PROCEDURE rrd_mpi_io_real_3d
195       MODULE PROCEDURE rrd_mpi_io_real_3d_soil
196    END INTERFACE rrd_mpi_io
197
198    INTERFACE rrd_mpi_io_global_array
199       MODULE PROCEDURE rrd_mpi_io_global_array_int_1d
200       MODULE PROCEDURE rrd_mpi_io_global_array_real_1d
201       MODULE PROCEDURE rrd_mpi_io_global_array_real_2d
202       MODULE PROCEDURE rrd_mpi_io_global_array_real_3d
203       MODULE PROCEDURE rrd_mpi_io_global_array_real_4d
204    END INTERFACE rrd_mpi_io_global_array
205
206    INTERFACE rrd_mpi_io_surface
207       MODULE PROCEDURE rrd_mpi_io_surface
208       MODULE PROCEDURE rrd_mpi_io_surface_2d
209    END INTERFACE rrd_mpi_io_surface
210
211    INTERFACE rd_mpi_io_surface_filetypes
212       MODULE PROCEDURE rd_mpi_io_surface_filetypes
213    END INTERFACE rd_mpi_io_surface_filetypes
214
215    INTERFACE wrd_mpi_io
216       MODULE PROCEDURE wrd_mpi_io_char
217       MODULE PROCEDURE wrd_mpi_io_int
218       MODULE PROCEDURE wrd_mpi_io_int_2d
219       MODULE PROCEDURE wrd_mpi_io_logical
220       MODULE PROCEDURE wrd_mpi_io_real
221       MODULE PROCEDURE wrd_mpi_io_real_2d
222       MODULE PROCEDURE wrd_mpi_io_real_3d
223       MODULE PROCEDURE wrd_mpi_io_real_3d_soil
224    END INTERFACE wrd_mpi_io
225
226    INTERFACE wrd_mpi_io_global_array
227       MODULE PROCEDURE wrd_mpi_io_global_array_int_1d
228       MODULE PROCEDURE wrd_mpi_io_global_array_real_1d
229       MODULE PROCEDURE wrd_mpi_io_global_array_real_2d
230       MODULE PROCEDURE wrd_mpi_io_global_array_real_3d
231       MODULE PROCEDURE wrd_mpi_io_global_array_real_4d
232    END INTERFACE wrd_mpi_io_global_array
233
234    INTERFACE wrd_mpi_io_surface
235       MODULE PROCEDURE wrd_mpi_io_surface
236       MODULE PROCEDURE wrd_mpi_io_surface_2d
237    END INTERFACE wrd_mpi_io_surface
238
239    PUBLIC  rd_mpi_io_check_array, rd_mpi_io_close, rd_mpi_io_open, rrd_mpi_io,                    &
240            rrd_mpi_io_global_array, rrd_mpi_io_surface, rd_mpi_io_surface_filetypes, wrd_mpi_io,  &
241            wrd_mpi_io_global_array, wrd_mpi_io_surface
242
243
244 CONTAINS
245
246
247!--------------------------------------------------------------------------------------------------!
248! Description:
249! ------------
250!> Open restart file for read or write with MPI-IO
251!--------------------------------------------------------------------------------------------------!
[4498]252 SUBROUTINE rd_mpi_io_open( action, file_name )
[4495]253
254    IMPLICIT NONE
255
256    CHARACTER(LEN=*), INTENT(IN)  ::  action                           !<
257    CHARACTER(LEN=*), INTENT(IN)  ::  file_name                        !<
258
259    INTEGER(iwp)                  ::  i                                !<
260    INTEGER(iwp)                  ::  gh_size                          !<
261
262    INTEGER(KIND=rd_offset_kind)  ::  offset                           !<
263
264#if defined( __parallel )
265    INTEGER, DIMENSION(rd_status_size) ::  status                      !<
266#endif
267
268#if ! defined( __parallel )
269    TYPE(C_PTR)                   ::  buf_ptr                          !<
270#endif
271
272
273    offset = 0
274
275    rd_flag = ( TRIM( action ) == 'READ'  .OR. TRIM( action ) == 'read'  )
276    wr_flag = ( TRIM( action ) == 'WRITE' .OR. TRIM( action ) == 'write' )
277
278!
279!-- Create subarrays and file types
280    filetypes_created = .FALSE.
281
282!
283!-- In case of read it is not known yet if data include total domain. Filetypes will be created
284!-- further below.
[4498]285    IF ( wr_flag)  THEN
[4495]286       CALL rs_mpi_io_create_filetypes
287       filetypes_created = .TRUE.
288    ENDIF
289
290!
291!-- Open file for MPI-IO
292#if defined( __parallel )
293    IF ( rd_flag )  THEN
294       CALL MPI_FILE_OPEN( comm2d, TRIM( file_name ), MPI_MODE_RDONLY, MPI_INFO_NULL, fh, ierr )
295       WRITE (9,*) 'Open MPI-IO restart file for read  ==> ', TRIM( file_name )
296    ELSEIF ( wr_flag )  THEN
297       CALL MPI_FILE_OPEN( comm2d, TRIM( file_name ), MPI_MODE_CREATE+MPI_MODE_WRONLY,             &
298                           MPI_INFO_NULL, fh, ierr )
299       WRITE (9,*) 'Open MPI-IO restart file for write ==> ', TRIM( file_name )
300    ELSE
301       CALL rs_mpi_io_error( 1 )
302    ENDIF
303#else
304    IF ( rd_flag )  THEN
305       fh = posix_open( TRIM( file_name ), .TRUE. )
306       WRITE (9,*) 'Open sequential restart file for read  ==> ', TRIM( file_name ), ' ', fh
307    ELSEIF ( wr_flag )  THEN
308       fh = posix_open( TRIM( file_name ), .FALSE. )
309       WRITE (9,*) 'Open sequential restart file for write ==> ', TRIM(file_name), ' ', fh
310    ELSE
311       CALL rs_mpi_io_error( 1 )
312    ENDIF
313
314    IF ( fh < 0 )  CALL rs_mpi_io_error( 6 )
315#endif
316
317    array_position  = 65536          !> Start offset for writing 2-D and 3.D arrays at 64 k
318    header_position = 0
319
320    header_int_index   = 1
321    header_char_index  = 1
322    header_real_index   = 1
323    header_arr_index   = 1
324
325    int_names    = ' '
326    int_values   = 0
327    text_lines   = ' '
328    real_names   = ' '
329    real_values  = 0.0
330    array_names  = ' '
331    array_offset = 0
332
333    int_names(1)     = 'nx'
334    int_values(1)    = nx
335    int_names(2)     = 'ny'
336    int_values(2)    = ny
337    int_names(3)     = 'nz'
338    int_values(3)    = nz
339    header_int_index = header_int_index+3
340
341    DO   i = 1, max_nr_arrays
342       array_offset(i) = 0
343       array_names(i)  = ' '
344    ENDDO
345
346    gh_size = STORAGE_SIZE( tgh ) / 8
347
348    IF ( rd_flag )  THEN
349!
350!--    File is open for read.
351#if defined( __parallel )
352!--    Set the default view
353       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
354!
355!--    Read the file header size
356       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
357       CALL MPI_FILE_READ( fh, tgh, gh_size, MPI_BYTE, status, ierr )
358#else
359       CALL posix_lseek( fh, header_position )
360       buf_ptr = C_LOC( tgh )
361       CALL posix_read( fh, buf_ptr, gh_size )
362#endif
363       header_position = header_position + gh_size
364
365       include_total_domain_boundaries = ( tgh%i_outer_bound == 1 )
366
367!
[4498]368!--    File types depend on if boundaries of the total domain is included in data. This has been
369!--    checked with the previous statement.
370       CALL rs_mpi_io_create_filetypes
371       filetypes_created = .TRUE.
[4495]372
373#if defined( __parallel )
374!
375!--    Read INTEGER values
376       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
377       CALL MPI_FILE_READ( fh, int_names, SIZE( int_names ) * 32, MPI_CHAR, status, ierr )
378       header_position = header_position + SIZE( int_names ) * 32
379
380       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
381       CALL MPI_FILE_READ (fh, int_values, SIZE( int_values ), MPI_INT, status, ierr )
382       header_position = header_position + SIZE( int_values ) * iwp
383!
384!--    Character entries
385       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
386       CALL MPI_FILE_READ( fh, text_lines, SIZE( text_lines ) * 128, MPI_CHAR, status, ierr )
387       header_position = header_position+size(text_lines) * 128
388!
389!--    REAL values
390       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
391       CALL MPI_FILE_READ( fh, real_names, SIZE( real_names ) * 32, MPI_CHAR, status, ierr )
392       header_position = header_position + SIZE( real_names ) * 32
393
394       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
395       CALL MPI_FILE_READ( fh, real_values, SIZE( real_values ), MPI_REAL, status, ierr )
396       header_position = header_position + SIZE( real_values ) * wp
397!
398!--    2d- and 3d-array headers
399       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
400       CALL MPI_FILE_READ( fh, array_names, SIZE( array_names ) * 32, MPI_CHAR, status, ierr )
401       header_position = header_position + SIZE( array_names ) * 32
402
403       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
404       CALL MPI_FILE_READ( fh, array_offset, SIZE( array_offset ) * MPI_OFFSET_KIND, MPI_BYTE,     &
405                           status,ierr )   ! there is no I*8 datatype in Fortran
406       header_position = header_position + SIZE( array_offset ) * rd_offset_kind
407#else
408       CALL posix_lseek( fh, header_position )
409       CALL posix_read( fh, int_names )
410       header_position = header_position + SIZE( int_names ) * 32
411
412       CALL posix_lseek( fh, header_position )
413       CALL posix_read( fh, int_values, SIZE( int_values ) )
414       header_position = header_position + SIZE( int_values ) * iwp
415!
416!--    Character entries
417       CALL posix_lseek( fh, header_position )
418       CALL posix_read( fh, text_lines )
419       header_position = header_position + SIZE( text_lines ) * 128
420!
421!--    REAL values
422       CALL posix_lseek( fh, header_position )
423       CALL posix_read( fh, real_names )
424       header_position = header_position + SIZE( real_names ) * 32
425
426       CALL posix_lseek( fh, header_position )
427       CALL posix_read( fh, real_values, SIZE( real_values ) )
428       header_position = header_position + SIZE( real_values ) * wp
429!
430!--    2d- and 3d-array headers
431       CALL posix_lseek( fh, header_position )
432       CALL posix_read( fh, array_names )
433       header_position = header_position + SIZE( array_names ) * 32
434
435       CALL posix_lseek( fh, header_position )
436       CALL posix_read( fh, array_offset, SIZE( array_offset ) ) ! there is no I*8 datatype in Fortran
437       header_position = header_position + SIZE( array_offset ) * rd_offset_kind
438#endif
439       IF ( debug_level >= 2 )  THEN
440          WRITE (9,*) 'header positio after array metadata  ', header_position
441       ENDIF
442
443       IF ( print_header_now )  CALL rs_mpi_io_print_header
444
445    ENDIF
446
447 END SUBROUTINE rd_mpi_io_open
448
449
450!--------------------------------------------------------------------------------------------------!
451! Description:
452! ------------
453!> Check, if array exists in restart file
454!--------------------------------------------------------------------------------------------------!
455 SUBROUTINE rd_mpi_io_check_array( name, found )
456
457    IMPLICIT NONE
458
459    CHARACTER(LEN=*), INTENT(IN) ::  name  !<
460
461    INTEGER(iwp) ::  i  !<
462
463    LOGICAl      ::  found  !<
464
465
466    DO  i = 1, tgh%nr_arrays
467       IF ( TRIM( array_names(i) ) == TRIM( name ) )  THEN
468          array_position = array_offset(i)
469          found = .TRUE.
470          RETURN
471       ENDIF
472    ENDDO
473
474    found = .FALSE.
475
476 END SUBROUTINE rd_mpi_io_check_array
477
478
479
480!--------------------------------------------------------------------------------------------------!
481! Description:
482! ------------
483!> Read INTEGER with MPI-IO
484!--------------------------------------------------------------------------------------------------!
485 SUBROUTINE rrd_mpi_io_int( name, value, found )
486
487    IMPLICIT NONE
488
489    CHARACTER(LEN=*), INTENT(IN)   :: name
490
491    INTEGER(iwp)                   ::  i
492    INTEGER(KIND=iwp), INTENT(OUT) ::  value
493
494    LOGICAL                        ::  lo_found
495    LOGICAL, INTENT(OUT), OPTIONAL ::  found
496
497
498    lo_found = .FALSE.
499    value = 0
500
501    DO  i = 1, tgh%nr_int
502       IF ( TRIM(int_names(i)) == TRIM( name ) )  THEN
503          IF ( debug_level >= 2 )  WRITE(9,*) 'INTEGER variable found ', name
504          value = int_values(i)
505          lo_found = .TRUE.
506          EXIT
507       ENDIF
508    ENDDO
509
510    IF ( PRESENT( found ) )  THEN
511       found = lo_found
512       RETURN
513    ENDIF
514
515    IF ( .NOT. lo_found )  THEN
516       WRITE(9,*)  'INTEGER not found ', name
517       CALL rs_mpi_io_error( 3 )
518    ENDIF
519
520 END SUBROUTINE rrd_mpi_io_int
521
522
523
524!--------------------------------------------------------------------------------------------------!
525! Description:
526! ------------
527!> Read REAL with MPI-IO
528!--------------------------------------------------------------------------------------------------!
529 SUBROUTINE rrd_mpi_io_real( name, value, found )
530
531    IMPLICIT NONE
532
533    CHARACTER(LEN=*), INTENT(IN)   ::  name
534
535    INTEGER(iwp)                   ::  i
536
537    LOGICAL                        ::  lo_found
538    LOGICAL, INTENT(OUT), OPTIONAL ::  found
539
540    REAL(KIND=wp), INTENT(OUT)     ::  value
541
542
543    lo_found = .FALSE.
544    value = 0.0
545
546    DO  i = 1, tgh%nr_real
547       IF ( TRIM(real_names(i)) == TRIM( name ) )  THEN
548          IF ( debug_level >= 2 )  WRITE(9,*) 'REAL variable found ', name
549          value = real_values(i)
550          lo_found = .TRUE.
551          EXIT
552       ENDIF
553    ENDDO
554
555    IF ( PRESENT( found ) )  THEN
556       found = lo_found
557       RETURN
558    ENDIF
559
560    IF ( .NOT. lo_found )  THEN
561       WRITE(9,*) 'REAL value not found ', name
562       CALL rs_mpi_io_error(3)
563    ENDIF
564
565 END SUBROUTINE rrd_mpi_io_real
566
567
568
569!--------------------------------------------------------------------------------------------------!
570! Description:
571! ------------
572!> Read 2d-real array with MPI-IO
573!--------------------------------------------------------------------------------------------------!
574 SUBROUTINE rrd_mpi_io_real_2d( name, data )
575
576    IMPLICIT NONE
577
578    CHARACTER(LEN=*), INTENT(IN)       ::  name
579
580#if defined( __parallel )
581    INTEGER, DIMENSION(rd_status_size) ::  status
582#endif
583    INTEGER(iwp)                       ::  i
584
585    LOGICAL                            ::  found
586
587    REAL(wp), INTENT(INOUT), DIMENSION(nysg:nyng,nxlg:nxrg) ::  data
588
589    REAL(KIND=wp), DIMENSION(lb%nxl:lb%nxr,lb%nys:lb%nyn)   ::  array_2d
590
591
592    found = .FALSE.
593
594    DO  i = 1, tgh%nr_arrays
595       IF ( TRIM(array_names(i)) == TRIM( name ) )  THEN
596          array_position = array_offset(i)
597          found = .TRUE.
598          EXIT
599       ENDIF
600    ENDDO
601
602     IF ( found )  THEN
603#if defined( __parallel )
604        CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL, ierr )
605        CALL MPI_FILE_READ_ALL( fh, array_2d, SIZE( array_2d), MPI_REAL, status, ierr )
606#else
607        CALL posix_lseek( fh, array_position )
608        CALL posix_read( fh, array_2d, SIZE( array_2d ) )
609#endif
610
611        IF ( include_total_domain_boundaries)  THEN
612           DO  i = lb%nxl, lb%nxr
613              data(lb%nys-nbgp:lb%nyn-nbgp,i-nbgp) = array_2d(i,lb%nys:lb%nyn)
614           ENDDO
615           IF ( debug_level >= 2)  WRITE(9,*) 'r2f_ob ', TRIM(name),' ', SUM( data(nys:nyn,nxl:nxr) )
616        ELSE
617           DO  i = nxl, nxr
618              data(nys:nyn,i) = array_2d(i,nys:nyn)
619           ENDDO
620           IF ( debug_level >= 2) WRITE(9,*) 'r2f ', TRIM( name ),' ', SUM( data)
621        ENDIF
622
623        CALL exchange_horiz_2d( data )
624
625     ELSE
626        WRITE(9,*) 'array_2D not found ', name
627        CALL rs_mpi_io_error( 2 )
628     ENDIF
629
630 END SUBROUTINE rrd_mpi_io_real_2d
631
632
633
634!--------------------------------------------------------------------------------------------------!
635! Description:
636! ------------
637!> Read 2d-INTEGER array with MPI-IO
638!--------------------------------------------------------------------------------------------------!
639 SUBROUTINE rrd_mpi_io_int_2d( name, data )
640
641    IMPLICIT NONE
642
643    CHARACTER(LEN=*), INTENT(IN)        ::  name
644
645    INTEGER(iwp)                        ::  i
646    INTEGER(iwp)                        ::  j
647
648#if defined( __parallel )
649    INTEGER, DIMENSION(rd_status_size)  ::  status
650#endif
651
652    INTEGER, DIMENSION(nxl:nxr,nys:nyn) ::  array_2d
653
654    INTEGER(KIND=iwp), INTENT(INOUT), DIMENSION(:,:) ::  data
655
656    LOGICAL                             ::  found
657
658
659    found = .FALSE.
660
661    DO  i = 1, tgh%nr_arrays
662       IF ( TRIM(array_names(i)) == TRIM( name ) )  THEN
663          array_position = array_offset(i)
664          found = .TRUE.
665          EXIT
666       ENDIF
667    ENDDO
668
669    IF ( found )  THEN
670
671       IF ( ( nxr - nxl + 1 + 2*nbgp ) == SIZE( data, 2 ) )  THEN
672!
673!--       Output array with Halos.
674!--       ATTENTION: INTEGER array with ghost boundaries are not implemented yet. This kind of array
675!--                  would be dimensioned in the caller subroutine like this:
676!--                  INTEGER, DIMENSION(nysg:nyng,nxlg:nxrg)::  data
677          CALL rs_mpi_io_error( 2 )
678
679       ELSEIF ( (nxr-nxl+1) == SIZE( data, 2 ) )  THEN
680!
681!--       INTEGER input array without Halos.
682!--       This kind of array is dimensioned in the caller subroutine
683!--       INTEGER, DIMENSION(nys:nyn,nxl:nxr) ::  data
684
685#if defined( __parallel )
686          CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native',            &
687                                  MPI_INFO_NULL, ierr )
688          CALL MPI_FILE_READ_ALL( fh, array_2d, SIZE( array_2d), MPI_INTEGER, status, ierr )
689#else
690          CALL posix_lseek( fh, array_position )
691          CALL posix_read( fh, array_2d, SIZE( array_2d ) )
692#endif
693
694          DO  j = nys, nyn
695             DO  i = nxl, nxr
696                data(j-nys+1,i-nxl+1) = array_2d(i,j)
697             ENDDO
698          ENDDO
699
700          IF ( debug_level >= 2 )  WRITE(9,*) 'r2i ', TRIM( name ),' ', SUM( array_2d )
701
702       ELSE
703          WRITE (9,*) '### rrd_mpi_io_int_2d  array: ', TRIM( name )
704          CALL rs_mpi_io_error( 4 )
705       ENDIF
706
707    ELSE
708
709       WRITE(9,*) 'array_2D not found ', name
710       CALL rs_mpi_io_error( 2 )
711
712    ENDIF
713
714 END SUBROUTINE rrd_mpi_io_int_2d
715
716
717
718!--------------------------------------------------------------------------------------------------!
719! Description:
720! ------------
721!> Read 2d-REAL array with MPI-IO
722!--------------------------------------------------------------------------------------------------!
723 SUBROUTINE rrd_mpi_io_real_3d( name, data )
724
725    IMPLICIT NONE
726
727    CHARACTER(LEN=*), INTENT(IN)       ::  name
728
729    INTEGER(iwp)                       ::  i
730
731#if defined( __parallel )
732    INTEGER, DIMENSION(rd_status_size) ::  status
733#endif
734
735    LOGICAL                            ::  found
736
737    REAL(KIND=wp), DIMENSION(nzb:nzt+1,lb%nxl:lb%nxr,lb%nys:lb%nyn)   ::  array_3d
738
739    REAL(wp), INTENT(INOUT), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  data
740
741
742    found = .FALSE.
743
744    DO  i = 1, tgh%nr_arrays
745       IF ( TRIM(array_names(i)) == TRIM( name ) )  THEN
746          array_position = array_offset(i)
747          found = .TRUE.
748          EXIT
749       ENDIF
750    ENDDO
751
752    IF ( found )  THEN
753#if defined( __parallel )
754       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL, ierr )
755       CALL MPI_FILE_READ_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
756#else
757       CALL posix_lseek( fh, array_position )
758       CALL posix_read(fh, array_3d, SIZE( array_3d ) )
759#endif
760       IF ( include_total_domain_boundaries)  THEN
761          DO  i = lb%nxl, lb%nxr
762             data(:,lb%nys-nbgp:lb%nyn-nbgp,i-nbgp) = array_3d(:,i,lb%nys:lb%nyn)
763          ENDDO
764          IF ( debug_level >= 2 )  WRITE(9,*) 'r3f_ob ', TRIM( name ),' ', SUM( data(:,nys:nyn,nxl:nxr) )
765       ELSE
766          DO  i = nxl, nxr
767             data(:,nys:nyn,i) = array_3d(:,i,nys:nyn)
768          ENDDO
769          IF ( debug_level >= 2 )  WRITE(9,*) 'r3f ', TRIM( name ),' ', SUM( data )
770       ENDIF
771
772       CALL exchange_horiz( data, nbgp )
773
774    ELSE
775       WRITE(9,*)  'array_3D not found ', name
776       CALL rs_mpi_io_error(2)
777    ENDIF
778
779 END SUBROUTINE rrd_mpi_io_real_3d
780
781
782
783!--------------------------------------------------------------------------------------------------!
784! Description:
785! ------------
786!> Read 3d-REAL soil array with MPI-IO
787!> nzb_soil, nzt_soil are located in the module land_surface_model_mod. Since Fortran does not allow
788!> cross referencing of module variables, it is required to pass these variables as arguments.
789!--------------------------------------------------------------------------------------------------!
790 SUBROUTINE rrd_mpi_io_real_3d_soil( name, data, nzb_soil, nzt_soil )
791
792    IMPLICIT NONE
793
794    CHARACTER(LEN=*), INTENT(IN)       ::  name
795
796    INTEGER(iwp)                       ::  i
797    INTEGER, INTENT(IN)                ::  nzb_soil
798    INTEGER, INTENT(IN)                ::  nzt_soil
799
800#if defined( __parallel )
801    INTEGER, DIMENSION(rd_status_size) ::  status
802#endif
803
804    LOGICAL                            ::  found
805
806    REAL(KIND=wp), DIMENSION(nzb_soil:nzt_soil,lb%nxl:lb%nxr,lb%nys:lb%nyn)   ::  array_3d
807
808    REAL(wp), INTENT(INOUT), DIMENSION(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ::  data
809
810
811    found = .FALSE.
812
813    DO  i = 1, tgh%nr_arrays
814       IF ( TRIM(array_names(i)) == TRIM( name ) )  THEN
815          array_position = array_offset(i)
816          found = .TRUE.
817          EXIT
818       ENDIF
819    ENDDO
820
821    IF ( found )  THEN
822#if defined( __parallel )
823       CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
824       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3dsoil, 'native', MPI_INFO_NULL,   &
825                               ierr )
826       CALL MPI_FILE_READ_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
827       CALL MPI_TYPE_FREE( ft_3dsoil, ierr )
828#else
829       CALL posix_lseek( fh, array_position )
830       CALL posix_read( fh, array_3d, SIZE( array_3d ) )
831#endif
832       IF ( include_total_domain_boundaries )  THEN
833          DO  i = lb%nxl, lb%nxr
834             data(:,lb%nys-nbgp:lb%nyn-nbgp,i-nbgp) = array_3d(:,i,lb%nys:lb%nyn)
835          ENDDO
836          IF ( debug_level >= 2 )  WRITE(9,*) 'r3f_ob_soil ', TRIM( name ),' ', SUM( data(:,nys:nyn,nxl:nxr) )
837       ELSE
838          DO  i = nxl, nxr
839             data(:,nys:nyn,i) = array_3d(:,i,nys:nyn)
840          ENDDO
841          IF ( debug_level >= 2 )  WRITE(9,*) 'r3f_soil ', TRIM( name ),' ', SUM( array_3d )
842       ENDIF
843
844    ELSE
845       WRITE(9,*)  'array_3D not found ', name
846       CALL rs_mpi_io_error( 2 )
847    ENDIF
848
849 END SUBROUTINE rrd_mpi_io_real_3d_soil
850
851
852
853!--------------------------------------------------------------------------------------------------!
854! Description:
855! ------------
856!> Read CHARACTER with MPI-IO
857!--------------------------------------------------------------------------------------------------!
858 SUBROUTINE rrd_mpi_io_char( name, text, found )
859
860    IMPLICIT NONE
861
862    CHARACTER(LEN=*), INTENT(IN)   ::  name
863    CHARACTER(LEN=*), INTENT(OUT)  ::  text
864    CHARACTER(LEN=128)             ::  lo_line
865
866    INTEGER(iwp)                   ::  i
867
868    LOGICAL, INTENT(OUT), OPTIONAL ::  found
869    LOGICAL                        ::  lo_found
870
871
872    lo_found = .FALSE.
873    text = ' '
874
875    DO  i = 1, tgh%nr_char
876       lo_line = text_lines(i)
877       IF ( lo_line(1:32) == name )  THEN
878          IF ( debug_level >= 2 )  WRITE(9,*)  'Character variable found ==> ', lo_line(1:32)
879          text = lo_line(33:)
880          lo_found = .TRUE.
881          EXIT
882       ENDIF
883    ENDDO
884
885    IF ( PRESENT( found ) )  THEN
886       found = lo_found
887       RETURN
888    ENDIF
889
890    IF ( .NOT. lo_found )  THEN
891       WRITE(9,*)  'Character variable not found ', name
892         CALL rs_mpi_io_error( 3 )
893    ENDIF
894
895 END SUBROUTINE rrd_mpi_io_char
896
897
898
899!--------------------------------------------------------------------------------------------------!
900! Description:
901! ------------
902!> Read LOGICAL with MPI-IO
903!--------------------------------------------------------------------------------------------------!
904 SUBROUTINE rrd_mpi_io_logical( name, value )
905
906    IMPLICIT NONE
907
908    CHARACTER(LEN=*), INTENT(IN) ::  name
909
910    INTEGER(iwp)                 ::  logical_as_integer
911
912    LOGICAL, INTENT(OUT)         ::  value
913
914
915    CALL rrd_mpi_io_int( name, logical_as_integer )
916    value = ( logical_as_integer == 1 )
917
918 END SUBROUTINE rrd_mpi_io_logical
919
920
921
922!--------------------------------------------------------------------------------------------------!
923! Description:
924! ------------
925!> Write INTEGER with MPI-IO
926!--------------------------------------------------------------------------------------------------!
927 SUBROUTINE wrd_mpi_io_int( name, value )
928
929    IMPLICIT NONE
930
931    CHARACTER(LEN=*), INTENT(IN)  ::  name
932
933    INTEGER(KIND=iwp), INTENT(IN) ::  value
934
935
936    int_names(header_int_index)  = name
937    int_values(header_int_index) = value
938    header_int_index = header_int_index + 1
939
940 END SUBROUTINE wrd_mpi_io_int
941
942
943
944 SUBROUTINE wrd_mpi_io_real( name, value )
945
946    IMPLICIT NONE
947
948    CHARACTER(LEN=*), INTENT(IN) ::  name
949
950    REAL(wp), INTENT(IN)         ::  value
951
952
953    real_names(header_real_index)  = name
954    real_values(header_real_index) = value
955    header_real_index = header_real_index + 1
956
957 END SUBROUTINE wrd_mpi_io_real
958
959
960
961!--------------------------------------------------------------------------------------------------!
962! Description:
963! ------------
964!> Write 2d-REAL array with MPI-IO
965!--------------------------------------------------------------------------------------------------!
966 SUBROUTINE wrd_mpi_io_real_2d( name, data )
967
968    IMPLICIT NONE
969
970    CHARACTER(LEN=*), INTENT(IN)       ::  name
971
972    INTEGER(iwp)                       ::  i
973
974#if defined( __parallel )
975    INTEGER, DIMENSION(rd_status_size) ::  status
976#endif
977
978    REAL(KIND=wp), DIMENSION(lb%nxl:lb%nxr,lb%nys:lb%nyn)  :: array_2d
979
980    REAL(wp), INTENT(IN), DIMENSION(nysg:nyng,nxlg:nxrg)    :: data
981
982
983    array_names(header_arr_index)  = name
984    array_offset(header_arr_index) = array_position
985    header_arr_index = header_arr_index + 1
986
987    IF ( include_total_domain_boundaries )  THEN
988!
989!--    Prepare Output with outer boundaries
990       DO  i = lb%nxl, lb%nxr
991          array_2d(i,lb%nys:lb%nyn) = data(lb%nys-nbgp:lb%nyn-nbgp,i-nbgp)
992       ENDDO
993       IF ( debug_level >= 2 )  WRITE(9,*)  'w2f_ob ', TRIM( name ), ' ',  SUM( data(nys:nyn,nxl:nxr) )
994    ELSE
995!
996!--    Prepare Output without outer boundaries
997       DO  i = nxl,nxr
998          array_2d(i,lb%nys:lb%nyn) = data(nys:nyn,i)
999       ENDDO
1000       IF ( debug_level >= 2 )  WRITE(9,*)  'w2f ', TRIM( name ),' ', SUM( array_2d )
1001    ENDIF
1002
1003#if defined( __parallel )
1004    CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL, ierr )
1005    CALL MPI_FILE_WRITE_ALL( fh, array_2d, SIZE( array_2d), MPI_REAL, status, ierr )
1006#else
1007    CALL posix_lseek( fh, array_position )
1008    CALL posix_write( fh, array_2d, SIZE( array_2d ) )
1009#endif
1010!
1011!-- Type conversion required, otherwise rigth hand side brackets are calculated assuming 4 byte INT.
1012!-- Maybe a compiler problem.
1013    array_position = array_position + ( INT( lb%ny, KIND=rd_offset_kind ) + 1 ) *                  &
1014                                      ( INT( lb%nx, KIND=rd_offset_kind ) + 1 ) * wp
1015
1016 END SUBROUTINE wrd_mpi_io_real_2d
1017
1018
1019
1020!--------------------------------------------------------------------------------------------------!
1021! Description:
1022! ------------
1023!> Write 2d-INTEGER array with MPI-IO
1024!--------------------------------------------------------------------------------------------------!
1025 SUBROUTINE wrd_mpi_io_int_2d( name, data, ar_found )
1026
1027    IMPLICIT NONE
1028
1029    CHARACTER(LEN=*), INTENT(IN)                  ::  name
1030
1031    INTEGER(iwp)                                  ::  i
1032    INTEGER(iwp)                                  ::  j
1033
1034#if defined( __parallel )
1035    INTEGER, DIMENSION(rd_status_size)            ::  status
1036#endif
1037    INTEGER(KIND=iwp), INTENT(IN), DIMENSION(:,:) ::  data
1038    INTEGER, DIMENSION(nxl:nxr,nys:nyn)           ::  array_2d
1039
1040    LOGICAl, OPTIONAL                             ::  ar_found
1041
1042
1043    array_names(header_arr_index)  = name
1044    array_offset(header_arr_index) = array_position
1045    header_arr_index = header_arr_index + 1
1046
1047    IF ( ( nxr-nxl + 1 + 2 * nbgp ) == SIZE( data, 2 ) )  THEN
1048!
1049!--    Integer arrays with ghost layers are not implemented yet. These kind of arrays would be
1050!--    dimensioned in the caller subroutine as
1051!--    INTEGER, DIMENSION(nysg:nyng,nxlg:nxrg) ::  data
1052       WRITE (9,*) '### wrd_mpi_io_int_2d  IF  array: ', TRIM( name )
1053       CALL rs_mpi_io_error( 4 )
1054
1055    ELSEIF ( ( nxr-nxl+1 ) == SIZE( data, 2 ) )  THEN
1056!
1057!--    INTEGER input array without ghost layers.
1058!--    This kind of array is dimensioned in the caller subroutine as
1059!--    INTEGER, DIMENSION(nys:nyn,nxl:nxr) ::  data
1060       DO  j = nys, nyn
1061          DO  i = nxl, nxr
1062             array_2d(i,j) = data(j-nys+1,i-nxl+1)
1063          ENDDO
1064       ENDDO
1065       IF ( debug_level >= 2 )  WRITE(9,*) 'w2i ', TRIM( name ), ' ', SUM( array_2d )
1066#if defined( __parallel )
1067       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native', MPI_INFO_NULL,&
1068                               ierr )
1069       CALL MPI_FILE_WRITE_ALL( fh, array_2d, SIZE( array_2d ), MPI_INTEGER, status, ierr )
1070#else
1071       CALL posix_lseek( fh, array_position )
1072       CALL posix_write( fh, array_2d, SIZE( array_2d ) )
1073#endif
1074!
1075!--    Type conversion required, otherwise rigth hand side brackets are calculated assuming 4 byte
1076!--    INT. Maybe a compiler problem.
1077       array_position = array_position + INT( (ny+1), KIND=rd_offset_kind ) *                      &
1078                                         INT( (nx+1), KIND=rd_offset_kind ) * 4
1079
1080    ELSE
1081       WRITE (9,*) '### wrd_mpi_io_int_2d  array: ', TRIM( name )
1082       CALL rs_mpi_io_error( 4 )
1083    ENDIF
1084
1085    IF ( PRESENT( ar_found ) )  ar_found = .TRUE.
1086
1087 END SUBROUTINE wrd_mpi_io_int_2d
1088
1089
1090
1091!--------------------------------------------------------------------------------------------------!
1092! Description:
1093! ------------
1094!> Write 3d-REAL array with MPI-IO
1095!--------------------------------------------------------------------------------------------------!
1096 SUBROUTINE wrd_mpi_io_real_3d( name, data )
1097
1098    IMPLICIT NONE
1099
1100    CHARACTER(LEN=*), INTENT(IN)       ::  name
1101
1102    INTEGER(iwp)                       ::  i
1103#if defined( __parallel )
1104    INTEGER, DIMENSION(rd_status_size) ::  status
1105#endif
1106    REAL(wp), INTENT(IN), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg)  ::  data
1107
1108    REAL(KIND=wp), DIMENSION(nzb:nzt+1,lb%nxl:lb%nxr,lb%nys:lb%nyn) ::  array_3d
1109
1110
1111    array_names(header_arr_index)  = name
1112    array_offset(header_arr_index) = array_position
1113    header_arr_index = header_arr_index + 1
1114
1115    IF ( include_total_domain_boundaries )  THEN
1116!
1117!--    Prepare output of 3d-REAL-array with ghost layers.
1118!--    In the virtual PE grid, the first dimension is PEs along x, and the second along y.
1119!--    For MPI-IO it is recommended to have the index order of the array in the same way, i.e.
1120!--    the first dimension should be along x and the second along y.
1121!--    For this reason, the original PALM data need to be swaped.
1122       DO  i = lb%nxl, lb%nxr
1123          array_3d(:,i,lb%nys:lb%nyn) = data(:,lb%nys-nbgp:lb%nyn-nbgp,i-nbgp)
1124       ENDDO
1125       IF ( debug_level >= 2 )  WRITE(9,*) 'w3f_ob ', TRIM( name ),' ', SUM( data(:,nys:nyn,nxl:nxr) )
1126    ELSE
1127!
1128!--    Prepare output of 3d-REAL-array without ghost layers
1129       DO  i = nxl, nxr
1130           array_3d(:,i,lb%nys:lb%nyn) = data(:,nys:nyn,i)
1131       ENDDO
1132       IF ( debug_level >= 2 )  WRITE(9,*)  'w3f ', TRIM( name ),' ', SUM( array_3d )
1133    ENDIF
1134#if defined( __parallel )
1135    CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL, ierr )
1136    CALL MPI_FILE_WRITE_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
1137#else
1138    CALL posix_lseek( fh, array_position )
1139    CALL posix_write( fh, array_3d, SIZE( array_3d ) )
1140#endif
1141!
1142!-- Type conversion required, otherwise rigth hand side brackets are calculated assuming 4 byte INT.
1143!-- Maybe a compiler problem.
1144    array_position = array_position + INT(    (nz+2), KIND=rd_offset_kind ) *                      &
1145                                      INT( (lb%ny+1), KIND=rd_offset_kind ) *                      &
1146                                      INT( (lb%nx+1), KIND=rd_offset_kind ) * wp
1147
1148 END SUBROUTINE wrd_mpi_io_real_3d
1149
1150
1151
1152!--------------------------------------------------------------------------------------------------!
1153! Description:
1154! ------------
1155!> Write 3d-REAL soil array with MPI-IO.
1156!> nzb_soil, nzt_soil are located in the module land_surface_model_mod. Since Fortran does not allow
1157!> cross referencing of module variables, it is required to pass these variables as arguments.
1158!--------------------------------------------------------------------------------------------------!
1159 SUBROUTINE wrd_mpi_io_real_3d_soil( name, data, nzb_soil, nzt_soil )
1160
1161    IMPLICIT NONE
1162
1163    CHARACTER(LEN=*), INTENT(IN)       ::  name
1164
1165    INTEGER(iwp)                       ::  i
1166    INTEGER, INTENT(IN)                ::  nzb_soil
1167    INTEGER, INTENT(IN)                ::  nzt_soil
1168
1169#if defined( __parallel )
1170    INTEGER, DIMENSION(rd_status_size) ::  status
1171#endif
1172
1173    REAL(wp), INTENT(IN), DIMENSION(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg)  ::  data
1174
1175    REAL(KIND=wp), DIMENSION(nzb_soil:nzt_soil,lb%nxl:lb%nxr,lb%nys:lb%nyn) ::  array_3d
1176
1177
1178    array_names(header_arr_index)  = name
1179    array_offset(header_arr_index) = array_position
1180    header_arr_index = header_arr_index + 1
1181
1182    IF ( include_total_domain_boundaries)  THEN
1183!
1184!--    Prepare output of 3d-REAL-array with ghost layers.
1185!--    In the virtual PE grid, the first dimension is PEs along x, and the second along y.
1186!--    For MPI-IO it is recommended to have the index order of the array in the same way, i.e.
1187!--    the first dimension should be along x and the second along y.
1188!--    For this reason, the original PALM data need to be swaped.
1189       DO  i = lb%nxl, lb%nxr
1190          array_3d(:,i,lb%nys:lb%nyn) = data(:,lb%nys-nbgp:lb%nyn-nbgp,i-nbgp)
1191       ENDDO
1192       IF ( debug_level >= 2 )  WRITE(9,*) 'w3f_ob_soil ', TRIM( name ), ' ', SUM( data(:,nys:nyn,nxl:nxr) )
1193    ELSE
1194!
1195!--    Prepare output of 3d-REAL-array without ghost layers
1196       DO  i = nxl, nxr
1197          array_3d(:,i,lb%nys:lb%nyn) = data(:,nys:nyn,i)
1198       ENDDO
1199       IF ( debug_level >= 2 )  WRITE(9,*) 'w3f_soil ', TRIM( name ), ' ', SUM( array_3d )
1200    ENDIF
1201#if defined( __parallel )
1202    CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
1203    CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3dsoil, 'native', MPI_INFO_NULL, ierr )
1204    CALL MPI_FILE_WRITE_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
1205    CALL MPI_TYPE_FREE( ft_3dsoil, ierr )
1206#else
1207    CALL posix_lseek( fh, array_position )
1208    CALL posix_write( fh, array_3d, SIZE( array_3d ) )
1209#endif
1210!
1211!-- Type conversion required, otherwise rigth hand side brackets are calculated assuming 4 byte INT.
1212!-- Maybe a compiler problem.
1213    array_position = array_position + INT( (nzt_soil-nzb_soil+1), KIND=rd_offset_kind ) *          &
1214                                      INT( (lb%ny+1),             KIND=rd_offset_kind ) *          &
1215                                      INT( (lb%nx+1),             KIND=rd_offset_kind ) * wp
1216
1217 END SUBROUTINE wrd_mpi_io_real_3d_soil
1218
1219
1220
1221!--------------------------------------------------------------------------------------------------!
1222! Description:
1223! ------------
1224!> Write CHARATCTER with MPI-IO
1225!--------------------------------------------------------------------------------------------------!
1226 SUBROUTINE wrd_mpi_io_char( name, text )
1227
1228    IMPLICIT NONE
1229
1230    CHARACTER(LEN=128)           ::  lo_line
1231    CHARACTER(LEN=*), INTENT(IN) ::  name
1232    CHARACTER(LEN=*), INTENT(IN) ::  text
1233
1234
1235    lo_line      = name
1236    lo_line(33:) = text
1237    text_lines(header_char_index) = lo_line
1238    header_char_index = header_char_index + 1
1239
1240 END SUBROUTINE wrd_mpi_io_char
1241
1242
1243
1244!--------------------------------------------------------------------------------------------------!
1245! Description:
1246! ------------
1247!> Write LOGICAL with MPI-IO
1248!--------------------------------------------------------------------------------------------------!
1249 SUBROUTINE wrd_mpi_io_logical( name, value )
1250
1251    IMPLICIT NONE
1252
1253    CHARACTER(LEN=*), INTENT(IN) ::  name
1254
1255    INTEGER(iwp)                 ::  logical_as_integer
1256
1257    LOGICAL, INTENT(IN)          ::  value
1258
1259
1260    IF ( value )  THEN
1261       logical_as_integer = 1
1262    ELSE
1263       logical_as_integer = 0
1264    ENDIF
1265
1266    CALL wrd_mpi_io_int( name, logical_as_integer )
1267
1268 END SUBROUTINE wrd_mpi_io_logical
1269
1270
1271
1272!--------------------------------------------------------------------------------------------------!
1273! Description:
1274! ------------
1275!> Read 1d-REAL global array with MPI-IO.
1276!> Array contains identical data on all PEs.
1277!--------------------------------------------------------------------------------------------------!
1278 SUBROUTINE rrd_mpi_io_global_array_real_1d( name, data )
1279
1280    IMPLICIT NONE
1281
1282    CHARACTER(LEN=*), INTENT(IN)       ::  name
1283
1284    INTEGER(iwp)                       ::  i
1285    INTEGER(KIND=rd_offset_kind)       ::  offset
1286
1287#if defined( __parallel )
1288    INTEGER, DIMENSION(rd_status_size) ::  status
1289#endif
1290
1291    LOGICAL                            ::  found
1292
1293    REAL(KIND=wp), INTENT(INOUT), DIMENSION(:) ::  data
1294
1295
1296    offset = 0
1297    found  = .FALSE.
1298
1299    DO  i = 1, tgh%nr_arrays
1300       IF ( TRIM(array_names(i)) == TRIM( name ) )  THEN
1301          array_position = array_offset(i)
1302          found = .TRUE.
1303          EXIT
1304       ENDIF
1305    ENDDO
1306
1307    IF ( found )  THEN
1308!
1309!--    Set default view
1310#if defined( __parallel )
1311       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
1312       CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
1313       CALL MPI_FILE_READ_ALL( fh, data, SIZE( data ), MPI_REAL, status, ierr )
1314#else
1315       CALL posix_lseek( fh, array_position )
1316       CALL posix_read( fh, data, SIZE( data ) )
1317#endif
1318       IF ( debug_level >= 2) WRITE(9,*) 'rr1f ',name,' ', SUM( data)
1319    ELSE
1320       WRITE(9,*)  'replicated array_1D not found ', name
1321       CALL rs_mpi_io_error( 2 )
1322    ENDIF
1323
1324 END SUBROUTINE rrd_mpi_io_global_array_real_1d
1325
1326
1327
1328!--------------------------------------------------------------------------------------------------!
1329! Description:
1330! ------------
1331!> Read 2d-REAL global array with MPI-IO.
1332!> Array contains identical data on all PEs.
1333!--------------------------------------------------------------------------------------------------!
1334 SUBROUTINE rrd_mpi_io_global_array_real_2d( name, data )
1335
1336    IMPLICIT NONE
1337
1338    CHARACTER(LEN=*), INTENT(IN)                      ::  name
1339
1340    INTEGER, DIMENSION(1)                             ::  bufshape
1341
1342    REAL(KIND=wp), INTENT(IN), DIMENSION(:,:), TARGET ::  data
1343    REAL(KIND=wp), POINTER, DIMENSION(:)              ::  buf
1344
1345    TYPE(C_PTR)                                       ::  c_data
1346
1347
1348    c_data = C_LOC( data )
1349    bufshape(1) = SIZE( data )
1350    CALL C_F_POINTER( c_data, buf, bufshape )
1351
1352    CALL rrd_mpi_io_global_array_real_1d( name, buf )
1353    IF ( debug_level >= 2 )  WRITE(9,*) 'rr2f ', TRIM( name ), ' ', bufshape(1), SUM( data )
1354
1355 END SUBROUTINE rrd_mpi_io_global_array_real_2d
1356
1357
1358
1359!--------------------------------------------------------------------------------------------------!
1360! Description:
1361! ------------
1362!> Read 3d-REAL global array with MPI-IO.
1363!> Array contains identical data on all PEs.
1364!--------------------------------------------------------------------------------------------------!
1365 SUBROUTINE rrd_mpi_io_global_array_real_3d( name, data )
1366
1367    IMPLICIT NONE
1368
1369    CHARACTER(LEN=*), INTENT(IN)                        ::  name
1370
1371    INTEGER, DIMENSION(1)                               ::  bufshape
1372
1373    REAL(KIND=wp), INTENT(IN), DIMENSION(:,:,:), TARGET ::  data
1374    REAL(KIND=wp), POINTER, DIMENSION(:)                ::  buf
1375
1376    TYPE(C_PTR)                                         ::  c_data
1377
1378
1379    c_data = C_LOC( data )
1380    bufshape(1) = SIZE( data )
1381    CALL C_F_POINTER( c_data, buf, bufshape )
1382
1383    CALL rrd_mpi_io_global_array_real_1d( name, buf )
1384
1385    IF ( debug_level >= 2 )  WRITE(9,*) 'rr3f ', TRIM( name ), ' ', bufshape(1), SUM( data )
1386
1387 END SUBROUTINE rrd_mpi_io_global_array_real_3d
1388
1389
1390
1391!--------------------------------------------------------------------------------------------------!
1392! Description:
1393! ------------
1394!> Read 4d-REAL global array with MPI-IO.
1395!> Array contains identical data on all PEs.
1396!--------------------------------------------------------------------------------------------------!
1397 SUBROUTINE rrd_mpi_io_global_array_real_4d( name, data )
1398
1399    IMPLICIT NONE
1400
1401    CHARACTER(LEN=*), INTENT(IN)                          ::  name
1402
1403    INTEGER, DIMENSION(1)                                 ::  bufshape
1404
1405    REAL(KIND=wp), INTENT(IN), DIMENSION(:,:,:,:), TARGET ::  data
1406    REAL(KIND=wp), POINTER, DIMENSION(:)                  ::  buf
1407
1408    TYPE(C_PTR)                                           ::  c_data
1409
1410
1411    c_data = C_LOC( data )
1412    bufshape(1) = SIZE( data)
1413    CALL C_F_POINTER( c_data, buf, bufshape )
1414
1415    CALL rrd_mpi_io_global_array_real_1d( name, buf )
1416    IF ( debug_level >= 2 )  WRITE(9,*) 'rr4f ', TRIM( name ), ' ', bufshape(1), SUM( data )
1417
1418 END SUBROUTINE rrd_mpi_io_global_array_real_4d
1419
1420
1421
1422!--------------------------------------------------------------------------------------------------!
1423! Description:
1424! ------------
1425!> Read 1d-INTEGER global array with MPI-IO.
1426!> Array contains identical data on all PEs.
1427!--------------------------------------------------------------------------------------------------!
1428 SUBROUTINE rrd_mpi_io_global_array_int_1d( name, data, ar_found )
1429
1430    IMPLICIT NONE
1431
1432    CHARACTER(LEN=*), INTENT(IN)                   ::  name
1433
1434    INTEGER(iwp)                                   ::  i
1435    INTEGER(KIND=rd_offset_kind)                   ::  offset
1436
1437#if defined( __parallel )
1438    INTEGER, DIMENSION(rd_status_size)             ::  status
1439#endif
1440    INTEGER(KIND=iwp), INTENT(INOUT), DIMENSION(:) ::  data
1441
1442    LOGICAl, OPTIONAL                              ::  ar_found
1443    LOGICAL                                        ::  found
1444
1445
1446    offset = 0
1447    found  = .FALSE.
1448
1449    DO  i = 1, tgh%nr_arrays
1450       IF ( TRIM(array_names(i)) == TRIM( name ) )  THEN
1451          array_position = array_offset(i)
1452          found = .TRUE.
1453          EXIT
1454       ENDIF
1455    ENDDO
1456
1457    IF ( found )  THEN
1458!
1459!--    Set default view
1460#if defined( __parallel )
1461       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
1462       CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
1463       CALL MPI_FILE_READ_ALL( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
1464#else
1465       CALL posix_lseek( fh, array_position )
1466       CALL posix_read( fh, data, SIZE( data ) )
1467#endif
1468    ELSE
1469       IF ( PRESENT( ar_found ) )  THEN
1470          ar_found =.FALSE.
1471          RETURN
1472       ELSE
1473          WRITE (9,*) '### rrd_mpi_io_global_array_int_1d ', TRIM( name )
1474          CALL rs_mpi_io_error( 4 )
1475          WRITE(9,*)  'replicated array_1D not found ', name
1476          CALL rs_mpi_io_error( 2 )
1477       ENDIF
1478    ENDIF
1479
1480    IF ( PRESENT( ar_found ) )  ar_found =.TRUE.
1481
1482 END SUBROUTINE rrd_mpi_io_global_array_int_1d
1483
1484
1485
1486!--------------------------------------------------------------------------------------------------!
1487! Description:
1488! ------------
1489!> Write 1d-REAL global array with MPI-IO.
1490!> Array contains identical data on all PEs.
1491!--------------------------------------------------------------------------------------------------!
1492 SUBROUTINE wrd_mpi_io_global_array_real_1d( name, data )
1493
1494    IMPLICIT NONE
1495
1496    CHARACTER(LEN=*), INTENT(IN)            ::  name
1497
1498    INTEGER(KIND=rd_offset_kind)            ::  offset
1499
1500#if defined( __parallel )
1501    INTEGER, DIMENSION(rd_status_size)      ::  status
1502#endif
1503
1504    REAL(KIND=wp), INTENT(IN), DIMENSION(:) ::  data
1505
1506
1507    offset = 0
1508
1509    array_names(header_arr_index)  = name
1510    array_offset(header_arr_index) = array_position
1511    header_arr_index = header_arr_index + 1
1512
1513    IF ( debug_level >= 2 )  WRITE(9,*)  'wr1f ', TRIM( name ), ' ', SUM( data )
1514!
1515!--    Set default view
1516#if defined( __parallel )
1517       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
1518!
1519!--    Only PE 0 writes replicated data
1520       IF ( myid == 0 )  THEN
1521          CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
1522          CALL MPI_FILE_WRITE( fh, data, SIZE( data), MPI_REAL, status, ierr )
1523       ENDIF
1524#else
1525       CALL posix_lseek( fh, array_position )
1526       CALL posix_write( fh, data, SIZE( data ) )
1527#endif
1528       array_position = array_position + SIZE( data ) * wp
1529
1530 END SUBROUTINE wrd_mpi_io_global_array_real_1d
1531
1532
1533
1534!--------------------------------------------------------------------------------------------------!
1535! Description:
1536! ------------
1537!> Write 2d-REAL global array with MPI-IO.
1538!> Array contains identical data on all PEs.
1539!--------------------------------------------------------------------------------------------------!
1540 SUBROUTINE wrd_mpi_io_global_array_real_2d( name, data )
1541
1542    IMPLICIT NONE
1543
1544    CHARACTER(LEN=*), INTENT(IN)                      ::  name
1545
1546    INTEGER, DIMENSION(1)                             ::  bufshape
1547
1548    REAL(KIND=wp), POINTER, DIMENSION(:)              ::  buf
1549    REAL(KIND=wp), INTENT(IN), DIMENSION(:,:), TARGET ::  data
1550
1551    TYPE(C_PTR)                                       ::  c_data
1552
1553
1554    c_data = C_LOC( data )
1555    bufshape(1) = SIZE( data)
1556    CALL C_F_POINTER( c_data, buf, bufshape )
1557
1558    IF ( debug_level >= 2 )  WRITE(9,*)  'wr2f ', TRIM( name ), ' ', bufshape(1), SUM( data )
1559
1560    CALL wrd_mpi_io_global_array_real_1d( name, buf )
1561
1562 END SUBROUTINE wrd_mpi_io_global_array_real_2d
1563
1564
1565
1566!--------------------------------------------------------------------------------------------------!
1567! Description:
1568! ------------
1569!> Write 3d-REAL global array with MPI-IO.
1570!> Array contains identical data on all PEs.
1571!--------------------------------------------------------------------------------------------------!
1572 SUBROUTINE wrd_mpi_io_global_array_real_3d( name, data )
1573
1574    IMPLICIT NONE
1575
1576    CHARACTER(LEN=*), INTENT(IN)                        ::  name
1577
1578    INTEGER, DIMENSION(1)                               ::  bufshape
1579
1580    REAL(KIND=wp), POINTER, DIMENSION(:)                ::  buf
1581    REAL(KIND=wp), INTENT(IN), DIMENSION(:,:,:), TARGET ::  data
1582
1583    TYPE(C_PTR)                                         ::  c_data
1584
1585
1586    c_data = C_LOC( data )
1587    bufshape(1) = SIZE( data )
1588    CALL C_F_POINTER( c_data, buf, bufshape )
1589
1590    IF ( debug_level >= 2 )  WRITE(9,*)  'wr3f ', TRIM( name ), ' ', bufshape(1), SUM( data )
1591
1592    CALL wrd_mpi_io_global_array_real_1d( name, buf )
1593
1594 END SUBROUTINE wrd_mpi_io_global_array_real_3d
1595
1596
1597
1598!--------------------------------------------------------------------------------------------------!
1599! Description:
1600! ------------
1601!> Write 4d-REAL global array with MPI-IO.
1602!> Array contains identical data on all PEs.
1603!--------------------------------------------------------------------------------------------------!
1604 SUBROUTINE wrd_mpi_io_global_array_real_4d( name, data )
1605
1606    IMPLICIT NONE
1607
1608    CHARACTER(LEN=*), INTENT(IN)                          ::  name
1609
1610    INTEGER, DIMENSION(1)                                 ::  bufshape
1611
1612    REAL(KIND=wp), POINTER, DIMENSION(:)                  ::  buf
1613    REAL(KIND=wp), INTENT(IN), DIMENSION(:,:,:,:), TARGET ::  data
1614
1615    TYPE(C_PTR)                                           ::  c_data
1616
1617
1618    c_data = C_LOC( data )
1619    bufshape(1) = SIZE( data)
1620    CALL C_F_POINTER( c_data, buf, bufshape )
1621
1622    IF ( debug_level >= 2 )  WRITE(9,*) 'wr4f ', TRIM( name ), ' ', bufshape(1), SUM( data )
1623
1624    CALL wrd_mpi_io_global_array_real_1d( name, buf )
1625
1626 END SUBROUTINE wrd_mpi_io_global_array_real_4d
1627
1628
1629
1630!--------------------------------------------------------------------------------------------------!
1631! Description:
1632! ------------
1633!> Write 1d-INTEGER global array with MPI-IO.
1634!> Array contains identical data on all PEs.
1635!--------------------------------------------------------------------------------------------------!
1636 SUBROUTINE wrd_mpi_io_global_array_int_1d( name, data )
1637
1638    IMPLICIT NONE
1639
1640    CHARACTER(LEN=*), INTENT(IN)                ::  name
1641
1642    INTEGER(KIND=rd_offset_kind)                ::  offset
1643
1644    INTEGER(KIND=iwp), INTENT(IN), DIMENSION(:) ::  data
1645#if defined( __parallel )
1646    INTEGER, DIMENSION(rd_status_size)          ::  status
1647#endif
1648
1649    offset = 0
1650    array_names(header_arr_index)  = name
1651    array_offset(header_arr_index) = array_position
1652    header_arr_index = header_arr_index + 1
1653
1654    IF ( debug_level >= 2 )  WRITE(9,*)  'wr1i ', TRIM( name ), ' ', SUM( data )
1655!
1656!-- Set default view
1657#if defined( __parallel )
1658    CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
1659!
1660!-- Only PE 0 writes replicated data
1661    IF ( myid == 0 )  THEN                        !
1662       CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
1663       CALL MPI_FILE_WRITE( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
1664    ENDIF
1665#else
1666    CALL posix_lseek( fh, array_position )
1667    CALL posix_write( fh, data, SIZE( data ) )
1668#endif
1669    array_position = array_position + SIZE( data ) * 4
1670
1671 END SUBROUTINE wrd_mpi_io_global_array_int_1d
1672
1673
1674
1675!--------------------------------------------------------------------------------------------------!
1676! Description:
1677! ------------
1678!> Read 1d-REAL surface data array with MPI-IO.
1679!--------------------------------------------------------------------------------------------------!
1680 SUBROUTINE rrd_mpi_io_surface( name, data, first_index )
1681
1682    IMPLICIT NONE
1683
1684    CHARACTER(LEN=*), INTENT(IN) ::  name
1685
1686    INTEGER(KIND=rd_offset_kind) ::  disp          !< displacement of actual indices
1687    INTEGER(KIND=rd_offset_kind) ::  disp_f        !< displacement in file
1688    INTEGER(KIND=rd_offset_kind) ::  disp_n        !< displacement of next column
1689    INTEGER(iwp), OPTIONAL       ::  first_index
1690
1691    INTEGER(iwp)                 ::  i
1692    INTEGER(iwp)                 ::  i_f
1693    INTEGER(iwp)                 ::  j
1694    INTEGER(iwp)                 ::  j_f
1695    INTEGER(iwp)                 ::  lo_first_index
1696    INTEGER(iwp)                 ::  nr_bytes
1697    INTEGER(iwp)                 ::  nr_bytes_f
1698    INTEGER(iwp)                 ::  nr_words
1699#if defined( __parallel )
1700    INTEGER, DIMENSION(rd_status_size)  ::  status
1701#endif
1702
1703    LOGICAL                             ::  found
1704
1705    REAL(wp), INTENT(OUT), DIMENSION(:) ::  data
1706
1707
1708    found = .FALSE.
1709    lo_first_index = 1
1710
1711    IF ( MAXVAL( m_global_start ) == -1 )   RETURN   ! nothing to do on this PE
1712
1713    IF ( PRESENT( first_index ) )   THEN
1714       lo_first_index = first_index
1715    ENDIF
1716
1717    DO  i = 1, tgh%nr_arrays
1718        IF ( TRIM( array_names(i) ) == TRIM( name ) )  THEN
1719           array_position = array_offset(i) + ( lo_first_index - 1 ) *                             &
1720                                              total_number_of_surface_values * wp
1721           found = .TRUE.
1722           EXIT
1723        ENDIF
1724    ENDDO
1725
1726    disp   = -1
1727    disp_f = -1
1728    disp_n = -1
1729    IF ( found )  THEN
1730
1731       DO  i = nxl, nxr
1732          DO  j = nys, nyn
1733
1734             IF ( m_global_start(j,i) > 0 )  THEN
1735                disp     = array_position+(m_global_start(j,i)-1) * wp
1736                nr_words = m_end_index(j,i)-m_start_index(j,i)+1
1737                nr_bytes = nr_words * wp
1738             ENDIF
1739             IF ( disp >= 0  .AND.  disp_f == -1 )  THEN   ! first Entry
1740                disp_f     = disp
1741                nr_bytes_f = 0
1742                i_f = i
1743                j_f = j
1744             ENDIF
1745             IF ( j == nyn  .AND.  i == nxr )  THEN        ! last Entry
1746                disp_n = -1
1747                IF (  nr_bytes > 0 )  THEN
1748                   nr_bytes_f = nr_bytes_f+nr_bytes
1749                ENDIF
1750             ELSEIF ( j == nyn )  THEN                     ! next x
1751                IF ( m_global_start(nys,i+1) > 0  .AND.  disp > 0 )  THEN
1752                   disp_n = array_position + ( m_global_start(nys,i+1) - 1 ) * wp
1753                ELSE
1754                   CYCLE
1755                ENDIF
1756             ELSE
1757                IF ( m_global_start(j+1,i) > 0  .AND.  disp > 0 )  THEN
1758                   disp_n = array_position + ( m_global_start(j+1,i) - 1 ) * wp
1759                ELSE
1760                   CYCLE
1761                ENDIF
1762             ENDIF
1763
1764
1765             IF ( disp + nr_bytes == disp_n )  THEN        ! contiguous block
1766                nr_bytes_f = nr_bytes_f + nr_bytes
1767             ELSE                                          ! read
1768#if defined( __parallel )
[4500]1769                IF ( debug_level >= 2 )  WRITE(9,'(a,8i8)') 'read block ', j, i, j_f, i_f, &
1770                                                            m_start_index(j_f,i_f), nr_bytes_f, disp_f
[4495]1771                CALL MPI_FILE_SEEK( fh, disp_f, MPI_SEEK_SET, ierr )
1772                nr_words = nr_bytes_f / wp
1773                CALL MPI_FILE_READ( fh, data(m_start_index(j_f,i_f)), nr_words, MPI_REAL, status, ierr )
1774#else
1775                CALL posix_lseek( fh, disp_f )
[4497]1776!                CALL posix_read( fh, data(m_start_index(j_f:,i_f:)), nr_bytes_f )
[4495]1777#endif
1778                disp_f     = disp
1779                nr_bytes_f = nr_bytes
1780                i_f = i
1781                j_f = j
1782             ENDIF
1783
1784          ENDDO
1785       ENDDO
1786
1787    ELSE
1788       WRITE(9,*) 'surface array not found ', name
1789       CALL rs_mpi_io_error( 2 )
1790    ENDIF
1791
1792      IF ( lo_first_index == 1 )  THEN
1793         IF ( debug_level >= 2 .AND. nr_val > 0 )  WRITE(9,*)  'r_surf ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
1794      ELSE
[4500]1795         IF ( debug_level >= 2 .AND. nr_val > 0 )  WRITE(9,*)  'r_surf_next ', TRIM( name ), ' ', &
1796                                                               lo_first_index,nr_val, SUM( data(1:nr_val) )
[4495]1797      ENDIF
1798
1799 END SUBROUTINE rrd_mpi_io_surface
1800
1801
1802
1803!--------------------------------------------------------------------------------------------------!
1804! Description:
1805! ------------
1806!> Read 2d-REAL surface data array with MPI-IO.
1807!> These consist of multiple 1d-REAL surface data arrays.
1808!--------------------------------------------------------------------------------------------------!
1809 SUBROUTINE rrd_mpi_io_surface_2d( name, data )
1810
1811    IMPLICIT NONE
1812
1813    CHARACTER(LEN=*), INTENT(IN)          ::  name
1814
1815    INTEGER(iwp)                          ::  i
1816
1817    REAL(wp), INTENT(OUT), DIMENSION(:,:) ::  data
1818    REAL(wp), DIMENSION(SIZE( data,2))    ::  tmp
1819
1820
1821    DO  i = 1, SIZE( data,1)
1822       CALL rrd_mpi_io_surface( name, tmp, first_index = i )
1823       data(i,:) = tmp
1824    ENDDO
1825
1826!
1827!-- Comment from Klaus Ketelsen (September 2018)
1828!-- The intention of the following loop was let the compiler do the copying on return.
1829!-- In my understanding is it standard conform to pass the second dimension to a 1d-
1830!-- array inside a subroutine and the compiler is responsible to generate code for
1831!-- copying. Acually this works fine for INTENT(IN) variables (wrd_mpi_io_surface_2d).
1832!-- For INTENT(OUT) like in this case the code works on pgi compiler. But both, the Intel 16
1833!-- and the Cray compiler show wrong answers using this loop.
1834!-- That is the reason why the above auxiliary array tmp was introduced.
1835!    DO  i = 1, SIZE(  data,1)
1836!       CALL rrd_mpi_io_surface( name, data(i,:), first_index = i )
1837!    ENDDO
1838
1839 END SUBROUTINE rrd_mpi_io_surface_2d
1840
1841
1842
1843!--------------------------------------------------------------------------------------------------!
1844! Description:
1845! ------------
1846!> Write 1d-REAL surface data array with MPI-IO.
1847!--------------------------------------------------------------------------------------------------!
1848 SUBROUTINE wrd_mpi_io_surface( name, data, first_index )
1849
1850    IMPLICIT NONE
1851
1852    CHARACTER(LEN=*), INTENT(IN)       ::  name
1853
1854#if defined( __parallel )
1855    INTEGER(KIND=rd_offset_kind)       ::  disp
1856#endif
1857    INTEGER(iwp), OPTIONAL             ::  first_index
1858    INTEGER(iwp)                       ::  lo_first_index
1859    INTEGER(KIND=rd_offset_kind)       ::  offset
1860
1861#if defined( __parallel )
1862    INTEGER, DIMENSION(rd_status_size) ::  status
1863#endif
1864
1865    REAL(wp), INTENT(IN), DIMENSION(:) ::  data
1866
1867
1868    offset = 0
1869    lo_first_index = 1
1870
1871    IF ( PRESENT(first_index) )  THEN
1872       lo_first_index = first_index
1873    ENDIF
1874!
1875!-- In case of 2d-data, name is writen only once
1876    IF ( lo_first_index == 1 )  THEN
1877       array_names(header_arr_index)  = name
1878       array_offset(header_arr_index) = array_position
1879       header_arr_index = header_arr_index + 1
1880    ENDIF
1881#if defined( __parallel )
1882    IF ( all_pes_write )  THEN
1883       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_surf, 'native', MPI_INFO_NULL, ierr )
1884       CALL MPI_FILE_WRITE_ALL( fh, data, nr_val, MPI_REAL, status, ierr )
1885    ELSE
1886       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
1887       IF ( nr_val > 0 )  THEN
1888          disp = array_position + 8 * ( glo_start - 1 )
1889          CALL MPI_FILE_SEEK( fh, disp, MPI_SEEK_SET, ierr )
1890          CALL MPI_FILE_WRITE( fh, data, nr_val, MPI_REAL, status, ierr )
1891       ENDIF
1892    ENDIF
1893#else
1894    CALL posix_lseek( fh, array_position )
1895    CALL posix_write( fh, data, nr_val )
1896#endif
1897    array_position = array_position + total_number_of_surface_values * wp
1898
1899    IF ( lo_first_index == 1 )  THEN
1900       IF ( debug_level >= 2 .AND. nr_val  > 0 )  WRITE(9,*) 'w_surf ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
1901    ELSE
[4500]1902       IF ( debug_level >= 2 .AND. nr_val  > 0 ) WRITE(9,*) 'w_surf_n ', TRIM( name ), ' ', &
1903                                                            lo_first_index, nr_val, SUM( data(1:nr_val) )
[4495]1904    ENDIF
1905
1906 END SUBROUTINE wrd_mpi_io_surface
1907
1908
1909
1910!--------------------------------------------------------------------------------------------------!
1911! Description:
1912! ------------
1913!> Read 2d-REAL surface data array with MPI-IO.
1914!> These consist of multiple 1d-REAL surface data arrays.
1915!--------------------------------------------------------------------------------------------------!
1916 SUBROUTINE wrd_mpi_io_surface_2d( name, data )
1917
1918    IMPLICIT NONE
1919
1920    CHARACTER(LEN=*), INTENT(IN)         ::  name
1921
1922    INTEGER(iwp)                         ::  i
1923
1924    REAL(wp), INTENT(IN), DIMENSION(:,:) ::  data
1925
1926
1927    DO  i = 1, SIZE( data,1)
1928       CALL wrd_mpi_io_surface( name, data(i,:), first_index = i )
1929    ENDDO
1930
1931 END SUBROUTINE wrd_mpi_io_surface_2d
1932
1933
1934
1935!--------------------------------------------------------------------------------------------------!
1936! Description:
1937! ------------
1938!> Close restart file for MPI-IO
1939!--------------------------------------------------------------------------------------------------!
1940 SUBROUTINE rd_mpi_io_close
1941
1942    IMPLICIT NONE
1943
1944    INTEGER(iwp)                       ::  gh_size
1945    INTEGER(KIND=rd_offset_kind)       ::  offset
1946#if defined( __parallel )
1947    INTEGER, DIMENSION(rd_status_size) ::  status
1948#endif
1949
1950#if ! defined( __parallel )
1951    TYPE(C_PTR)                        ::  buf_ptr
1952#endif
1953
1954
1955    offset = 0
1956
1957    IF ( wr_flag )  THEN
1958
1959       tgh%nr_int    = header_int_index - 1
1960       tgh%nr_char   = header_char_index - 1
1961       tgh%nr_real   = header_real_index - 1
1962       tgh%nr_arrays = header_arr_index - 1
1963       tgh%total_nx  = lb%nx + 1
1964       tgh%total_ny  = lb%ny + 1
1965       IF ( include_total_domain_boundaries )  THEN   ! not sure, if LOGICAL interpretation is the same on all compilers,
1966          tgh%i_outer_bound = 1        ! therefore store as INTEGER in general header
1967       ELSE
1968          tgh%i_outer_bound = 0
1969       ENDIF
1970!
1971!--    Check for big/little endian format. This check is currently not used, and could be removed
1972!--    if we can assume little endian as the default on all machines.
1973       CALL rs_mpi_io_check_endian( tgh%endian )
1974
1975!
1976!--    Set default view
1977#if defined( __parallel )
1978       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
1979#endif
1980!
1981!--    Write header file
1982       gh_size = storage_size(tgh) / 8
1983       IF ( myid == 0 )  THEN
1984#if defined( __parallel )
1985          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
1986          CALL MPI_FILE_WRITE( fh, tgh, gh_size, MPI_INT, status, ierr )
1987          header_position = header_position + gh_size
1988!
1989!--       INTEGER values
1990          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
1991          CALL MPI_FILE_WRITE( fh, int_names, SIZE( int_names )*32, MPI_CHAR, status, ierr )
1992          header_position = header_position + SIZE( int_names ) * 32
1993
1994          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
1995          CALL MPI_FILE_WRITE( fh, int_values, SIZE( int_values ), MPI_INT, status, ierr )
1996          header_position = header_position + SIZE( int_values ) * iwp
1997!
1998!--       Character entries
1999          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
2000          CALL MPI_FILE_WRITE( fh, text_lines, SIZE( text_lines )*128, MPI_CHAR, status, ierr )
2001          header_position = header_position + SIZE( text_lines ) * 128
2002!
2003!---      REAL values
2004          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
2005          CALL MPI_FILE_WRITE( fh, real_names, SIZE( real_names )*32, MPI_CHAR, status, ierr )
2006          header_position = header_position + SIZE( real_names ) * 32
2007
2008          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
2009          CALL MPI_FILE_WRITE( fh, real_values, SIZE( real_values ), MPI_REAL, status, ierr )
2010          header_position = header_position + SIZE( real_values ) * wp
2011!
2012!--       2d- and 3d- distributed array headers, all replicated array headers
2013          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
2014          CALL MPI_FILE_WRITE( fh, array_names, SIZE( array_names )*32, MPI_CHAR, status, ierr )
2015          header_position = header_position + SIZE( array_names ) * 32
2016
2017          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
2018          CALL MPI_FILE_WRITE( fh, array_offset, SIZE( array_offset )*MPI_OFFSET_KIND, MPI_BYTE, status, ierr )   !There is no I*8 datatype in FORTRAN
2019          header_position = header_position + SIZE( array_offset ) * rd_offset_kind
2020#else
2021          CALL posix_lseek( fh, header_position )
2022          buf_ptr = C_LOC( tgh )
2023          CALL posix_write( fh, buf_ptr, gh_size )
2024          header_position = header_position + gh_size
2025!
2026!--       INTEGER values
2027          CALL posix_lseek( fh, header_position )
2028          CALL posix_write( fh, int_names )
2029          header_position = header_position + SIZE( int_names ) * 32
2030
2031          CALL posix_lseek( fh, header_position )
2032          CALL posix_write( fh, int_values, SIZE( int_values ) )
2033          header_position = header_position + SIZE( int_values ) * iwp
2034!
2035!--       Character entries
2036          CALL posix_lseek( fh, header_position )
2037          CALL posix_write( fh, text_lines )
2038          header_position = header_position + SIZE( text_lines ) * 128
2039!
2040!--       REAL values
2041          CALL posix_lseek( fh, header_position )
2042          CALL posix_write( fh, real_names )
2043          header_position = header_position + SIZE( real_names ) * 32
2044
2045          CALL posix_lseek( fh, header_position )
2046          CALL posix_write( fh, real_values, SIZE( real_values ) )
2047          header_position = header_position + SIZE( real_values ) * wp
2048!
2049!--       2d- and 3d-distributed array headers, all replicated array headers
2050          CALL posix_lseek( fh, header_position )
2051          CALL posix_write( fh, array_names )
2052          header_position = header_position + SIZE( array_names ) * 32
2053
2054          CALL posix_lseek( fh, header_position )
2055          CALL posix_write( fh, array_offset, SIZE( array_offset ) )
2056          header_position = header_position + SIZE( array_offset ) * rd_offset_kind
2057#endif
2058          IF ( debug_level >= 2 )  THEN
2059             WRITE(9,*)  'header position after arrays  ', header_position, gh_size
2060          ENDIF
2061
2062          IF ( print_header_now )  CALL rs_mpi_io_print_header
2063       ENDIF
2064
2065    ENDIF
2066
2067!
2068!-- Free file types
2069    CALL rs_mpi_io_free_filetypes
2070
2071!
2072!-- Close MPI-IO file
2073#if defined( __parallel )
2074    CALL MPI_FILE_CLOSE( fh, ierr )
2075#else
2076    CALL posix_close( fh )
2077#endif
2078
2079    mb_processed = array_position / ( 1024.0_dp * 1024.0_dp )
2080
2081 END SUBROUTINE rd_mpi_io_close
2082
2083
2084
2085!--------------------------------------------------------------------------------------------------!
2086! Description:
2087! ------------
2088!> This subroutine prepares a filetype and some variables for the I/O of surface data.
2089!> Whenever a new set of start_index and end_index is used, rd_mpi_io_surface_filetypes has to be
2090!> called. A main feature of this subroutine is computing the global start indices of the 1d- and
2091!> 2d- surface arrays.
2092!--------------------------------------------------------------------------------------------------!
2093 SUBROUTINE rd_mpi_io_surface_filetypes( start_index, end_index, data_to_write, global_start )
2094
2095    IMPLICIT NONE
2096
2097    INTEGER(iwp)                          ::  i            !<  loop index
2098    INTEGER(KIND=rd_offset_kind)          ::  offset
2099
2100    INTEGER(iwp), DIMENSION(1)            ::  dims1
2101    INTEGER(iwp), DIMENSION(1)            ::  lize1
2102    INTEGER(iwp), DIMENSION(1)            ::  start1
2103    INTEGER(iwp), DIMENSION(0:numprocs-1) ::  lo_nr_val    !< local number of values in x and y direction
2104    INTEGER(iwp), DIMENSION(0:numprocs-1) ::  all_nr_val   !< number of values for all PEs
2105
2106    INTEGER, INTENT(IN), DIMENSION(nys:nyn,nxl:nxr)  ::  end_index
2107    INTEGER, INTENT(OUT), DIMENSION(nys:nyn,nxl:nxr) ::  global_start
2108    INTEGER, INTENT(IN), DIMENSION(nys:nyn,nxl:nxr)  ::  start_index
2109
2110    LOGICAL, INTENT(OUT)                  ::  data_to_write  !< returns, if surface data have to be written
2111
2112
2113    offset = 0
2114    lo_nr_val= 0
2115    lo_nr_val(myid) = MAXVAL( end_index )
2116#if defined( __parallel )
2117    CALL MPI_ALLREDUCE( lo_nr_val, all_nr_val, numprocs, MPI_INTEGER, MPI_SUM, comm2d, ierr )
2118    IF ( ft_surf /= -1 )  THEN
2119       CALL MPI_TYPE_FREE( ft_surf, ierr )    ! if set, free last surface filetype
2120    ENDIF
2121#else
2122    all_nr_val(myid) = lo_nr_val(myid)
2123#endif
2124    nr_val = lo_nr_val(myid)
2125
2126    total_number_of_surface_values = 0
2127    DO  i = 0, numprocs-1
2128       IF ( i == myid )  THEN
2129          glo_start = total_number_of_surface_values + 1
2130       ENDIF
2131       total_number_of_surface_values = total_number_of_surface_values + all_nr_val(i)
2132    ENDDO
2133
2134!
2135!-- Actions during read
2136    IF ( rd_flag )  THEN
2137       IF ( .NOT. ALLOCATED( m_start_index )  )  ALLOCATE( m_start_index(nys:nyn,nxl:nxr)  )
2138       IF ( .NOT. ALLOCATED( m_end_index )    )  ALLOCATE( m_end_index(nys:nyn,nxl:nxr)    )
2139       IF ( .NOT. ALLOCATED( m_global_start ) )  ALLOCATE( m_global_start(nys:nyn,nxl:nxr) )
2140!
2141!--    Save arrays for later reading
2142       m_start_index  = start_index
2143       m_end_index    = end_index
2144       m_global_start = global_start
2145       nr_val = MAXVAL( end_index )
2146
2147#if defined( __parallel )
2148       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
2149#endif
2150    ENDIF
2151
2152!
2153!-- Actions during write
2154    IF ( wr_flag )  THEN
2155!
2156!--    Create surface filetype
2157       ft_surf      = -1
2158       global_start = start_index + glo_start - 1
2159
2160       WHERE ( end_index < start_index )
2161          global_start = -1
2162       ENDWHERE
2163
2164!
2165!--    Check, if surface data exist on this PE
2166       data_to_write = .TRUE.
2167       IF ( total_number_of_surface_values == 0 )  THEN
2168          data_to_write = .FALSE.
2169          RETURN
2170       ENDIF
2171
2172       all_pes_write = ( MINVAL( all_nr_val ) > 0 )
2173
2174       IF ( all_pes_write )  THEN
2175          dims1(1)  = total_number_of_surface_values
2176          lize1(1)  = nr_val
2177          start1(1) = glo_start-1
2178
2179#if defined( __parallel )
2180          IF ( total_number_of_surface_values > 0 )  THEN
2181              CALL MPI_TYPE_CREATE_SUBARRAY( 1, dims1, lize1, start1, MPI_ORDER_FORTRAN, MPI_REAL, ft_surf, ierr )
2182              CALL MPI_TYPE_COMMIT( ft_surf, ierr )
2183          ENDIF
2184#endif
2185       ENDIF
2186    ENDIF
2187
2188 END SUBROUTINE rd_mpi_io_surface_filetypes
2189
2190
2191
2192!--------------------------------------------------------------------------------------------------!
2193! Description:
2194! ------------
2195!> This subroutine creates file types to access 2d-/3d-REAL arrays and 2d-INTEGER arrays
2196!> distributed in blocks among processes to a single file that contains the global arrays.
2197!--------------------------------------------------------------------------------------------------!
2198  SUBROUTINE rs_mpi_io_create_filetypes
2199
2200    IMPLICIT NONE
2201
2202    INTEGER, DIMENSION(2) ::  dims2
2203    INTEGER, DIMENSION(2) ::  lize2
2204    INTEGER, DIMENSION(2) ::  start2
2205
2206    INTEGER, DIMENSION(3) ::  dims3
2207    INTEGER, DIMENSION(3) ::  lize3
2208    INTEGER, DIMENSION(3) ::  start3
2209
2210
2211!
2212!-- Decision, if storage with or without ghost layers.
2213!-- Please note that the indexing of the global array always starts at 0, even in Fortran.
2214!-- Therefore the PE local indices have to be shifted by nbgp in the case with ghost layers.
2215    IF ( include_total_domain_boundaries )  THEN
2216
2217       lb%nxl = nxl + nbgp
2218       lb%nxr = nxr + nbgp
2219       lb%nys = nys + nbgp
2220       lb%nyn = nyn + nbgp
2221       lb%nnx = nnx
2222       lb%nny = nny
2223       lb%nx  = nx + 2 * nbgp
2224       lb%ny  = ny + 2 * nbgp
2225       IF ( myidx == 0 )  THEN
2226          lb%nxl = lb%nxl - nbgp
2227          lb%nnx = lb%nnx + nbgp
2228       ENDIF
2229       IF ( myidx == npex-1  .OR.  npex == -1 )  THEN   ! npex == 1 if -D__parallel not set
2230          lb%nxr = lb%nxr + nbgp
2231          lb%nnx = lb%nnx + nbgp
2232       ENDIF
2233       IF ( myidy == 0 )  THEN
2234          lb%nys = lb%nys - nbgp
2235          lb%nny = lb%nny + nbgp
2236       ENDIF
2237       IF ( myidy == npey-1  .OR.  npey == -1 )  THEN   ! npey == 1 if -D__parallel not set
2238          lb%nyn = lb%nyn + nbgp
2239          lb%nny = lb%nny + nbgp
2240       ENDIF
2241
2242    ELSE
2243
2244       lb%nxl = nxl
2245       lb%nxr = nxr
2246       lb%nys = nys
2247       lb%nyn = nyn
2248       lb%nnx = nnx
2249       lb%nny = nny
2250       lb%nx  = nx
2251       lb%ny  = ny
2252
2253    ENDIF
2254
2255!
2256!-- Create filetype for 2d-REAL array with ghost layers around the total domain
2257    dims2(1)  = lb%nx + 1
2258    dims2(2)  = lb%ny + 1
2259
2260    lize2(1)  = lb%nnx
2261    lize2(2)  = lb%nny
2262
2263    start2(1) = lb%nxl
2264    start2(2) = lb%nys
2265
2266#if defined( __parallel )
2267    CALL MPI_TYPE_CREATE_SUBARRAY( 2, dims2, lize2, start2, MPI_ORDER_FORTRAN, MPI_REAL, ft_2d, ierr )
2268    CALL MPI_TYPE_COMMIT( ft_2d, ierr )
2269#endif
2270!
2271!-- Create filetype for 2d-INTEGER array without ghost layers around the total domain
2272    dims2(1)  = nx + 1
2273    dims2(2)  = ny + 1
2274
2275    lize2(1)  = nnx
2276    lize2(2)  = nny
2277
2278    start2(1) = nxl
2279    start2(2) = nys
2280
2281#if defined( __parallel )
2282    CALL MPI_TYPE_CREATE_SUBARRAY( 2, dims2, lize2, start2, MPI_ORDER_FORTRAN, MPI_INTEGER, ft_2di_nb, ierr )
2283    CALL MPI_TYPE_COMMIT( ft_2di_nb, ierr )
2284#endif
2285!
2286!-- Create filetype for 3d-REAL array
2287    dims3(1)  = nz + 2
2288    dims3(2)  = lb%nx + 1
2289    dims3(3)  = lb%ny + 1
2290
2291    lize3(1)  = dims3(1)
2292    lize3(2)  = lb%nnx
2293    lize3(3)  = lb%nny
2294
2295    start3(1) = nzb
2296    start3(2) = lb%nxl
2297    start3(3) = lb%nys
2298
2299#if defined( __parallel )
2300    CALL MPI_TYPE_CREATE_SUBARRAY( 3, dims3, lize3, start3, MPI_ORDER_FORTRAN, MPI_REAL, ft_3d, ierr )
2301    CALL MPI_TYPE_COMMIT( ft_3d, ierr )
2302#endif
2303
2304 END SUBROUTINE rs_mpi_io_create_filetypes
2305
2306
2307
2308!--------------------------------------------------------------------------------------------------!
2309! Description:
2310! ------------
2311!> This subroutine creates file types to access 3d-soil arrays
2312!> distributed in blocks among processes to a single file that contains the global arrays.
[4496]2313!> It is not required for the serial mode.
[4495]2314!--------------------------------------------------------------------------------------------------!
2315#if defined( __parallel )
2316 SUBROUTINE rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
2317
2318    IMPLICIT NONE
2319
2320    INTEGER, INTENT(IN)  ::  nzb_soil
2321    INTEGER, INTENT(IN)  ::  nzt_soil
2322
2323    INTEGER, DIMENSION(3) ::  dims3
2324    INTEGER, DIMENSION(3) ::  lize3
2325    INTEGER, DIMENSION(3) ::  start3
2326
2327
2328!
2329!-- Create filetype for 3d-soil array
2330    dims3(1)  = nzt_soil - nzb_soil + 1
2331    dims3(2)  = lb%nx + 1
2332    dims3(3)  = lb%ny + 1
2333
2334    lize3(1)  = dims3(1)
2335    lize3(2)  = lb%nnx
2336    lize3(3)  = lb%nny
2337
2338    start3(1) = nzb_soil
2339    start3(2) = lb%nxl
2340    start3(3) = lb%nys
2341
[4496]2342    CALL MPI_TYPE_CREATE_SUBARRAY( 3, dims3, lize3, start3, MPI_ORDER_FORTRAN, MPI_REAL,           &
2343                                   ft_3dsoil, ierr )
[4495]2344    CALL MPI_TYPE_COMMIT( ft_3dsoil, ierr )
2345
2346 END SUBROUTINE rd_mpi_io_create_filetypes_3dsoil
2347#endif
2348
2349
2350
2351!--------------------------------------------------------------------------------------------------!
2352! Description:
2353! ------------
2354!> Free all file types that have been created for MPI-IO.
2355!--------------------------------------------------------------------------------------------------!
2356 SUBROUTINE rs_mpi_io_free_filetypes
2357
2358    IMPLICIT NONE
2359
2360
2361#if defined( __parallel )
2362    IF ( filetypes_created )  THEN
2363       CALL MPI_TYPE_FREE( ft_2d, ierr )
2364       CALL MPI_TYPE_FREE( ft_2di_nb, ierr )
2365       CALL MPI_TYPE_FREE( ft_3d, ierr )
2366    ENDIF
2367!
2368!-- Free last surface filetype
2369    IF ( ft_surf /= -1 )  THEN 
2370       CALL MPI_TYPE_FREE( ft_surf, ierr )
2371    ENDIF
2372#endif
2373
2374 END SUBROUTINE rs_mpi_io_free_filetypes
2375
2376
2377
2378!--------------------------------------------------------------------------------------------------!
2379! Description:
2380! ------------
2381!> Print the restart data file header (MPI-IO format) for debugging.
2382!--------------------------------------------------------------------------------------------------!
2383 SUBROUTINE rs_mpi_io_print_header
2384
2385    IMPLICIT NONE
2386
2387    INTEGER(iwp) ::  i
2388
2389
2390    IF ( debug_level >= 1 )  THEN
2391 
2392       WRITE (9,*)  ' '
2393       WRITE (9,*)  ' CHARACTER header values ', tgh%nr_char
2394       WRITE (9,*)  ' '
2395       DO  i = 1, tgh%nr_char
2396          WRITE(9,*)  text_lines(i)(1:80)
2397       ENDDO
2398
2399       WRITE (9,*)  ' '
2400       WRITE (9,*) ' INTEGER header values ', tgh%nr_int
2401       WRITE (9,*)  ' '
2402       DO  i = 1, tgh%nr_int
2403          WRITE(9,*)  'INTERGER value:   ', int_names(i), ' ', int_values(i)
2404       ENDDO
2405
2406       WRITE (9,*)  ' '
2407       WRITE (9,*)  ' REAL header values ', tgh%nr_real
2408       WRITE (9,*)  ' '
2409       DO  i = 1, tgh%nr_real
2410          WRITE(9,*) 'REAL     value:   ', real_names(i), ' ', real_values(i)
2411       ENDDO
2412
2413       WRITE (9,*)  ' '
2414       WRITE (9,*)  ' Header entries with Offset ',tgh%nr_arrays
2415       WRITE (9,*)  '                Name                                  Offset '
2416       DO  i = 1, tgh%nr_arrays
2417          WRITE(9,'(a,1x,a30,1x,i16)') 'Header entiy:   ', array_names(i), array_offset(i)
2418       ENDDO
2419       WRITE (9,*)  ' '
2420    ENDIF
2421
2422    print_header_now = .FALSE.
2423
2424 END SUBROUTINE rs_mpi_io_print_header
2425
2426
2427
2428!--------------------------------------------------------------------------------------------------!
2429! Description:
2430! ------------
2431!> Print error messages for reading/writing restart data with MPI-IO
2432!--------------------------------------------------------------------------------------------------!
2433 SUBROUTINE rs_mpi_io_error( error_number )
2434
2435    IMPLICIT NONE
2436
2437    INTEGER, INTENT(IN) ::  error_number
2438
2439    IF ( myid == 0)  THEN
2440
2441       SELECT CASE (error_number)
2442 
2443          CASE ( 1 )
2444             WRITE(6,*)  'illegal action while opening restart File'
2445          CASE ( 2 )
2446             WRITE(6,*)  'data array not found in restart File'
2447          CASE ( 3 )
2448             WRITE(6,*)  'INTEGER or REAL value not found in restart File'
2449          CASE ( 4 )
2450             WRITE(6,*)  'Arrays only array with nbgp Halos or without halos legal'
2451          CASE ( 5 )
2452             WRITE(6,*)  'outer boundary in model and outer boundary in restart file do not match'
2453          CASE ( 6 )
2454             WRITE(6,*)  'posix IO: ERROR Opening Restart File'
2455          CASE DEFAULT
2456             WRITE(6,*)  'rs_mpi_io_error: illegal error number: ',error_number
2457
2458       END SELECT
2459
2460    ENDIF
2461#if defined( __parallel )
2462    CALL MPI_BARRIER( comm2d, ierr )
2463    CALL MPI_ABORT( comm2d, 1, ierr )
2464#else
2465    CALL ABORT
2466#endif
2467
2468 END SUBROUTINE rs_mpi_io_error
2469
2470
2471
2472!--------------------------------------------------------------------------------------------------!
2473! Description:
2474! ------------
2475!> Check if big/little endian data format is used.
2476!> An int*4 pointer is set to a int*8 variable, the int*8 is set to 1, and then it is checked, if
2477!> the first 4 bytes of the pointer are equal 1 (little endian) or not.
2478!--------------------------------------------------------------------------------------------------!
2479 SUBROUTINE rs_mpi_io_check_endian( i_endian )
2480
2481    IMPLICIT NONE
2482
2483    INTEGER, INTENT(out)                   ::  i_endian
2484    INTEGER(KIND=8), TARGET                ::  int8
2485
2486    INTEGER, DIMENSION(1)                  ::  bufshape
2487    INTEGER(KIND=4), POINTER, DIMENSION(:) ::  int4
2488
2489    TYPE(C_PTR)                            ::  ptr
2490
2491
2492    ptr = C_LOC( int8 )
2493    bufshape(1) = 2
2494    CALL C_F_POINTER( ptr, int4, bufshape )
2495
2496    int8 = 1
2497
2498    IF ( int4(1) == 1 )  THEN
2499       i_endian = 1    ! little endian
2500    ELSE
2501       i_endian = 2    ! big endian
2502    ENDIF
2503
2504 END SUBROUTINE rs_mpi_io_check_endian
2505
2506 END MODULE restart_data_mpi_io_mod
Note: See TracBrowser for help on using the repository browser.