Home

Context Navigation

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

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