Home

Context Navigation

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

Last change on this file since 4534 was 4534, checked in by raasch, 4 years ago
I/O on reduced number of cores added (using shared memory MPI)
Property svn:keywords set to `Id`
File size: 93.2 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |