Home

Context Navigation

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

Last change on this file since 4618 was 4617, checked in by raasch, 4 years ago
cyclic fill mode implemented for MPI-IO, check, if boundary conditions in the prerun are both set to cyclic
Property svn:keywords set to `Id`
File size: 125.4 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 4617 2020-07-22 09:48:50Z raasch $
27	! Cyclic fill mode implemented
28	!
29	! 4598 2020-07-10 10:13:23Z suehring
30	! Bugfix in treatment of 3D soil arrays
31	!
32	! 4591 2020-07-06 15:56:08Z raasch
33	! File re-formatted to follow the PALM coding standard
34	!
35	! 4539 2020-05-18 14:05:17Z raasch
36	! Checks added, if index limits in header are exceeded
37	! Bugfix in rrd_mpi_io_int_2d
38	!
39	! 4536 2020-05-17 17:24:13Z raasch
40	! Messages and debug output converted to PALM routines
41	!
42	! 4534 2020-05-14 18:35:22Z raasch
43	! I/O on reduced number of cores added (using shared memory MPI)
44	!
45	! 4500 2020-04-17 10:12:45Z suehring
46	! Fix too long lines
47	!
48	! 4498 2020-04-15 14:26:31Z raasch
49	! Bugfix for creation of filetypes, argument removed from rd_mpi_io_open
50	!
51	! 4497 2020-04-15 10:20:51Z raasch
52	! Last bugfix deactivated because of compile problems
53	!
54	! 4496 2020-04-15 08:37:26Z raasch
55	! Problem with posix read arguments for surface data fixed
56	!
57	! 4495 2020-04-13 20:11:20Z raasch
58	! Initial version (K. Ketelsen), adjusted to PALM formatting standards (s. Raasch)
59	!
60	!
61	!
62	! Description:
63	! ------------
64	!> Routines for restart data handling using MPI-IO.
65	!--------------------------------------------------------------------------------------------------!
66	MODULE restart_data_mpi_io_mod
67
68	#if defined( __parallel )
69	#if defined( __mpifh )
70	INCLUDE "mpif.h"
71	#else
72	USE MPI
73	#endif
74	#else
75	USE posix_interface, &
76	ONLY: posix_close, &
77	posix_lseek, &
78	posix_open, &
79	posix_read, &
80	posix_write
81	#endif
82
83	USE, INTRINSIC :: ISO_C_BINDING
84
85	USE control_parameters, &
86	ONLY: debug_output, &
87	debug_string, &
88	include_total_domain_boundaries, &
89	message_string, &
90	restart_data_format_input, &
91	restart_data_format_output, &
92	restart_file_size
93
94	USE exchange_horiz_mod, &
95	ONLY: exchange_horiz, &
96	exchange_horiz_2d
97
98	USE indices, &
99	ONLY: nbgp, &
100	nnx, &
101	nny, &
102	nx, &
103	nxl, &
104	nxlg, &
105	nx_on_file, &
106	nxr, &
107	nxrg, &
108	ny, &
109	nyn, &
110	nyng, &
111	ny_on_file, &
112	nys, &
113	nysg, &
114	nz, &
115	nzb, &
116	nzt
117
118	USE kinds
119
120	USE pegrid, &
121	ONLY: comm1dx, &
122	comm1dy, &
123	comm2d, &
124	communicator_configurations, &
125	myid, &
126	myidx, &
127	myidy, &
128	npex, &
129	npey, &
130	numprocs, &
131	pdims
132
133	USE shared_memory_io_mod, &
134	ONLY: domain_decomposition_grid_features, &
135	sm_class
136
137
138	IMPLICIT NONE
139
140	CHARACTER(LEN=128) :: io_file_name !> internal variable to communicate filename between different subroutines
141
142	#if defined( __parallel )
143	INTEGER(iwp) :: ierr !< error status of MPI-calls
144	INTEGER(iwp), PARAMETER :: rd_offset_kind = MPI_OFFSET_KIND !< Adress or Offset kind
145	INTEGER(iwp), PARAMETER :: rd_status_size = MPI_STATUS_SIZE !<
146	#else
147	INTEGER(iwp), PARAMETER :: rd_offset_kind = C_SIZE_T !<
148	INTEGER(iwp), PARAMETER :: rd_status_size = 1 !< Not required in sequential mode
149	#endif
150
151	INTEGER(iwp) :: debug_level = 1 !< TODO: replace with standard debug output steering
152
153	INTEGER(iwp) :: comm_io !< Communicator for MPI-IO
154	INTEGER(iwp) :: fh !< MPI-IO file handle
155	#if defined( __parallel )
156	INTEGER(iwp) :: fhs = -1 !< MPI-IO file handle to open file with comm2d always
157	#endif
158	INTEGER(iwp) :: ft_surf = -1 !< MPI filetype surface data
159	#if defined( __parallel )
160	INTEGER(iwp) :: ft_2di_nb !< MPI filetype 2D array INTEGER no outer boundary
161	INTEGER(iwp) :: ft_2d !< MPI filetype 2D array REAL with outer boundaries
162	INTEGER(iwp) :: ft_3d !< MPI filetype 3D array REAL with outer boundaries
163	INTEGER(iwp) :: ft_3dsoil !< MPI filetype for 3d-soil array
164	#endif
165	INTEGER(iwp) :: glo_start !< global start index on this PE
166	#if defined( __parallel )
167	INTEGER(iwp) :: local_start !<
168	#endif
169	INTEGER(iwp) :: nr_iope !<
170	INTEGER(iwp) :: nr_val !< local number of values in x and y direction
171	#if defined( __parallel )
172	INTEGER(iwp) :: win_2di !<
173	INTEGER(iwp) :: win_2dr !<
174	INTEGER(iwp) :: win_3dr !<
175	INTEGER(iwp) :: win_3ds !<
176	INTEGER(iwp) :: win_surf = -1 !<
177	#endif
178	INTEGER(iwp) :: total_number_of_surface_values !< total number of values for one variable
179
180	INTEGER(KIND=rd_offset_kind) :: array_position !<
181	INTEGER(KIND=rd_offset_kind) :: header_position !<
182
183	INTEGER(iwp), DIMENSION(:,:), POINTER, CONTIGUOUS :: array_2di !<
184
185	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: m_end_index !<
186	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: m_global_start !<
187	INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: m_start_index !<
188
189
190	LOGICAL :: all_pes_write !< all PEs have data to write
191	LOGICAL :: filetypes_created !<
192	LOGICAL :: io_on_limited_cores_per_node !< switch to shared memory MPI-IO
193	LOGICAL :: rd_flag !< file is opened for read
194	LOGICAL :: wr_flag !< file is opened for write
195
196	#if defined( __parallel )
197	REAL(wp), DIMENSION(:), POINTER, CONTIGUOUS :: array_1d !<
198	#endif
199	REAL(wp), DIMENSION(:,:), POINTER, CONTIGUOUS :: array_2d !<
200	REAL(wp), DIMENSION(:,:,:), POINTER, CONTIGUOUS :: array_3d !<
201	REAL(wp), DIMENSION(:,:,:), POINTER, CONTIGUOUS :: array_3d_soil !<
202
203	!
204	!-- Variable to store the grid features (index bounds) of the temporary arrays that are used
205	!-- to read and write the restart data. They differ depending on if the outer boundary of the
206	!-- total domain is contained in the restart data or not. iog stands for IO-grid.
207	TYPE(domain_decomposition_grid_features) :: iog !<
208
209	!
210	!-- General Header (first 32 byte in restart file)
211	TYPE general_header
212	INTEGER(iwp) :: endian !< little endian (1) or big endian (2) internal format
213	INTEGER(iwp) :: i_outer_bound !< if 1, outer boundaries are stored in restart file
214	INTEGER(iwp) :: nr_arrays !< number of arrays in restart files
215	INTEGER(iwp) :: nr_char !< number of Text strings entries in header
216	INTEGER(iwp) :: nr_int !< number of INTEGER entries in header
217	INTEGER(iwp) :: nr_real !< number of REAL entries in header
218	INTEGER(iwp) :: total_nx !< total number of points in x-direction
219	INTEGER(iwp) :: total_ny !< total number of points in y-direction
220	END TYPE general_header
221
222	TYPE(general_header), TARGET :: tgh !<
223
224	TYPE(sm_class) :: sm_io !<
225
226	!
227	!-- Declaration of varibales for file header section
228	INTEGER(KIND=rd_offset_kind) :: header_int_index
229	INTEGER, PARAMETER :: max_nr_int=256
230	CHARACTER(LEN=32), DIMENSION(max_nr_int) :: int_names
231	INTEGER(KIND=iwp), DIMENSION(max_nr_int) :: int_values
232
233	INTEGER(KIND=rd_offset_kind) :: header_char_index
234	INTEGER, PARAMETER :: max_nr_char=128
235	CHARACTER(LEN=128), DIMENSION(max_nr_char) :: text_lines
236
237	INTEGER(KIND=rd_offset_kind) :: header_real_index
238	INTEGER, PARAMETER :: max_nr_real=256
239	CHARACTER(LEN=32), DIMENSION(max_nr_real) :: real_names
240	REAL(KIND=wp), DIMENSION(max_nr_real) :: real_values
241
242	INTEGER(KIND=rd_offset_kind) :: header_array_index
243	INTEGER, PARAMETER :: max_nr_arrays=600
244	CHARACTER(LEN=32), DIMENSION(max_nr_arrays) :: array_names
245	INTEGER(KIND=rd_offset_kind), DIMENSION(max_nr_arrays) :: array_offset
246
247	!
248	!-- Variables to handle the cyclic fill initialization mode
249	INTEGER :: comm_cyclic_fill !< communicator for cyclic fill PEs
250	INTEGER :: rmawin_2di !< RMA window 2d INTEGER
251	INTEGER :: rmawin_2d !< RMA window 2d REAL
252	INTEGER :: rmawin_3d !< RMA window 3d
253
254	INTEGER(iwp), ALLOCATABLE, DIMENSION(:,:) :: remote_pe
255	INTEGER(iwp), ALLOCATABLE, DIMENSION(:,:) :: remote_pe_s
256	INTEGER(iwp), ALLOCATABLE, DIMENSION(:,:) :: rma_offset
257	INTEGER(iwp), ALLOCATABLE, DIMENSION(:,:) :: rma_offset_s
258	INTEGER(iwp), ALLOCATABLE, DIMENSION(:,:) :: rmabuf_2di
259
260	LOGICAL :: cyclic_fill_mode !< arrays are filled cyclically with data from prerun
261	LOGICAL :: pe_active_for_read = .TRUE. !< this PE is active for reading data from prerun or
262	!< restart run. For restarts all PEs are active.
263
264	REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: rmabuf_2d
265	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: rmabuf_3d
266
267	TYPE(domain_decomposition_grid_features) :: mainrun_grid !< grid variables for the main run
268	TYPE(domain_decomposition_grid_features) :: prerun_grid !< grid variables for the prerun
269
270
271	SAVE
272
273	PRIVATE
274
275	PUBLIC restart_file_size, total_number_of_surface_values
276
277	!
278	!-- PALM interfaces
279	INTERFACE rd_mpi_io_check_array
280	MODULE PROCEDURE rd_mpi_io_check_array
281	END INTERFACE rd_mpi_io_check_array
282
283	INTERFACE rd_mpi_io_close
284	MODULE PROCEDURE rd_mpi_io_close
285	END INTERFACE rd_mpi_io_close
286
287	INTERFACE rd_mpi_io_open
288	MODULE PROCEDURE rd_mpi_io_open
289	END INTERFACE rd_mpi_io_open
290
291	INTERFACE rrd_mpi_io
292	MODULE PROCEDURE rrd_mpi_io_char
293	MODULE PROCEDURE rrd_mpi_io_int
294	MODULE PROCEDURE rrd_mpi_io_int_2d
295	MODULE PROCEDURE rrd_mpi_io_logical
296	MODULE PROCEDURE rrd_mpi_io_real
297	MODULE PROCEDURE rrd_mpi_io_real_2d
298	MODULE PROCEDURE rrd_mpi_io_real_3d
299	MODULE PROCEDURE rrd_mpi_io_real_3d_soil
300	END INTERFACE rrd_mpi_io
301
302	INTERFACE rrd_mpi_io_global_array
303	MODULE PROCEDURE rrd_mpi_io_global_array_int_1d
304	MODULE PROCEDURE rrd_mpi_io_global_array_real_1d
305	MODULE PROCEDURE rrd_mpi_io_global_array_real_2d
306	MODULE PROCEDURE rrd_mpi_io_global_array_real_3d
307	MODULE PROCEDURE rrd_mpi_io_global_array_real_4d
308	END INTERFACE rrd_mpi_io_global_array
309
310	INTERFACE rrd_mpi_io_surface
311	MODULE PROCEDURE rrd_mpi_io_surface
312	MODULE PROCEDURE rrd_mpi_io_surface_2d
313	END INTERFACE rrd_mpi_io_surface
314
315	INTERFACE rd_mpi_io_surface_filetypes
316	MODULE PROCEDURE rd_mpi_io_surface_filetypes
317	END INTERFACE rd_mpi_io_surface_filetypes
318
319	INTERFACE wrd_mpi_io
320	MODULE PROCEDURE wrd_mpi_io_char
321	MODULE PROCEDURE wrd_mpi_io_int
322	MODULE PROCEDURE wrd_mpi_io_int_2d
323	MODULE PROCEDURE wrd_mpi_io_logical
324	MODULE PROCEDURE wrd_mpi_io_real
325	MODULE PROCEDURE wrd_mpi_io_real_2d
326	MODULE PROCEDURE wrd_mpi_io_real_3d
327	MODULE PROCEDURE wrd_mpi_io_real_3d_soil
328	END INTERFACE wrd_mpi_io
329
330	INTERFACE wrd_mpi_io_global_array
331	MODULE PROCEDURE wrd_mpi_io_global_array_int_1d
332	MODULE PROCEDURE wrd_mpi_io_global_array_real_1d
333	MODULE PROCEDURE wrd_mpi_io_global_array_real_2d
334	MODULE PROCEDURE wrd_mpi_io_global_array_real_3d
335	MODULE PROCEDURE wrd_mpi_io_global_array_real_4d
336	END INTERFACE wrd_mpi_io_global_array
337
338	INTERFACE wrd_mpi_io_surface
339	MODULE PROCEDURE wrd_mpi_io_surface
340	MODULE PROCEDURE wrd_mpi_io_surface_2d
341	END INTERFACE wrd_mpi_io_surface
342
343	PUBLIC rd_mpi_io_check_array, &
344	rd_mpi_io_close, &
345	rd_mpi_io_open, &
346	rrd_mpi_io, &
347	rrd_mpi_io_global_array, &
348	rrd_mpi_io_surface, &
349	rd_mpi_io_surface_filetypes, &
350	wrd_mpi_io, &
351	wrd_mpi_io_global_array, &
352	wrd_mpi_io_surface
353
354
355	CONTAINS
356
357
358	!--------------------------------------------------------------------------------------------------!
359	! Description:
360	! ------------
361	!> Open restart file for read or write with MPI-IO
362	!--------------------------------------------------------------------------------------------------!
363	SUBROUTINE rd_mpi_io_open( action, file_name, open_for_global_io_only )
364
365	IMPLICIT NONE
366
367	CHARACTER(LEN=*), INTENT(IN) :: action !<
368	CHARACTER(LEN=*), INTENT(IN) :: file_name !<
369
370	INTEGER(iwp) :: i !<
371	INTEGER(iwp) :: gh_size !<
372
373	INTEGER(KIND=rd_offset_kind) :: offset !<
374
375	#if defined( __parallel )
376	INTEGER, DIMENSION(rd_status_size) :: status !<
377	#endif
378
379	LOGICAL, INTENT(IN), OPTIONAL :: open_for_global_io_only !<
380	LOGICAL :: set_filetype !<
381
382	#if ! defined( __parallel )
383	TYPE(C_PTR) :: buf_ptr !<
384	#endif
385
386	! write(9,*) 'Here is rd_mpi_io_open',nx,nx_on_file,ny,ny_on_file,TRIM(action) !kk may become Debug Output
387
388	offset = 0
389	io_on_limited_cores_per_node = .FALSE.
390
391	rd_flag = ( TRIM( action ) == 'READ' .OR. TRIM( action ) == 'read' )
392	wr_flag = ( TRIM( action ) == 'WRITE' .OR. TRIM( action ) == 'write' )
393
394	IF ( .NOT. ( rd_flag .OR. wr_flag ) ) THEN
395	message_string = 'illegal action "' // TRIM( action ) // '" for opening restart files'
396	CALL message( 'restart_data_mpi_io_mod', 'PA0720', 1, 2, 0, 6, 0 )
397	ENDIF
398	!
399	!-- Store name of I/O file to communicate it internally within this module.
400	io_file_name = file_name
401	!
402	!-- Setup for IO on a limited number of threads per node (using shared memory MPI)
403	IF ( rd_flag ) THEN
404	set_filetype = .TRUE.
405	IF ( TRIM( restart_data_format_input ) == 'mpi_shared_memory' ) THEN
406	io_on_limited_cores_per_node = .TRUE.
407	ENDIF
408	ENDIF
409
410	IF ( TRIM( restart_data_format_output ) == 'mpi_shared_memory' .AND. wr_flag ) THEN
411	io_on_limited_cores_per_node = .TRUE.
412	ENDIF
413	!
414	!-- Shared memory MPI is not used for reading of global data
415	IF ( PRESENT( open_for_global_io_only ) .AND. rd_flag ) THEN
416	IF ( open_for_global_io_only ) THEN
417	io_on_limited_cores_per_node = .FALSE.
418	set_filetype = .FALSE.
419	ENDIF
420	ENDIF
421
422	!
423	!-- Determine, if prerun data shall be read and mapped cyclically to the mainrun arrays.
424	!-- In cyclic fill mode only a subset of the PEs will read.
425	cyclic_fill_mode = .FALSE.
426	pe_active_for_read = .TRUE.
427
428	IF ( rd_flag .AND. .NOT. PRESENT( open_for_global_io_only ) .AND. &
429	nx_on_file < nx .AND. ny_on_file < ny ) THEN
430	cyclic_fill_mode = .TRUE.
431	CALL setup_cyclic_fill
432	!
433	!-- Shared memory IO on limited cores is not allowed for cyclic fill mode
434	CALL sm_io%sm_init_comm( .FALSE. ) !
435	ELSE
436	CALL sm_io%sm_init_comm( io_on_limited_cores_per_node )
437	ENDIF
438
439	!
440	!-- TODO: add a more detailed meaningful comment about what is happening here
441	!-- activate model grid
442	IF( cyclic_fill_mode .AND. .NOT. pe_active_for_read ) THEN
443	CALL mainrun_grid%activate_grid_from_this_class()
444	RETURN
445	ENDIF
446
447
448	!
449	!-- Set communicator to be used. If all cores are doing I/O, comm2d is used as usual.
450	IF( sm_io%is_sm_active() ) THEN
451	comm_io = sm_io%comm_io
452	ELSEIF ( cyclic_fill_mode ) THEN
453	comm_io = comm_cyclic_fill
454	ELSE
455	comm_io = comm2d
456	ENDIF
457
458	!
459	!-- Create subarrays and file types
460	filetypes_created = .FALSE.
461
462	!
463	!-- In case of read it is not known yet if data include total domain. Filetypes will be created
464	!-- further below.
465	IF ( wr_flag ) THEN
466	CALL rd_mpi_io_create_filetypes
467	filetypes_created = .TRUE.
468	ENDIF
469
470	!
471	!-- Open file for MPI-IO
472	#if defined( __parallel )
473	IF ( sm_io%iam_io_pe ) THEN
474
475	IF ( rd_flag ) THEN
476
477	IF ( debug_output ) THEN
478	WRITE( debug_string, * ) 'open joint restart file "' // TRIM( io_file_name ) // &
479	'" for read with MPI-IO'
480	CALL debug_message( debug_string, 'start' )
481	ENDIF
482
483	CALL MPI_FILE_OPEN( comm_io, TRIM( io_file_name ), MPI_MODE_RDONLY, MPI_INFO_NULL, fh, &
484	ierr )
485
486	IF ( ierr /= 0 ) THEN
487	message_string = 'error opening restart file "' // TRIM( io_file_name ) // &
488	'" for reading with MPI-IO'
489	CALL message( 'rrd_mpi_io_open', 'PA0727', 3, 2, 0, 6, 0 )
490	ENDIF
491
492	IF ( debug_output ) THEN
493	WRITE( debug_string, * ) 'open joint restart file "' // TRIM( io_file_name ) // &
494	'" for read with MPI-IO'
495	CALL debug_message( debug_string, 'end' )
496	ENDIF
497
498	ELSEIF ( wr_flag ) THEN
499
500	IF ( debug_output ) THEN
501	WRITE( debug_string, * ) 'open joint restart file "' // TRIM( io_file_name ) // &
502	'" for write with MPI-IO'
503	CALL debug_message( debug_string, 'start' )
504	ENDIF
505
506	CALL MPI_FILE_OPEN( comm_io, TRIM( io_file_name ), MPI_MODE_CREATE+MPI_MODE_WRONLY, &
507	MPI_INFO_NULL, fh, ierr )
508
509	IF ( ierr /= 0 ) THEN
510	message_string = 'error opening restart file "' // TRIM( io_file_name ) // &
511	'" for writing with MPI-IO'
512	CALL message( 'rrd_mpi_io_open', 'PA0728', 3, 2, 0, 6, 0 )
513	ENDIF
514
515	IF ( debug_output ) THEN
516	WRITE( debug_string, * ) 'open joint restart file "' // TRIM( io_file_name ) // &
517	'" for write with MPI-IO'
518	CALL debug_message( debug_string, 'end' )
519	ENDIF
520
521	ENDIF
522
523	ENDIF
524	#else
525	IF ( rd_flag ) THEN
526
527	IF ( debug_output ) THEN
528	WRITE( debug_string, * ) 'open restart file "' // TRIM( io_file_name ) // &
529	'" for read in serial mode (posix)'
530	CALL debug_message( debug_string, 'start' )
531	ENDIF
532
533	fh = posix_open( TRIM( io_file_name ), .TRUE. )
534
535	IF ( debug_output ) THEN
536	WRITE( debug_string, * ) 'open restart file "' // TRIM( io_file_name ) // &
537	'" for read in serial mode (posix)'
538	CALL debug_message( debug_string, 'end' )
539	ENDIF
540
541	ELSEIF ( wr_flag ) THEN
542
543	IF ( debug_output ) THEN
544	WRITE( debug_string, * ) 'open restart file "' // TRIM( io_file_name ) // &
545	'" for write in serial mode (posix)'
546	CALL debug_message( debug_string, 'start' )
547	ENDIF
548
549	fh = posix_open( TRIM( io_file_name ), .FALSE. )
550
551	IF ( debug_output ) THEN
552	WRITE( debug_string, * ) 'open restart file "' // TRIM( io_file_name ) // &
553	'" for write in serial mode (posix)'
554	CALL debug_message( debug_string, 'end' )
555	ENDIF
556
557	ENDIF
558
559	IF ( fh < 0 ) THEN
560	message_string = 'error opening restart file for posix I/O'
561	CALL message( 'restart_data_mpi_io_mod', 'PA0721', 1, 2, 0, 6, 0 )
562	ENDIF
563	#endif
564
565	array_position = 65536 !> Start offset for writing 2-D and 3.D arrays at 64 k
566	header_position = 0
567
568	header_int_index = 1
569	header_char_index = 1
570	header_real_index = 1
571	header_array_index = 1
572
573	int_names = ' '
574	int_values = 0
575	text_lines = ' '
576	real_names = ' '
577	real_values = 0.0
578	array_names = ' '
579	array_offset = 0
580
581	int_names(1) = 'nx'
582	int_values(1) = nx
583	int_names(2) = 'ny'
584	int_values(2) = ny
585	int_names(3) = 'nz'
586	int_values(3) = nz
587	header_int_index = header_int_index+3
588
589	DO i = 1, max_nr_arrays
590	array_offset(i) = 0
591	array_names(i) = ' '
592	ENDDO
593
594	gh_size = STORAGE_SIZE( tgh ) / 8
595
596	IF ( rd_flag ) THEN
597	IF ( sm_io%iam_io_pe ) THEN
598	!
599	!-- File is open for reading
600	#if defined( __parallel )
601	!-- Set the default view
602	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
603	!
604	!-- Read the file header size
605	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
606	CALL MPI_FILE_READ( fh, tgh, gh_size, MPI_BYTE, status, ierr )
607	#else
608	CALL posix_lseek( fh, header_position )
609	buf_ptr = C_LOC( tgh )
610	CALL posix_read( fh, buf_ptr, gh_size )
611	#endif
612	ENDIF
613	#if defined( __parallel )
614	IF ( sm_io%is_sm_active() ) THEN
615	CALL MPI_BCAST( tgh, gh_size, MPI_BYTE, 0, sm_io%comm_shared, ierr )
616	ENDIF
617	#endif
618	header_position = header_position + gh_size
619
620	include_total_domain_boundaries = ( tgh%i_outer_bound == 1 )
621
622	!
623	!-- File types depend on boundaries of the total domain being included in data. This has been
624	!-- checked with the previous statement.
625	IF ( set_filetype ) THEN
626	CALL rd_mpi_io_create_filetypes
627	filetypes_created = .TRUE.
628	ENDIF
629
630	IF ( sm_io%iam_io_pe ) THEN
631	#if defined( __parallel )
632	!
633	!-- Read INTEGER values
634	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
635	CALL MPI_FILE_READ( fh, int_names, SIZE( int_names ) * 32, MPI_CHAR, status, ierr )
636	header_position = header_position + SIZE( int_names ) * 32
637
638	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
639	CALL MPI_FILE_READ (fh, int_values, SIZE( int_values ), MPI_INT, status, ierr )
640	header_position = header_position + SIZE( int_values ) * iwp
641	!
642	!-- Character entries
643	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
644	CALL MPI_FILE_READ( fh, text_lines, SIZE( text_lines ) * 128, MPI_CHAR, status, ierr )
645	header_position = header_position + SIZE ( text_lines ) * 128
646	!
647	!-- REAL values
648	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
649	CALL MPI_FILE_READ( fh, real_names, SIZE( real_names ) * 32, MPI_CHAR, status, ierr )
650	header_position = header_position + SIZE( real_names ) * 32
651
652	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
653	CALL MPI_FILE_READ( fh, real_values, SIZE( real_values ), MPI_REAL, status, ierr )
654	header_position = header_position + SIZE( real_values ) * wp
655	!
656	!-- 2d- and 3d-array headers
657	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
658	CALL MPI_FILE_READ( fh, array_names, SIZE( array_names ) * 32, MPI_CHAR, status, ierr )
659	header_position = header_position + SIZE( array_names ) * 32
660
661	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
662	CALL MPI_FILE_READ( fh, array_offset, SIZE( array_offset ) * MPI_OFFSET_KIND, MPI_BYTE, &
663	status,ierr ) ! there is no I*8 datatype in Fortran
664	header_position = header_position + SIZE( array_offset ) * rd_offset_kind
665	#else
666	CALL posix_lseek( fh, header_position )
667	CALL posix_read( fh, int_names )
668	header_position = header_position + SIZE( int_names ) * 32
669
670	CALL posix_lseek( fh, header_position )
671	CALL posix_read( fh, int_values, SIZE( int_values ) )
672	header_position = header_position + SIZE( int_values ) * iwp
673	!
674	!-- Character entries
675	CALL posix_lseek( fh, header_position )
676	CALL posix_read( fh, text_lines )
677	header_position = header_position + SIZE( text_lines ) * 128
678	!
679	!-- REAL values
680	CALL posix_lseek( fh, header_position )
681	CALL posix_read( fh, real_names )
682	header_position = header_position + SIZE( real_names ) * 32
683
684	CALL posix_lseek( fh, header_position )
685	CALL posix_read( fh, real_values, SIZE( real_values ) )
686	header_position = header_position + SIZE( real_values ) * wp
687	!
688	!-- 2d- and 3d-array headers
689	CALL posix_lseek( fh, header_position )
690	CALL posix_read( fh, array_names )
691	header_position = header_position + SIZE( array_names ) * 32
692
693	CALL posix_lseek( fh, header_position )
694	CALL posix_read( fh, array_offset, SIZE( array_offset ) ) ! there is no I*8 datatype in Fortran
695	header_position = header_position + SIZE( array_offset ) * rd_offset_kind
696	#endif
697	IF ( debug_output ) CALL rd_mpi_io_print_header
698
699	ENDIF
700
701	#if defined( __parallel )
702	!
703	!-- Broadcast header to all remaining cores that are not involved in I/O
704	IF ( sm_io%is_sm_active() ) THEN
705	!
706	!-- Not sure, that it is possible to broadcast CHARACTER array in one MPI_Bcast call
707	DO i = 1, SIZE( int_names )
708	CALL MPI_BCAST( int_names(i), 32, MPI_CHARACTER, 0, sm_io%comm_shared, ierr )
709	ENDDO
710	CALL MPI_BCAST( int_values, SIZE( int_values ), MPI_INTEGER, 0, sm_io%comm_shared, ierr )
711
712	DO i = 1, SIZE( text_lines )
713	CALL MPI_BCAST( text_lines(i), 128, MPI_CHARACTER, 0, sm_io%comm_shared, ierr )
714	ENDDO
715
716	DO i = 1, SIZE( real_names )
717	CALL MPI_BCAST( real_names(i), 32, MPI_CHARACTER, 0, sm_io%comm_shared, ierr )
718	ENDDO
719	CALL MPI_BCAST( real_values, SIZE( real_values ), MPI_REAL, 0, sm_io%comm_shared, ierr )
720
721	DO i = 1, SIZE( array_names )
722	CALL MPI_BCAST( array_names(i), 32, MPI_CHARACTER, 0, sm_io%comm_shared, ierr )
723	ENDDO
724	CALL MPI_BCAST( array_offset, SIZE( array_offset )*8, MPI_BYTE, 0, sm_io%comm_shared, &
725	ierr ) ! there is no I8 datatype in Fortran (array_offset is I8!)
726
727	CALL MPI_BCAST( header_position, rd_offset_kind, MPI_BYTE, 0, sm_io%comm_shared, ierr )
728
729	ENDIF
730	#endif
731
732
733	ENDIF
734
735	!
736	!-- TODO: describe in more detail what is happening here
737	!-- activate model grid
738	IF ( cyclic_fill_mode ) CALL mainrun_grid%activate_grid_from_this_class()
739
740	CONTAINS
741
742	SUBROUTINE setup_cyclic_fill
743
744	IMPLICIT NONE
745
746	INTEGER :: color !< used to set the IO PEs for MPI_COMM_SPLIT
747	INTEGER :: ierr !<
748	INTEGER(iwp) :: i !<
749	INTEGER(iwp) :: j !<
750	INTEGER(KIND=MPI_ADDRESS_KIND) :: winsize !< size of RMA window
751
752	!
753	!-- TODO: describe in more detail what is done here and why it is done
754	!-- save grid of main run
755	CALL mainrun_grid%save_grid_into_this_class()
756
757	ALLOCATE( remote_pe(0:nx_on_file,0:ny_on_file) )
758	ALLOCATE( remote_pe_s(0:nx_on_file,0:ny_on_file) )
759	ALLOCATE( rma_offset(0:nx_on_file,0:ny_on_file) )
760	ALLOCATE( rma_offset_s(0:nx_on_file,0:ny_on_file) )
761
762	remote_pe_s = 0
763	rma_offset_s = 0
764	!
765	!-- Determine, if gridpoints of the prerun are located on this thread.
766	!-- Set the (cyclic) prerun grid.
767	nxr = MIN( nxr, nx_on_file )
768	IF ( nxl > nx_on_file ) THEN
769	nxl = -99
770	nxr = -99
771	nnx = 0
772	ELSE
773	nnx =nxr-nxl+1
774	ENDIF
775
776	nyn = MIN( nyn, ny_on_file )
777	IF ( nys > ny_on_file ) THEN
778	nys = -99
779	nyn = -99
780	nny = 0
781	ELSE
782	nny = nyn-nys+1
783	ENDIF
784
785	nx = nx_on_file
786	ny = ny_on_file
787	!
788	!-- Determine, if this thread is doing IO
789	IF ( nnx > 0 .AND. nny > 0 ) THEN
790	color = 1
791	pe_active_for_read = .TRUE.
792	remote_pe_s(nxl:nxr,nys:nyn) = myid ! myid from comm2d
793	DO j = nys, nyn
794	DO i = nxl, nxr
795	rma_offset_s(i,j) = ( j-nys ) + ( i-nxl ) * nny
796	ENDDO
797	ENDDO
798	ELSE
799	color = MPI_UNDEFINED
800	pe_active_for_read = .FALSE.
801	ENDIF
802
803	#if defined( __parallel )
804	CALL MPI_ALLREDUCE( remote_pe_s, remote_pe, SIZE(remote_pe_s), MPI_INTEGER, MPI_SUM, &
805	comm2d, ierr )
806	CALL MPI_ALLREDUCE( rma_offset_s, rma_offset, SIZE(rma_offset_s), MPI_INTEGER, MPI_SUM, &
807	comm2d, ierr )
808	CALL MPI_COMM_SPLIT( comm2d, color, 0, comm_cyclic_fill, ierr )
809
810	IF ( pe_active_for_read ) THEN
811	CALL MPI_COMM_SIZE( comm_cyclic_fill, numprocs, ierr )
812	CALL MPI_COMM_RANK( comm_cyclic_fill, myid, ierr )
813	ENDIF
814	#else
815	remote_pe = remote_pe_s
816	rma_offset = rma_offset_s
817	myid = 0
818	numprocs = 1
819	#endif
820	!
821	!-- Allocate 2d buffers as RMA window, accessible on all threads
822	IF ( pe_active_for_read ) THEN
823	ALLOCATE( rmabuf_2di(nys:nyn,nxl:nxr) )
824	ELSE
825	ALLOCATE( rmabuf_2di(1,1) )
826	ENDIF
827	winsize = SIZE( rmabuf_2di ) * iwp
828
829	#if defined( __parallel )
830	CALL MPI_WIN_CREATE( rmabuf_2di, winsize, iwp, MPI_INFO_NULL, comm2d, rmawin_2di, ierr )
831	CALL MPI_WIN_FENCE( 0, rmawin_2di, ierr )
832	#endif
833
834	IF ( pe_active_for_read ) THEN
835	ALLOCATE( rmabuf_2d(nys:nyn,nxl:nxr) )
836	ELSE
837	ALLOCATE( rmabuf_2d(1,1) )
838	ENDIF
839	winsize = SIZE( rmabuf_2d ) * wp
840
841	#if defined( __parallel )
842	CALL MPI_WIN_CREATE( rmabuf_2d, winsize, wp, MPI_INFO_NULL, comm2d, rmawin_2d, ierr )
843	CALL MPI_WIN_FENCE( 0, rmawin_2d, ierr )
844	#endif
845
846	!
847	!-- Allocate 3d buffer as RMA window, accessable on all threads
848	IF ( pe_active_for_read ) THEN
849	ALLOCATE( rmabuf_3d(nzb:nzt+1,nys:nyn,nxl:nxr) )
850	ELSE
851	ALLOCATE( rmabuf_3d(1,1,1) )
852	ENDIF
853	winsize = SIZE( rmabuf_3d ) * wp
854
855	#if defined( __parallel )
856	CALL MPI_WIN_CREATE( rmabuf_3d, winsize, wp, MPI_INFO_NULL, comm2d, rmawin_3d, ierr )
857	CALL MPI_WIN_FENCE( 0, rmawin_3d, ierr )
858	#endif
859
860	!
861	!-- TODO: comment in more detail, what is done here, and why
862	!-- save small grid
863	CALL prerun_grid%save_grid_into_this_class()
864	prerun_grid%comm2d = comm_cyclic_fill
865
866	DEALLOCATE( remote_pe_s, rma_offset_s )
867
868	END SUBROUTINE setup_cyclic_fill
869
870	END SUBROUTINE rd_mpi_io_open
871
872
873	!--------------------------------------------------------------------------------------------------!
874	! Description:
875	! ------------
876	!> Check, if array exists in restart file
877	!--------------------------------------------------------------------------------------------------!
878	SUBROUTINE rd_mpi_io_check_array( name, found )
879
880	IMPLICIT NONE
881
882	CHARACTER(LEN=*), INTENT(IN) :: name !<
883
884	INTEGER(iwp) :: i !<
885
886	LOGICAl :: found !<
887
888
889	DO i = 1, tgh%nr_arrays
890	IF ( TRIM( array_names(i) ) == TRIM( name ) ) THEN
891	array_position = array_offset(i)
892	found = .TRUE.
893	RETURN
894	ENDIF
895	ENDDO
896
897	found = .FALSE.
898
899	END SUBROUTINE rd_mpi_io_check_array
900
901
902
903	!--------------------------------------------------------------------------------------------------!
904	! Description:
905	! ------------
906	!> Read INTEGER with MPI-IO
907	!--------------------------------------------------------------------------------------------------!
908	SUBROUTINE rrd_mpi_io_int( name, value )
909
910	IMPLICIT NONE
911
912	CHARACTER(LEN=*), INTENT(IN) :: name !<
913
914	INTEGER(iwp) :: i !<
915	INTEGER(KIND=iwp), INTENT(OUT) :: value !<
916
917	LOGICAL :: found !<
918
919
920	found = .FALSE.
921	value = 0
922
923	DO i = 1, tgh%nr_int
924	IF ( TRIM(int_names(i)) == TRIM( name ) ) THEN
925	value = int_values(i)
926	found = .TRUE.
927	EXIT
928	ENDIF
929	ENDDO
930
931	IF ( .NOT. found ) THEN
932	message_string = 'INTEGER variable "' // TRIM( name ) // '" not found in restart file'
933	CALL message( 'rrd_mpi_io_int', 'PA0722', 3, 2, 0, 6, 0 )
934	ENDIF
935
936	END SUBROUTINE rrd_mpi_io_int
937
938
939
940	!--------------------------------------------------------------------------------------------------!
941	! Description:
942	! ------------
943	!> Read REAL with MPI-IO
944	!--------------------------------------------------------------------------------------------------!
945	SUBROUTINE rrd_mpi_io_real( name, value )
946
947	IMPLICIT NONE
948
949	CHARACTER(LEN=*), INTENT(IN) :: name !<
950
951	INTEGER(iwp) :: i !<
952
953	LOGICAL :: found !<
954
955	REAL(KIND=wp), INTENT(OUT) :: value !<
956
957
958	found = .FALSE.
959	value = 0.0
960
961	DO i = 1, tgh%nr_real
962	IF ( TRIM(real_names(i)) == TRIM( name ) ) THEN
963	value = real_values(i)
964	found = .TRUE.
965	EXIT
966	ENDIF
967	ENDDO
968
969	IF ( .NOT. found ) THEN
970	message_string = 'REAL variable "' // TRIM( name ) // '" not found in restart file'
971	CALL message( 'rrd_mpi_io_int', 'PA0722', 3, 2, 0, 6, 0 )
972	ENDIF
973
974	END SUBROUTINE rrd_mpi_io_real
975
976
977
978	!--------------------------------------------------------------------------------------------------!
979	! Description:
980	! ------------
981	!> Read 2d-real array with MPI-IO
982	!--------------------------------------------------------------------------------------------------!
983	SUBROUTINE rrd_mpi_io_real_2d( name, data )
984
985	IMPLICIT NONE
986
987	CHARACTER(LEN=*), INTENT(IN) :: name !<
988
989	#if defined( __parallel )
990	INTEGER, DIMENSION(rd_status_size) :: status !<
991	#endif
992	INTEGER(iwp) :: i !<
993
994	LOGICAL :: found !<
995
996	REAL(wp), INTENT(INOUT), DIMENSION(nysg:nyng,nxlg:nxrg) :: data !<
997
998
999	found = .FALSE.
1000
1001	DO i = 1, tgh%nr_arrays
1002	IF ( TRIM(array_names(i)) == TRIM( name ) ) THEN
1003	array_position = array_offset(i)
1004	found = .TRUE.
1005	EXIT
1006	ENDIF
1007	ENDDO
1008
1009	IF ( found ) THEN
1010
1011	IF ( cyclic_fill_mode ) THEN
1012
1013	CALL rrd_mpi_io_real_2d_cyclic_fill
1014
1015	ELSE
1016
1017	#if defined( __parallel )
1018	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited # of cores is inactive
1019	IF ( sm_io%iam_io_pe ) THEN
1020	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL, &
1021	ierr )
1022	CALL MPI_FILE_READ_ALL( fh, array_2d, SIZE( array_2d ), MPI_REAL, status, ierr )
1023	ENDIF
1024	CALL sm_io%sm_node_barrier()
1025	#else
1026	CALL posix_lseek( fh, array_position )
1027	CALL posix_read( fh, array_2d, SIZE( array_2d ) )
1028	#endif
1029
1030	IF ( include_total_domain_boundaries ) THEN
1031	DO i = iog%nxl, iog%nxr
1032	data(iog%nys-nbgp:iog%nyn-nbgp,i-nbgp) = array_2d(i,iog%nys:iog%nyn)
1033	ENDDO
1034	IF ( debug_level >= 2) THEN
1035	WRITE(9,*) 'r2f_ob ', TRIM(name),' ', SUM( data(nys:nyn,nxl:nxr) )
1036	ENDIF
1037	ELSE
1038	DO i = nxl, nxr
1039	data(nys:nyn,i) = array_2d(i,nys:nyn)
1040	ENDDO
1041	IF ( debug_level >= 2) THEN
1042	WRITE(9,*) 'r2f ', TRIM( name ),' ', SUM( data(nys:nyn,nxl:nxr) )
1043	ENDIF
1044	ENDIF
1045
1046	ENDIF
1047
1048	CALL exchange_horiz_2d( data )
1049
1050	ELSE
1051	message_string = '2d-REAL array "' // TRIM( name ) // '" not found in restart file'
1052	CALL message( 'rrd_mpi_io_int', 'PA0722', 3, 2, 0, 6, 0 )
1053	ENDIF
1054
1055
1056	CONTAINS
1057
1058	SUBROUTINE rrd_mpi_io_real_2d_cyclic_fill
1059
1060	IMPLICIT NONE
1061
1062	INTEGER(iwp) :: i !<
1063	INTEGER(iwp) :: ie !<
1064	INTEGER(iwp) :: ierr !<
1065	INTEGER(iwp) :: is !<
1066	INTEGER(iwp) :: i_remote !<
1067	INTEGER(iwp) :: j !<
1068	INTEGER(iwp) :: je !<
1069	INTEGER(iwp) :: js !<
1070	INTEGER(iwp) :: j_remote !<
1071	INTEGER(iwp) :: nval !<
1072	INTEGER(iwp) :: rem_pe !<
1073
1074	INTEGER(KIND=MPI_ADDRESS_KIND) :: rem_offs !<
1075
1076
1077	!kk write(9,*) 'Here is rma_cylic_fill_real_2d ',nxl,nxr,nys,nyn; FLUSH(9)
1078
1079	!
1080	!-- Reading 2d real array on prerun grid
1081	CALL prerun_grid%activate_grid_from_this_class()
1082
1083	IF ( pe_active_for_read ) THEN
1084
1085	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL, &
1086	ierr )
1087	CALL MPI_FILE_READ_ALL( fh, array_2d, SIZE( array_2d ), MPI_REAL, status, ierr )
1088
1089	DO i = nxl, nxr
1090	rmabuf_2d(nys:nyn,i) = array_2d(i,nys:nyn)
1091	ENDDO
1092	data(nys:nyn,nxl:nxr) = rmabuf_2d ! copy prerund data directly into output array data
1093	ENDIF
1094
1095	CALL mainrun_grid%activate_grid_from_this_class()
1096
1097	#if defined( __parallel )
1098	!
1099	!-- Close RMA window to allow remote access
1100	CALL MPI_WIN_FENCE( 0, rmawin_2d, ierr )
1101	#endif
1102
1103	!
1104	!-- TODO: describe in more detail what is happening in this IF/ELSE clause
1105	IF ( .NOT. pe_active_for_read ) THEN
1106
1107	is = nxl
1108	ie = nxr
1109	js = nys
1110	je = nyn
1111
1112	ELSE
1113	!
1114	!-- Extra get for cyclic data north of prerun data
1115	is = nxl
1116	ie = nxr
1117	js = prerun_grid%nys+1
1118	je = nyn
1119	DO i = is, ie
1120	DO j = js, je
1121	i_remote = MOD(i,nx_on_file+1)
1122	j_remote = MOD(j,ny_on_file+1)
1123	rem_pe = remote_pe(i_remote,j_remote)
1124	rem_offs = rma_offset(i_remote,j_remote)
1125	nval = 1
1126
1127	#if defined( __parallel )
1128	IF ( rem_pe /= myid ) THEN
1129	CALL MPI_GET( data(j,i), nval, MPI_REAL, rem_pe, rem_offs, nval, MPI_REAL, &
1130	rmawin_2d, ierr )
1131	ELSE
1132	data(j,i) = rmabuf_2d(j_remote,i_remote)
1133	ENDIF
1134	#else
1135	data(j,i) = array_2d(i_remote,j_remote)
1136	#endif
1137	ENDDO
1138	ENDDO
1139	!
1140	!-- Prepare setup for stripe right of prerun data
1141	is = prerun_grid%nxr+1
1142	ie = nxr
1143	js = nys
1144	je = nyn
1145
1146	ENDIF
1147
1148	DO i = is, ie
1149	DO j = js, je
1150	i_remote = MOD(i,nx_on_file+1)
1151	j_remote = MOD(j,ny_on_file+1)
1152	rem_pe = remote_pe(i_remote,j_remote)
1153	rem_offs = rma_offset(i_remote,j_remote)
1154	nval = 1
1155
1156	#if defined( __parallel )
1157	IF ( rem_pe /= myid ) THEN
1158	CALL MPI_GET( data(j,i), nval, MPI_REAL, rem_pe, rem_offs, nval, MPI_REAL, &
1159	rmawin_2d, ierr )
1160	ELSE
1161	data(j,i) = rmabuf_2d(j_remote,i_remote)
1162	ENDIF
1163	#else
1164	data(j,i) = array_2d(i_remote,j_remote)
1165	#endif
1166	ENDDO
1167	ENDDO
1168
1169	#if defined( __parallel )
1170	!
1171	!-- Reopen RMA window to allow filling
1172	CALL MPI_WIN_FENCE( 0, rmawin_2d, ierr )
1173	#endif
1174
1175	END SUBROUTINE rrd_mpi_io_real_2d_cyclic_fill
1176
1177	END SUBROUTINE rrd_mpi_io_real_2d
1178
1179
1180
1181	!--------------------------------------------------------------------------------------------------!
1182	! Description:
1183	! ------------
1184	!> Read 2d-INTEGER array with MPI-IO
1185	!--------------------------------------------------------------------------------------------------!
1186	SUBROUTINE rrd_mpi_io_int_2d( name, data )
1187
1188	IMPLICIT NONE
1189
1190	CHARACTER(LEN=*), INTENT(IN) :: name !<
1191
1192	INTEGER(iwp) :: i !<
1193	INTEGER(iwp) :: j !<
1194
1195	#if defined( __parallel )
1196	INTEGER, DIMENSION(rd_status_size) :: status !<
1197	#endif
1198
1199	INTEGER(KIND=iwp), INTENT(INOUT), DIMENSION(:,:) :: data !<
1200
1201	LOGICAL :: found !<
1202
1203
1204	found = .FALSE.
1205
1206	DO i = 1, tgh%nr_arrays
1207	IF ( TRIM(array_names(i)) == TRIM( name ) ) THEN
1208	array_position = array_offset(i)
1209	found = .TRUE.
1210	EXIT
1211	ENDIF
1212	ENDDO
1213
1214	IF ( found ) THEN
1215
1216	IF ( ( nxr - nxl + 1 + 2 * nbgp ) == SIZE( data, 2 ) ) THEN
1217	!
1218	!-- Output array with Halos.
1219	!-- ATTENTION: INTEGER arrays with ghost boundaries are not implemented yet. This kind of
1220	!-- array would be dimensioned in the caller subroutine like this:
1221	!-- INTEGER, DIMENSION(nysg:nyng,nxlg:nxrg):: data
1222	message_string = '2d-INTEGER array "' // TRIM( name ) // '" to be read from restart ' // &
1223	'file is defined with illegal dimensions in the PALM code'
1224	CALL message( 'rrd_mpi_io_int_2d', 'PA0723', 3, 2, 0, 6, 0 )
1225
1226	ELSEIF ( (nxr - nxl + 1) == SIZE( data, 2 ) ) THEN
1227	!
1228	!-- INTEGER input array without Halos.
1229	!-- This kind of array is dimensioned in the caller subroutine
1230	!-- INTEGER, DIMENSION(nys:nyn,nxl:nxr) :: data
1231	IF ( cyclic_fill_mode ) THEN
1232
1233	CALL rrd_mpi_io_int_2d_cyclic_fill
1234
1235	ELSE
1236
1237	#if defined( __parallel )
1238	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited # of cores is inactive
1239	IF ( sm_io%iam_io_pe ) THEN
1240	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native', &
1241	MPI_INFO_NULL, ierr )
1242	CALL MPI_FILE_READ_ALL( fh, array_2di, SIZE( array_2di ), MPI_INTEGER, status, &
1243	ierr )
1244	ENDIF
1245	CALL sm_io%sm_node_barrier()
1246	#else
1247	CALL posix_lseek( fh, array_position )
1248	CALL posix_read( fh, array_2di, SIZE( array_2di ) )
1249	#endif
1250	DO j = nys, nyn
1251	DO i = nxl, nxr
1252	data(j-nys+1,i-nxl+1) = array_2di(i,j)
1253	ENDDO
1254	ENDDO
1255
1256	ENDIF
1257
1258	ELSE
1259
1260	message_string = '2d-INTEGER array "' // TRIM( name ) // '" to be read from restart ' // &
1261	'file is defined with illegal dimensions in the PALM code'
1262	CALL message( 'rrd_mpi_io_int_2d', 'PA0723', 3, 2, 0, 6, 0 )
1263
1264	ENDIF
1265
1266	ELSE
1267
1268	message_string = '2d-INTEGER array "' // TRIM( name ) // '" not found in restart file'
1269	CALL message( 'rrd_mpi_io_int_2d', 'PA0722', 3, 2, 0, 6, 0 )
1270
1271	ENDIF
1272
1273
1274	CONTAINS
1275
1276	SUBROUTINE rrd_mpi_io_int_2d_cyclic_fill
1277
1278	IMPLICIT NONE
1279
1280	INTEGER(iwp) :: i !<
1281	INTEGER(iwp) :: ie !<
1282	INTEGER(iwp) :: ierr !<
1283	INTEGER(iwp) :: is !<
1284	INTEGER(iwp) :: i_remote !<
1285	INTEGER(iwp) :: j !<
1286	INTEGER(iwp) :: je !<
1287	INTEGER(iwp) :: js !<
1288	INTEGER(iwp) :: j_remote !<
1289	INTEGER(iwp) :: nval !<
1290	INTEGER(iwp) :: rem_pe !<
1291
1292	INTEGER(KIND=MPI_ADDRESS_KIND) :: rem_offs !<
1293
1294
1295	CALL prerun_grid%activate_grid_from_this_class()
1296
1297	IF ( pe_active_for_read ) THEN
1298	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native', &
1299	MPI_INFO_NULL, ierr )
1300	CALL MPI_FILE_READ_ALL( fh, array_2di, SIZE( array_2di ), MPI_INTEGER, status, ierr )
1301
1302	DO i = nxl, nxr
1303	rmabuf_2di(nys:nyn,i) = array_2di(i,nys:nyn)
1304	ENDDO
1305	data(1:nny,1:nnx) = rmabuf_2di
1306	ENDIF
1307
1308	CALL mainrun_grid%activate_grid_from_this_class()
1309
1310	#if defined( __parallel )
1311	!
1312	!-- Close RMA window to allow remote access
1313	CALL MPI_WIN_FENCE( 0, rmawin_2di, ierr )
1314	#endif
1315
1316	IF ( .NOT. pe_active_for_read ) THEN
1317
1318	is = nxl
1319	ie = nxr
1320	js = nys
1321	je = nyn
1322
1323	ELSE
1324
1325	is = nxl
1326	ie = nxr
1327	js = prerun_grid%nys+1
1328	je = nyn
1329	DO i = is, ie
1330	DO j = js, je
1331	i_remote = MOD(i,nx_on_file+1)
1332	j_remote = MOD(j,ny_on_file+1)
1333	rem_pe = remote_pe(i_remote,j_remote)
1334	rem_offs = rma_offset(i_remote,j_remote)
1335	nval = 1
1336
1337	#if defined( __parallel )
1338	IF ( rem_pe /= myid ) THEN
1339	CALL MPI_GET( data(j-nys+1,i-nxl+1), nval, MPI_INTEGER, rem_pe, rem_offs, nval, &
1340	MPI_INTEGER, rmawin_2di, ierr )
1341	ELSE
1342	data(j-nys+1,i-nxl+1) = rmabuf_2di(j_remote,i_remote)
1343	ENDIF
1344	#else
1345	data(j-nys+1,i-nxl+1) = array_2di(i_remote,j_remote)
1346	#endif
1347	ENDDO
1348	ENDDO
1349	is = prerun_grid%nxr+1
1350	ie = nxr
1351	js = nys
1352	je = nyn
1353
1354	ENDIF
1355
1356	DO i = is, ie
1357	DO j = js, je
1358	i_remote = MOD(i,nx_on_file+1)
1359	j_remote = MOD(j,ny_on_file+1)
1360	rem_pe = remote_pe(i_remote,j_remote)
1361	rem_offs = rma_offset(i_remote,j_remote)
1362	nval = 1
1363	#if defined( __parallel )
1364	IF ( rem_pe /= myid ) THEN
1365	CALL MPI_GET( data(j-nys+1,i-nxl+1), nval, MPI_INTEGER, rem_pe, rem_offs, nval, &
1366	MPI_INTEGER, rmawin_2di, ierr)
1367	ELSE
1368	data(j-nys+1,i-nxl+1) = rmabuf_2di(j_remote,i_remote)
1369	ENDIF
1370	#else
1371	data(j-nys+1,i-nxl+1) = array_2di(i_remote,j_remote)
1372	#endif
1373	ENDDO
1374	ENDDO
1375
1376	#if defined( __parallel )
1377	!
1378	!-- Reopen RMA window to allow filling
1379	CALL MPI_WIN_FENCE( 0, rmawin_2di, ierr )
1380	#endif
1381
1382	END SUBROUTINE rrd_mpi_io_int_2d_cyclic_fill
1383
1384	END SUBROUTINE rrd_mpi_io_int_2d
1385
1386
1387
1388	!--------------------------------------------------------------------------------------------------!
1389	! Description:
1390	! ------------
1391	!> Read 2d-REAL array with MPI-IO
1392	!--------------------------------------------------------------------------------------------------!
1393	SUBROUTINE rrd_mpi_io_real_3d( name, data )
1394
1395	IMPLICIT NONE
1396
1397	CHARACTER(LEN=*), INTENT(IN) :: name !<
1398
1399	INTEGER(iwp) :: i !<
1400	INTEGER(iwp) :: j !<
1401
1402	#if defined( __parallel )
1403	INTEGER, DIMENSION(rd_status_size) :: status !<
1404	#endif
1405
1406	LOGICAL :: found !<
1407
1408	REAL(wp), INTENT(INOUT), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: data !<
1409
1410
1411	found = .FALSE.
1412	data = -1.0
1413
1414	DO i = 1, tgh%nr_arrays
1415	IF ( TRIM(array_names(i)) == TRIM( name ) ) THEN
1416	array_position = array_offset(i)
1417	found = .TRUE.
1418	EXIT
1419	ENDIF
1420	ENDDO
1421
1422	IF ( found ) THEN
1423
1424	IF ( cyclic_fill_mode ) THEN
1425
1426	CALL rrd_mpi_io_real_3d_cyclic_fill
1427	!
1428	!-- Cyclic fill mode requires to use the "cyclic" communicator, in order to initialize
1429	!-- grid points at the outer boundaries (ghost layers) of the total domain. These points
1430	!-- are not contained in the prerun data, because the prerun used cyclic boundary conditions.
1431	CALL exchange_horiz( data, nbgp, alternative_communicator = 1 )
1432
1433	ELSE
1434	#if defined( __parallel )
1435	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited # of cores is inactive
1436	IF( sm_io%iam_io_pe ) THEN
1437	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL, &
1438	ierr )
1439	CALL MPI_FILE_READ_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
1440	ENDIF
1441	CALL sm_io%sm_node_barrier()
1442	#else
1443	CALL posix_lseek( fh, array_position )
1444	CALL posix_read(fh, array_3d, SIZE( array_3d ) )
1445	#endif
1446	IF ( include_total_domain_boundaries ) THEN
1447	DO i = iog%nxl, iog%nxr
1448	data(:,iog%nys-nbgp:iog%nyn-nbgp,i-nbgp) = array_3d(:,i,iog%nys:iog%nyn)
1449	ENDDO
1450	ELSE
1451	DO i = nxl, nxr
1452	data(:,nys:nyn,i) = array_3d(:,i,nys:nyn)
1453	ENDDO
1454	ENDIF
1455
1456	CALL exchange_horiz( data, nbgp )
1457
1458	ENDIF
1459
1460	ELSE
1461
1462	message_string = '3d-REAL array "' // TRIM( name ) // '" not found in restart file'
1463	CALL message( 'rrd_mpi_io_real_3d', 'PA0722', 3, 2, 0, 6, 0 )
1464
1465	ENDIF
1466
1467
1468	CONTAINS
1469
1470	SUBROUTINE rrd_mpi_io_real_3d_cyclic_fill
1471
1472	IMPLICIT NONE
1473
1474	INTEGER(iwp) :: i !<
1475	INTEGER(iwp) :: ie !<
1476	INTEGER(iwp) :: ierr !<
1477	INTEGER(iwp) :: is !<
1478	INTEGER(iwp) :: i_remote !<
1479	INTEGER(iwp) :: j !<
1480	INTEGER(iwp) :: je !<
1481	INTEGER(iwp) :: js !<
1482	INTEGER(iwp) :: j_remote !<
1483	INTEGER(iwp) :: nval !<
1484	INTEGER(iwp) :: rem_pe !<
1485
1486	INTEGER(KIND=MPI_ADDRESS_KIND) :: rem_offs !<
1487
1488
1489	CALL prerun_grid%activate_grid_from_this_class()
1490
1491	IF ( pe_active_for_read ) THEN
1492	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL, &
1493	ierr )
1494	CALL MPI_FILE_READ_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
1495
1496	DO i = nxl, nxr
1497	rmabuf_3d(:,nys:nyn,i) = array_3d(:,i,nys:nyn)
1498	ENDDO
1499	data(:,nys:nyn,nxl:nxr) = rmabuf_3d
1500	ENDIF
1501	CALL mainrun_grid%activate_grid_from_this_class ()
1502
1503	#if defined( __parallel )
1504	!
1505	!-- Close RMA window to allow remote access
1506	CALL MPI_WIN_FENCE( 0, rmawin_3d, ierr )
1507	#endif
1508
1509	IF ( .NOT. pe_active_for_read ) THEN
1510
1511	is = nxl
1512	ie = nxr
1513	js = nys
1514	je = nyn
1515
1516	ELSE
1517
1518	is = nxl
1519	ie = nxr
1520	js = prerun_grid%nys+1
1521	je = nyn
1522
1523	DO i = is, ie
1524	DO j = js, je
1525	i_remote = MOD(i,nx_on_file+1)
1526	j_remote = MOD(j,ny_on_file+1)
1527	rem_pe = remote_pe(i_remote,j_remote)
1528	rem_offs = rma_offset(i_remote,j_remote)*(nzt-nzb+2)
1529	nval = nzt-nzb+2
1530
1531	#if defined( __parallel )
1532	IF(rem_pe /= myid) THEN
1533	CALL MPI_GET( data(nzb,j,i), nval, MPI_REAL, rem_pe, rem_offs, nval, MPI_REAL, &
1534	rmawin_3d, ierr)
1535	ELSE
1536	data(:,j,i) = rmabuf_3d(:,j_remote,i_remote)
1537	ENDIF
1538	#else
1539	data(:,j,i) = array_3d(:,i_remote,j_remote)
1540	#endif
1541	ENDDO
1542	ENDDO
1543	is = prerun_grid%nxr+1
1544	ie = nxr
1545	js = nys
1546	je = nyn
1547
1548	ENDIF
1549
1550	DO i = is, ie
1551	DO j = js, je
1552	i_remote = MOD(i,nx_on_file+1)
1553	j_remote = MOD(j,ny_on_file+1)
1554	rem_pe = remote_pe(i_remote,j_remote)
1555	rem_offs = rma_offset(i_remote,j_remote) * ( nzt-nzb+2 )
1556	nval = nzt-nzb+2
1557
1558	#if defined( __parallel )
1559	IF ( rem_pe /= myid ) THEN
1560	CALL MPI_GET( data(nzb,j,i), nval, MPI_REAL, rem_pe, rem_offs, nval, MPI_REAL, &
1561	rmawin_3d, ierr)
1562	ELSE
1563	data(:,j,i) = rmabuf_3d(:,j_remote,i_remote)
1564	ENDIF
1565	#else
1566	data(:,j,i) = array_3d(:,i_remote,j_remote)
1567	#endif
1568	ENDDO
1569	ENDDO
1570
1571	#if defined( __parallel )
1572	!
1573	!-- Reopen RMA window to allow filling
1574	CALL MPI_WIN_FENCE( 0, rmawin_3d, ierr )
1575	#endif
1576
1577	END SUBROUTINE rrd_mpi_io_real_3d_cyclic_fill
1578
1579	END SUBROUTINE rrd_mpi_io_real_3d
1580
1581
1582
1583	!--------------------------------------------------------------------------------------------------!
1584	! Description:
1585	! ------------
1586	!> Read 3d-REAL soil array with MPI-IO
1587	!> nzb_soil, nzt_soil are located in the module land_surface_model_mod. Since Fortran does not allow
1588	!> cross referencing of module variables, it is required to pass these variables as arguments.
1589	!--------------------------------------------------------------------------------------------------!
1590	SUBROUTINE rrd_mpi_io_real_3d_soil( name, data, nzb_soil, nzt_soil )
1591
1592	IMPLICIT NONE
1593
1594	CHARACTER(LEN=*), INTENT(IN) :: name !<
1595
1596	INTEGER(iwp) :: i !<
1597	INTEGER, INTENT(IN) :: nzb_soil !<
1598	INTEGER, INTENT(IN) :: nzt_soil !<
1599
1600	#if defined( __parallel )
1601	INTEGER, DIMENSION(rd_status_size) :: status !<
1602	#endif
1603
1604	LOGICAL :: found !<
1605	INTEGER(iwp) :: ierr !<
1606
1607	REAL(wp), INTENT(INOUT), DIMENSION(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) :: data !<
1608
1609
1610	!
1611	!-- Prerun data is not allowed to contain soil information so far
1612	IF ( cyclic_fill_mode ) THEN
1613	message_string = 'prerun data is not allowed to contain soil information'
1614	CALL message( 'rrd_mpi_io_real_3d_soil', 'PA0729', 3, 2, -1, 6, 0 )
1615	ENDIF
1616
1617	found = .FALSE.
1618
1619	DO i = 1, tgh%nr_arrays
1620	IF ( TRIM(array_names(i)) == TRIM( name ) ) THEN
1621	array_position = array_offset(i)
1622	found = .TRUE.
1623	EXIT
1624	ENDIF
1625	ENDDO
1626
1627	IF ( found ) THEN
1628	#if defined( __parallel )
1629	CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
1630	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited number of cores is inactive
1631	IF ( sm_io%iam_io_pe ) THEN
1632	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3dsoil, 'native', &
1633	MPI_INFO_NULL, ierr )
1634	CALL MPI_FILE_READ_ALL( fh, array_3d_soil, SIZE( array_3d_soil ), MPI_REAL, status, ierr )
1635	CALL MPI_TYPE_FREE( ft_3dsoil, ierr )
1636	ENDIF
1637	CALL sm_io%sm_node_barrier()
1638	#else
1639	CALL posix_lseek( fh, array_position )
1640	CALL posix_read( fh, array_3d_soil, SIZE( array_3d_soil ) )
1641	#endif
1642	IF ( include_total_domain_boundaries ) THEN
1643	DO i = iog%nxl, iog%nxr
1644	data(:,iog%nys-nbgp:iog%nyn-nbgp,i-nbgp) = array_3d_soil(:,i,iog%nys:iog%nyn)
1645	ENDDO
1646	ELSE
1647	DO i = nxl, nxr
1648	data(:,nys:nyn,i) = array_3d_soil(:,i,nys:nyn)
1649	ENDDO
1650	ENDIF
1651
1652	ELSE
1653
1654	message_string = '3d-REAL soil array "' // TRIM( name ) // '" not found in restart file'
1655	CALL message( 'rrd_mpi_io_real_3d_soil', 'PA0722', 3, 2, 0, 6, 0 )
1656
1657	ENDIF
1658
1659	END SUBROUTINE rrd_mpi_io_real_3d_soil
1660
1661
1662
1663	!--------------------------------------------------------------------------------------------------!
1664	! Description:
1665	! ------------
1666	!> Read CHARACTER with MPI-IO
1667	!--------------------------------------------------------------------------------------------------!
1668	SUBROUTINE rrd_mpi_io_char( name, text )
1669
1670	IMPLICIT NONE
1671
1672	CHARACTER(LEN=*), INTENT(IN) :: name !<
1673	CHARACTER(LEN=*), INTENT(OUT) :: text !<
1674	CHARACTER(LEN=128) :: line !<
1675
1676	INTEGER(iwp) :: i !<
1677
1678	LOGICAL :: found !<
1679
1680
1681	found = .FALSE.
1682	text = ' '
1683
1684	DO i = 1, tgh%nr_char
1685	line = text_lines(i)
1686	IF ( TRIM( line(1:32) ) == TRIM( name ) ) THEN
1687	text = line(33:)
1688	found = .TRUE.
1689	EXIT
1690	ENDIF
1691	ENDDO
1692
1693	IF ( .NOT. found ) THEN
1694	message_string = 'CHARACTER variable "' // TRIM( name ) // '" not found in restart file'
1695	CALL message( 'rrd_mpi_io_char', 'PA0722', 3, 2, 0, 6, 0 )
1696	ENDIF
1697
1698	END SUBROUTINE rrd_mpi_io_char
1699
1700
1701
1702	!--------------------------------------------------------------------------------------------------!
1703	! Description:
1704	! ------------
1705	!> Read LOGICAL with MPI-IO
1706	!--------------------------------------------------------------------------------------------------!
1707	SUBROUTINE rrd_mpi_io_logical( name, value )
1708
1709	IMPLICIT NONE
1710
1711	CHARACTER(LEN=*), INTENT(IN) :: name !<
1712
1713	INTEGER(iwp) :: logical_as_integer !<
1714
1715	LOGICAL, INTENT(OUT) :: value !<
1716
1717
1718	CALL rrd_mpi_io_int( name, logical_as_integer )
1719	value = ( logical_as_integer == 1 )
1720
1721	END SUBROUTINE rrd_mpi_io_logical
1722
1723
1724
1725	!--------------------------------------------------------------------------------------------------!
1726	! Description:
1727	! ------------
1728	!> Write INTEGER with MPI-IO
1729	!--------------------------------------------------------------------------------------------------!
1730	SUBROUTINE wrd_mpi_io_int( name, value )
1731
1732	IMPLICIT NONE
1733
1734	CHARACTER(LEN=*), INTENT(IN) :: name !<
1735
1736	INTEGER(KIND=iwp), INTENT(IN) :: value !<
1737
1738
1739	IF ( header_int_index == max_nr_int ) THEN
1740	STOP '+++ maximum number of INTEGER entries in restart file header exceeded'
1741	ENDIF
1742
1743	int_names(header_int_index) = name
1744	int_values(header_int_index) = value
1745	header_int_index = header_int_index + 1
1746
1747	END SUBROUTINE wrd_mpi_io_int
1748
1749
1750	!--------------------------------------------------------------------------------------------------!
1751	! Description:
1752	! ------------
1753	!> To do: Description missing!
1754	!--------------------------------------------------------------------------------------------------!
1755	SUBROUTINE wrd_mpi_io_real( name, value )
1756
1757	IMPLICIT NONE
1758
1759	CHARACTER(LEN=*), INTENT(IN) :: name !<
1760
1761	REAL(wp), INTENT(IN) :: value !<
1762
1763
1764	IF ( header_real_index == max_nr_real ) THEN
1765	STOP '+++ maximum number of REAL entries in restart file header exceeded'
1766	ENDIF
1767
1768	real_names(header_real_index) = name
1769	real_values(header_real_index) = value
1770	header_real_index = header_real_index + 1
1771
1772	END SUBROUTINE wrd_mpi_io_real
1773
1774
1775
1776	!--------------------------------------------------------------------------------------------------!
1777	! Description:
1778	! ------------
1779	!> Write 2d-REAL array with MPI-IO
1780	!--------------------------------------------------------------------------------------------------!
1781	SUBROUTINE wrd_mpi_io_real_2d( name, data )
1782
1783	IMPLICIT NONE
1784
1785	CHARACTER(LEN=*), INTENT(IN) :: name !<
1786
1787	INTEGER(iwp) :: i !<
1788
1789	#if defined( __parallel )
1790	INTEGER, DIMENSION(rd_status_size) :: status !<
1791	#endif
1792
1793	REAL(wp), INTENT(IN), DIMENSION(nysg:nyng,nxlg:nxrg) :: data !<
1794
1795
1796	IF ( header_array_index == max_nr_arrays ) THEN
1797	STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
1798	ENDIF
1799
1800	array_names(header_array_index) = name
1801	array_offset(header_array_index) = array_position
1802	header_array_index = header_array_index + 1
1803
1804	IF ( include_total_domain_boundaries ) THEN
1805	!
1806	!-- Prepare output with outer boundaries
1807	DO i = iog%nxl, iog%nxr
1808	array_2d(i,iog%nys:iog%nyn) = data(iog%nys-nbgp:iog%nyn-nbgp,i-nbgp)
1809	ENDDO
1810
1811	ELSE
1812	!
1813	!-- Prepare output without outer boundaries
1814	DO i = nxl,nxr
1815	array_2d(i,iog%nys:iog%nyn) = data(nys:nyn,i)
1816	ENDDO
1817
1818	ENDIF
1819
1820	#if defined( __parallel )
1821	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited number of cores is inactive
1822	IF ( sm_io%iam_io_pe ) THEN
1823	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL, ierr )
1824	CALL MPI_FILE_WRITE_ALL( fh, array_2d, SIZE( array_2d), MPI_REAL, status, ierr )
1825	ENDIF
1826	CALL sm_io%sm_node_barrier()
1827	#else
1828	CALL posix_lseek( fh, array_position )
1829	CALL posix_write( fh, array_2d, SIZE( array_2d ) )
1830	#endif
1831	!
1832	!-- Type conversion required, otherwise right hand side brackets are calculated assuming 4 byte INT.
1833	!-- Maybe a compiler problem.
1834	array_position = array_position + ( INT( iog%ny, KIND=rd_offset_kind ) + 1 ) * &
1835	( INT( iog%nx, KIND=rd_offset_kind ) + 1 ) * wp
1836
1837	END SUBROUTINE wrd_mpi_io_real_2d
1838
1839
1840
1841	!--------------------------------------------------------------------------------------------------!
1842	! Description:
1843	! ------------
1844	!> Write 2d-INTEGER array with MPI-IO
1845	!--------------------------------------------------------------------------------------------------!
1846	SUBROUTINE wrd_mpi_io_int_2d( name, data )
1847
1848	IMPLICIT NONE
1849
1850	CHARACTER(LEN=*), INTENT(IN) :: name !<
1851
1852	INTEGER(iwp) :: i !<
1853	INTEGER(iwp) :: j !<
1854
1855	#if defined( __parallel )
1856	INTEGER, DIMENSION(rd_status_size) :: status !<
1857	#endif
1858	INTEGER(KIND=iwp), INTENT(IN), DIMENSION(:,:) :: data !<
1859
1860
1861	IF ( header_array_index == max_nr_arrays ) THEN
1862	STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
1863	ENDIF
1864
1865	array_names(header_array_index) = name
1866	array_offset(header_array_index) = array_position
1867	header_array_index = header_array_index + 1
1868
1869	IF ( ( nxr - nxl + 1 + 2 * nbgp ) == SIZE( data, 2 ) ) THEN
1870	!
1871	!-- Integer arrays with ghost layers are not implemented yet. These kind of arrays would be
1872	!-- dimensioned in the caller subroutine as
1873	!-- INTEGER, DIMENSION(nysg:nyng,nxlg:nxrg) :: data
1874	message_string = '2d-INTEGER array "' // TRIM( name ) // '" to be written to restart ' // &
1875	'file is defined with illegal dimensions in the PALM code'
1876	CALL message( 'wrd_mpi_io_int_2d', 'PA0723', 3, 2, 0, 6, 0 )
1877
1878	ELSEIF ( ( nxr-nxl+1 ) == SIZE( data, 2 ) ) THEN
1879	!
1880	!-- INTEGER input array without ghost layers.
1881	!-- This kind of array is dimensioned in the caller subroutine as
1882	!-- INTEGER, DIMENSION(nys:nyn,nxl:nxr) :: data
1883	DO j = nys, nyn
1884	DO i = nxl, nxr
1885	array_2di(i,j) = data(j-nys+1,i-nxl+1)
1886	ENDDO
1887	ENDDO
1888	#if defined( __parallel )
1889	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited number of cores is inactive
1890	IF ( sm_io%iam_io_pe ) THEN
1891	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native', &
1892	MPI_INFO_NULL, ierr )
1893	CALL MPI_FILE_WRITE_ALL( fh, array_2di, SIZE( array_2di ), MPI_INTEGER, status, ierr )
1894	ENDIF
1895	CALL sm_io%sm_node_barrier()
1896	#else
1897	CALL posix_lseek( fh, array_position )
1898	CALL posix_write( fh, array_2di, SIZE( array_2di ) )
1899	#endif
1900	!
1901	!-- Type conversion required, otherwise rigth hand side brackets are calculated assuming 4 byte
1902	!-- INT. Maybe a compiler problem.
1903	array_position = array_position + INT( (ny+1), KIND=rd_offset_kind ) * &
1904	INT( (nx+1), KIND=rd_offset_kind ) * 4
1905
1906	ELSE
1907
1908	message_string = '2d-INTEGER array "' // TRIM( name ) // '" to be written to restart ' // &
1909	'file is defined with illegal dimensions in the PALM code'
1910	CALL message( 'wrd_mpi_io_int_2d', 'PA0723', 3, 2, 0, 6, 0 )
1911
1912	ENDIF
1913
1914	END SUBROUTINE wrd_mpi_io_int_2d
1915
1916
1917
1918	!--------------------------------------------------------------------------------------------------!
1919	! Description:
1920	! ------------
1921	!> Write 3d-REAL array with MPI-IO
1922	!--------------------------------------------------------------------------------------------------!
1923	SUBROUTINE wrd_mpi_io_real_3d( name, data )
1924
1925	IMPLICIT NONE
1926
1927	CHARACTER(LEN=*), INTENT(IN) :: name !<
1928
1929	INTEGER(iwp) :: i !<
1930	#if defined( __parallel )
1931	INTEGER, DIMENSION(rd_status_size) :: status !<
1932	#endif
1933	REAL(wp), INTENT(IN), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: data !<
1934
1935
1936	IF ( header_array_index == max_nr_arrays ) THEN
1937	STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
1938	ENDIF
1939
1940	array_names(header_array_index) = name
1941	array_offset(header_array_index) = array_position
1942	header_array_index = header_array_index + 1
1943
1944	IF ( include_total_domain_boundaries ) THEN
1945	!
1946	!-- Prepare output of 3d-REAL-array with ghost layers. In the virtual PE grid, the first
1947	!-- dimension is PEs along x, and the second along y. For MPI-IO it is recommended to have the
1948	!-- index order of the array in the same way, i.e. the first dimension should be along x and the
1949	!-- second along y. For this reason, the original PALM data need to be swaped.
1950	DO i = iog%nxl, iog%nxr
1951	array_3d(:,i,iog%nys:iog%nyn) = data(:,iog%nys-nbgp:iog%nyn-nbgp,i-nbgp)
1952	ENDDO
1953
1954	ELSE
1955	!
1956	!-- Prepare output of 3d-REAL-array without ghost layers
1957	DO i = nxl, nxr
1958	array_3d(:,i,iog%nys:iog%nyn) = data(:,nys:nyn,i)
1959	ENDDO
1960
1961	ENDIF
1962	#if defined( __parallel )
1963	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited number of cores is inactive
1964	IF ( sm_io%iam_io_pe ) THEN
1965	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL, ierr )
1966	CALL MPI_FILE_WRITE_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
1967	ENDIF
1968	CALL sm_io%sm_node_barrier()
1969	#else
1970	CALL posix_lseek( fh, array_position )
1971	CALL posix_write( fh, array_3d, SIZE( array_3d ) )
1972	#endif
1973	!
1974	!-- Type conversion required, otherwise right hand side brackets are calculated assuming 4 byte INT.
1975	!-- Maybe a compiler problem.
1976	array_position = array_position + INT( (nz+2), KIND = rd_offset_kind ) * &
1977	INT( (iog%ny+1), KIND = rd_offset_kind ) * &
1978	INT( (iog%nx+1), KIND = rd_offset_kind ) * wp
1979
1980	END SUBROUTINE wrd_mpi_io_real_3d
1981
1982
1983
1984	!--------------------------------------------------------------------------------------------------!
1985	! Description:
1986	! ------------
1987	!> Write 3d-REAL soil array with MPI-IO.
1988	!> nzb_soil, nzt_soil are located in the module land_surface_model_mod. Since Fortran does not allow
1989	!> cross referencing of module variables, it is required to pass these variables as arguments.
1990	!--------------------------------------------------------------------------------------------------!
1991	SUBROUTINE wrd_mpi_io_real_3d_soil( name, data, nzb_soil, nzt_soil )
1992
1993	IMPLICIT NONE
1994
1995	CHARACTER(LEN=*), INTENT(IN) :: name !<
1996
1997	INTEGER(iwp) :: i !<
1998	INTEGER, INTENT(IN) :: nzb_soil !<
1999	INTEGER, INTENT(IN) :: nzt_soil !<
2000
2001	#if defined( __parallel )
2002	INTEGER, DIMENSION(rd_status_size) :: status !<
2003	#endif
2004
2005	REAL(wp), INTENT(IN), DIMENSION(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) :: data !<
2006
2007
2008	IF ( header_array_index == max_nr_arrays ) THEN
2009	STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
2010	ENDIF
2011
2012	array_names(header_array_index) = name
2013	array_offset(header_array_index) = array_position
2014	header_array_index = header_array_index + 1
2015
2016	#if defined( __parallel )
2017	CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
2018	#endif
2019
2020	IF ( include_total_domain_boundaries) THEN
2021	!
2022	!-- Prepare output of 3d-REAL-array with ghost layers. In the virtual PE grid, the first
2023	!-- dimension is PEs along x, and the second along y. For MPI-IO it is recommended to have the
2024	!-- index order of the array in the same way, i.e. the first dimension should be along x and the
2025	!-- second along y. For this reason, the original PALM data need to be swaped.
2026	DO i = iog%nxl, iog%nxr
2027	array_3d_soil(:,i,iog%nys:iog%nyn) = data(:,iog%nys-nbgp:iog%nyn-nbgp,i-nbgp)
2028	ENDDO
2029
2030	ELSE
2031	!
2032	!-- Prepare output of 3d-REAL-array without ghost layers
2033	DO i = nxl, nxr
2034	array_3d_soil(:,i,iog%nys:iog%nyn) = data(:,nys:nyn,i)
2035	ENDDO
2036
2037	ENDIF
2038	#if defined( __parallel )
2039	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited number of cores is inactive
2040	IF ( sm_io%iam_io_pe ) THEN
2041	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3dsoil, 'native', MPI_INFO_NULL, &
2042	ierr )
2043	CALL MPI_FILE_WRITE_ALL( fh, array_3d_soil, SIZE( array_3d_soil ), MPI_REAL, status, ierr )
2044	ENDIF
2045	CALL sm_io%sm_node_barrier()
2046	#else
2047	CALL posix_lseek( fh, array_position )
2048	CALL posix_write( fh, array_3d_soil, SIZE( array_3d_soil ) )
2049	#endif
2050	!
2051	!-- Type conversion required, otherwise right hand side brackets are calculated assuming 4 byte INT.
2052	!-- Maybe a compiler problem.
2053	array_position = array_position + INT( (nzt_soil-nzb_soil+1), KIND = rd_offset_kind ) * &
2054	INT( (iog%ny+1), KIND = rd_offset_kind ) * &
2055	INT( (iog%nx+1), KIND = rd_offset_kind ) * wp
2056
2057	END SUBROUTINE wrd_mpi_io_real_3d_soil
2058
2059
2060
2061	!--------------------------------------------------------------------------------------------------!
2062	! Description:
2063	! ------------
2064	!> Write CHARATCTER with MPI-IO
2065	!--------------------------------------------------------------------------------------------------!
2066	SUBROUTINE wrd_mpi_io_char( name, text )
2067
2068	IMPLICIT NONE
2069
2070	CHARACTER(LEN=128) :: lo_line !<
2071	CHARACTER(LEN=*), INTENT(IN) :: name !<
2072	CHARACTER(LEN=*), INTENT(IN) :: text !<
2073
2074
2075	IF ( header_char_index == max_nr_char ) THEN
2076	STOP '+++ maximum number of CHARACTER entries in restart file header exceeded'
2077	ENDIF
2078
2079	lo_line = name
2080	lo_line(33:) = text
2081	text_lines(header_char_index) = lo_line
2082	header_char_index = header_char_index + 1
2083
2084	END SUBROUTINE wrd_mpi_io_char
2085
2086
2087
2088	!--------------------------------------------------------------------------------------------------!
2089	! Description:
2090	! ------------
2091	!> Write LOGICAL with MPI-IO
2092	!--------------------------------------------------------------------------------------------------!
2093	SUBROUTINE wrd_mpi_io_logical( name, value )
2094
2095	IMPLICIT NONE
2096
2097	CHARACTER(LEN=*), INTENT(IN) :: name !<
2098
2099	INTEGER(iwp) :: logical_as_integer !<
2100
2101	LOGICAL, INTENT(IN) :: value !<
2102
2103
2104	IF ( value ) THEN
2105	logical_as_integer = 1
2106	ELSE
2107	logical_as_integer = 0
2108	ENDIF
2109
2110	CALL wrd_mpi_io_int( name, logical_as_integer )
2111
2112	END SUBROUTINE wrd_mpi_io_logical
2113
2114
2115
2116	!--------------------------------------------------------------------------------------------------!
2117	! Description:
2118	! ------------
2119	!> Read 1d-REAL global array with MPI-IO.
2120	!> Array contains identical data on all PEs.
2121	!--------------------------------------------------------------------------------------------------!
2122	SUBROUTINE rrd_mpi_io_global_array_real_1d( name, data )
2123
2124	IMPLICIT NONE
2125
2126	CHARACTER(LEN=*), INTENT(IN) :: name !<
2127
2128	INTEGER(iwp) :: i !<
2129	INTEGER(KIND=rd_offset_kind) :: offset !<
2130
2131	#if defined( __parallel )
2132	INTEGER, DIMENSION(rd_status_size) :: status !<
2133	#endif
2134
2135	LOGICAL :: found !<
2136
2137	REAL(KIND=wp), INTENT(INOUT), DIMENSION(:) :: data !<
2138
2139
2140	offset = 0
2141	found = .FALSE.
2142
2143	DO i = 1, tgh%nr_arrays
2144	IF ( TRIM(array_names(i)) == TRIM( name ) ) THEN
2145	array_position = array_offset(i)
2146	found = .TRUE.
2147	EXIT
2148	ENDIF
2149	ENDDO
2150
2151
2152	IF ( found ) THEN
2153
2154	!
2155	!-- Set default view
2156	#if defined( __parallel )
2157	IF ( cyclic_fill_mode ) THEN !kk This may be the general solution for all cases
2158	IF ( pe_active_for_read ) THEN
2159	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
2160	CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
2161	CALL MPI_FILE_READ_ALL( fh, data, SIZE( data ), MPI_REAL, status, ierr )
2162	ENDIF
2163	CALL MPI_BCAST( data, SIZE( data ), MPI_REAL, 0, comm2d, ierr )
2164	ELSE
2165	IF ( sm_io%iam_io_pe ) THEN
2166	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
2167	CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
2168	CALL MPI_FILE_READ_ALL( fh, data, SIZE( data ), MPI_REAL, status, ierr )
2169	ENDIF
2170	IF ( sm_io%is_sm_active() ) THEN
2171	CALL MPI_BCAST( data, SIZE( data ), MPI_REAL, 0, sm_io%comm_shared, ierr )
2172	ENDIF
2173	ENDIF
2174	#else
2175	CALL posix_lseek( fh, array_position )
2176	CALL posix_read( fh, data, SIZE( data ) )
2177	#endif
2178
2179	ELSE
2180
2181	message_string = '1d/2d/3d/4d-REAL global array "' // TRIM( name ) // '" not found in ' // &
2182	'restart file'
2183	CALL message( 'rrd_mpi_io_global_array_real_1d', 'PA0722', 3, 2, 0, 6, 0 )
2184
2185	ENDIF
2186
2187	END SUBROUTINE rrd_mpi_io_global_array_real_1d
2188
2189
2190
2191	!--------------------------------------------------------------------------------------------------!
2192	! Description:
2193	! ------------
2194	!> Read 2d-REAL global array with MPI-IO.
2195	!> Array contains identical data on all PEs.
2196	!--------------------------------------------------------------------------------------------------!
2197	SUBROUTINE rrd_mpi_io_global_array_real_2d( name, data )
2198
2199	IMPLICIT NONE
2200
2201	CHARACTER(LEN=*), INTENT(IN) :: name !<
2202
2203	INTEGER, DIMENSION(1) :: bufshape !<
2204
2205	REAL(KIND=wp), INTENT(IN), DIMENSION(:,:), TARGET :: data !<
2206	REAL(KIND=wp), POINTER, DIMENSION(:) :: buf !<
2207
2208	TYPE(C_PTR) :: c_data !<
2209
2210
2211	c_data = C_LOC( data )
2212	bufshape(1) = SIZE( data )
2213	CALL C_F_POINTER( c_data, buf, bufshape )
2214
2215	CALL rrd_mpi_io_global_array_real_1d( name, buf )
2216
2217	END SUBROUTINE rrd_mpi_io_global_array_real_2d
2218
2219
2220
2221	!--------------------------------------------------------------------------------------------------!
2222	! Description:
2223	! ------------
2224	!> Read 3d-REAL global array with MPI-IO.
2225	!> Array contains identical data on all PEs.
2226	!--------------------------------------------------------------------------------------------------!
2227	SUBROUTINE rrd_mpi_io_global_array_real_3d( name, data )
2228
2229	IMPLICIT NONE
2230
2231	CHARACTER(LEN=*), INTENT(IN) :: name !<
2232
2233	INTEGER, DIMENSION(1) :: bufshape !<
2234
2235	REAL(KIND=wp), INTENT(IN), DIMENSION(:,:,:), TARGET :: data !<
2236	REAL(KIND=wp), POINTER, DIMENSION(:) :: buf !<
2237
2238	TYPE(C_PTR) :: c_data !<
2239
2240
2241	c_data = C_LOC( data )
2242	bufshape(1) = SIZE( data )
2243	CALL C_F_POINTER( c_data, buf, bufshape )
2244
2245	CALL rrd_mpi_io_global_array_real_1d( name, buf )
2246
2247	END SUBROUTINE rrd_mpi_io_global_array_real_3d
2248
2249
2250
2251	!--------------------------------------------------------------------------------------------------!
2252	! Description:
2253	! ------------
2254	!> Read 4d-REAL global array with MPI-IO.
2255	!> Array contains identical data on all PEs.
2256	!--------------------------------------------------------------------------------------------------!
2257	SUBROUTINE rrd_mpi_io_global_array_real_4d( name, data )
2258
2259	IMPLICIT NONE
2260
2261	CHARACTER(LEN=*), INTENT(IN) :: name !<
2262
2263	INTEGER, DIMENSION(1) :: bufshape !<
2264
2265	REAL(KIND=wp), INTENT(IN), DIMENSION(:,:,:,:), TARGET :: data !<
2266	REAL(KIND=wp), POINTER, DIMENSION(:) :: buf !<
2267
2268	TYPE(C_PTR) :: c_data !<
2269
2270
2271	c_data = C_LOC( data )
2272	bufshape(1) = SIZE( data)
2273	CALL C_F_POINTER( c_data, buf, bufshape )
2274
2275	CALL rrd_mpi_io_global_array_real_1d( name, buf )
2276
2277	END SUBROUTINE rrd_mpi_io_global_array_real_4d
2278
2279
2280
2281	!--------------------------------------------------------------------------------------------------!
2282	! Description:
2283	! ------------
2284	!> Read 1d-INTEGER global array with MPI-IO.
2285	!> Array contains identical data on all PEs.
2286	!--------------------------------------------------------------------------------------------------!
2287	SUBROUTINE rrd_mpi_io_global_array_int_1d( name, data )
2288
2289	IMPLICIT NONE
2290
2291	CHARACTER(LEN=*), INTENT(IN) :: name !<
2292
2293	INTEGER(iwp) :: i !<
2294	INTEGER(KIND=rd_offset_kind) :: offset !<
2295
2296	#if defined( __parallel )
2297	INTEGER, DIMENSION(rd_status_size) :: status !<
2298	#endif
2299	INTEGER(KIND=iwp), INTENT(INOUT), DIMENSION(:) :: data !<
2300
2301	LOGICAL :: found !<
2302
2303
2304	offset = 0
2305	found = .FALSE.
2306
2307	DO i = 1, tgh%nr_arrays
2308	IF ( TRIM( array_names(i) ) == TRIM( name ) ) THEN
2309	array_position = array_offset(i)
2310	found = .TRUE.
2311	EXIT
2312	ENDIF
2313	ENDDO
2314
2315	IF ( found ) THEN
2316	!
2317	!-- Set default view
2318	#if defined( __parallel )
2319	IF ( cyclic_fill_mode ) THEN !kk This may be the general solution for all cases
2320	IF ( pe_active_for_read ) THEN
2321	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
2322	CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
2323	CALL MPI_FILE_READ_ALL( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
2324	ENDIF
2325	CALL MPI_BCAST( data, SIZE( data ), MPI_REAL, 0, comm2d, ierr )
2326	ELSE
2327	IF ( sm_io%iam_io_pe ) THEN
2328	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
2329	CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
2330	CALL MPI_FILE_READ_ALL( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
2331	ENDIF
2332	IF ( sm_io%is_sm_active() ) THEN
2333	CALL MPI_BCAST( data, SIZE( data ), MPI_INTEGER, 0, sm_io%comm_shared, ierr )
2334	ENDIF
2335	ENDIF
2336	#else
2337	CALL posix_lseek( fh, array_position )
2338	CALL posix_read( fh, data, SIZE( data ) )
2339	#endif
2340	ELSE
2341
2342	message_string = '1d-INTEGER global array "' // TRIM( name ) // '" not found in ' // &
2343	'restart file'
2344	CALL message( 'rrd_mpi_io_global_array_int_1d', 'PA0722', 3, 2, 0, 6, 0 )
2345
2346	ENDIF
2347
2348	END SUBROUTINE rrd_mpi_io_global_array_int_1d
2349
2350
2351
2352	!--------------------------------------------------------------------------------------------------!
2353	! Description:
2354	! ------------
2355	!> Write 1d-REAL global array with MPI-IO.
2356	!> Array contains identical data on all PEs.
2357	!--------------------------------------------------------------------------------------------------!
2358	SUBROUTINE wrd_mpi_io_global_array_real_1d( name, data )
2359
2360	IMPLICIT NONE
2361
2362	CHARACTER(LEN=*), INTENT(IN) :: name !<
2363
2364	INTEGER(KIND=rd_offset_kind) :: offset !<
2365
2366	#if defined( __parallel )
2367	INTEGER, DIMENSION(rd_status_size) :: status !<
2368	#endif
2369
2370	REAL(KIND=wp), INTENT(IN), DIMENSION(:) :: data !<
2371
2372
2373	offset = 0
2374
2375	IF ( header_array_index == max_nr_arrays ) THEN
2376	STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
2377	ENDIF
2378
2379	array_names(header_array_index) = name
2380	array_offset(header_array_index) = array_position
2381	header_array_index = header_array_index + 1
2382
2383	!
2384	!-- Set default view
2385	#if defined( __parallel )
2386	IF ( sm_io%iam_io_pe ) THEN
2387	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
2388	ENDIF
2389	!
2390	!-- Only PE 0 writes replicated data
2391	IF ( myid == 0 ) THEN
2392	CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
2393	CALL MPI_FILE_WRITE( fh, data, SIZE( data ), MPI_REAL, status, ierr )
2394	ENDIF
2395	#else
2396	CALL posix_lseek( fh, array_position )
2397	CALL posix_write( fh, data, SIZE( data ) )
2398	#endif
2399	array_position = array_position + SIZE( data ) * wp
2400
2401	END SUBROUTINE wrd_mpi_io_global_array_real_1d
2402
2403
2404
2405	!--------------------------------------------------------------------------------------------------!
2406	! Description:
2407	! ------------
2408	!> Write 2d-REAL global array with MPI-IO.
2409	!> Array contains identical data on all PEs.
2410	!--------------------------------------------------------------------------------------------------!
2411	SUBROUTINE wrd_mpi_io_global_array_real_2d( name, data )
2412
2413	IMPLICIT NONE
2414
2415	CHARACTER(LEN=*), INTENT(IN) :: name !<
2416
2417	INTEGER, DIMENSION(1) :: bufshape !<
2418
2419	REAL(KIND=wp), POINTER, DIMENSION(:) :: buf !<
2420	REAL(KIND=wp), INTENT(IN), DIMENSION(:,:), TARGET :: data !<
2421
2422	TYPE(C_PTR) :: c_data !<
2423
2424
2425	c_data = C_LOC( data )
2426	bufshape(1) = SIZE( data)
2427	CALL C_F_POINTER( c_data, buf, bufshape )
2428
2429	CALL wrd_mpi_io_global_array_real_1d( name, buf )
2430
2431	END SUBROUTINE wrd_mpi_io_global_array_real_2d
2432
2433
2434
2435	!--------------------------------------------------------------------------------------------------!
2436	! Description:
2437	! ------------
2438	!> Write 3d-REAL global array with MPI-IO.
2439	!> Array contains identical data on all PEs.
2440	!--------------------------------------------------------------------------------------------------!
2441	SUBROUTINE wrd_mpi_io_global_array_real_3d( name, data )
2442
2443	IMPLICIT NONE
2444
2445	CHARACTER(LEN=*), INTENT(IN) :: name !<
2446
2447	INTEGER, DIMENSION(1) :: bufshape !<
2448
2449	REAL(KIND=wp), POINTER, DIMENSION(:) :: buf !<
2450	REAL(KIND=wp), INTENT(IN), DIMENSION(:,:,:), TARGET :: data !<
2451
2452	TYPE(C_PTR) :: c_data !<
2453
2454
2455	c_data = C_LOC( data )
2456	bufshape(1) = SIZE( data )
2457	CALL C_F_POINTER( c_data, buf, bufshape )
2458
2459	CALL wrd_mpi_io_global_array_real_1d( name, buf )
2460
2461	END SUBROUTINE wrd_mpi_io_global_array_real_3d
2462
2463
2464
2465	!--------------------------------------------------------------------------------------------------!
2466	! Description:
2467	! ------------
2468	!> Write 4d-REAL global array with MPI-IO.
2469	!> Array contains identical data on all PEs.
2470	!--------------------------------------------------------------------------------------------------!
2471	SUBROUTINE wrd_mpi_io_global_array_real_4d( name, data )
2472
2473	IMPLICIT NONE
2474
2475	CHARACTER(LEN=*), INTENT(IN) :: name !<
2476
2477	INTEGER, DIMENSION(1) :: bufshape !<
2478
2479	REAL(KIND=wp), POINTER, DIMENSION(:) :: buf !<
2480	REAL(KIND=wp), INTENT(IN), DIMENSION(:,:,:,:), TARGET :: data !<
2481
2482	TYPE(C_PTR) :: c_data !<
2483
2484
2485	c_data = C_LOC( data )
2486	bufshape(1) = SIZE( data)
2487	CALL C_F_POINTER( c_data, buf, bufshape )
2488
2489	CALL wrd_mpi_io_global_array_real_1d( name, buf )
2490
2491	END SUBROUTINE wrd_mpi_io_global_array_real_4d
2492
2493
2494
2495	!--------------------------------------------------------------------------------------------------!
2496	! Description:
2497	! ------------
2498	!> Write 1d-INTEGER global array with MPI-IO.
2499	!> Array contains identical data on all PEs.
2500	!--------------------------------------------------------------------------------------------------!
2501	SUBROUTINE wrd_mpi_io_global_array_int_1d( name, data )
2502
2503	IMPLICIT NONE
2504
2505	CHARACTER(LEN=*), INTENT(IN) :: name !<
2506
2507	INTEGER(KIND=rd_offset_kind) :: offset !<
2508
2509	INTEGER(KIND=iwp), INTENT(IN), DIMENSION(:) :: data !<
2510	#if defined( __parallel )
2511	INTEGER, DIMENSION(rd_status_size) :: status !<
2512	#endif
2513
2514	IF ( header_array_index == max_nr_arrays ) THEN
2515	STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
2516	ENDIF
2517
2518	offset = 0
2519	array_names(header_array_index) = name
2520	array_offset(header_array_index) = array_position
2521	header_array_index = header_array_index + 1
2522
2523	!
2524	!-- Set default view
2525	#if defined( __parallel )
2526	IF ( sm_io%iam_io_pe ) THEN
2527	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
2528	ENDIF
2529	!
2530	!-- Only PE 0 writes replicated data
2531	IF ( myid == 0 ) THEN
2532	CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
2533	CALL MPI_FILE_WRITE( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
2534	ENDIF
2535	#else
2536	CALL posix_lseek( fh, array_position )
2537	CALL posix_write( fh, data, SIZE( data ) )
2538	#endif
2539	array_position = array_position + SIZE( data ) * 4
2540
2541	END SUBROUTINE wrd_mpi_io_global_array_int_1d
2542
2543
2544
2545	!--------------------------------------------------------------------------------------------------!
2546	! Description:
2547	! ------------
2548	!> Read 1d-REAL surface data array with MPI-IO.
2549	!--------------------------------------------------------------------------------------------------!
2550	SUBROUTINE rrd_mpi_io_surface( name, data, first_index )
2551
2552	IMPLICIT NONE
2553
2554	CHARACTER(LEN=*), INTENT(IN) :: name !<
2555
2556	INTEGER(KIND=rd_offset_kind) :: disp !< displacement of actual indices
2557	INTEGER(KIND=rd_offset_kind) :: disp_f !< displacement in file
2558	INTEGER(KIND=rd_offset_kind) :: disp_n !< displacement of next column
2559	INTEGER(iwp), OPTIONAL :: first_index !<
2560
2561	INTEGER(iwp) :: i !<
2562	INTEGER(iwp) :: i_f !<
2563	INTEGER(iwp) :: j !<
2564	INTEGER(iwp) :: j_f !<
2565	INTEGER(iwp) :: lo_first_index !<
2566	INTEGER(iwp) :: nr_bytes !<
2567	INTEGER(iwp) :: nr_bytes_f !<
2568	INTEGER(iwp) :: nr_words !<
2569	#if defined( __parallel )
2570	INTEGER, DIMENSION(rd_status_size) :: status !<
2571	#endif
2572
2573	LOGICAL :: found !<
2574
2575	REAL(wp), INTENT(OUT), DIMENSION(:) :: data !<
2576
2577
2578	found = .FALSE.
2579	lo_first_index = 1
2580
2581	IF ( PRESENT( first_index ) ) THEN
2582	lo_first_index = first_index
2583	ENDIF
2584
2585	DO i = 1, tgh%nr_arrays
2586	IF ( TRIM( array_names(i) ) == TRIM( name ) ) THEN
2587	array_position = array_offset(i) + ( lo_first_index - 1 ) * &
2588	total_number_of_surface_values * wp
2589	found = .TRUE.
2590	EXIT
2591	ENDIF
2592	ENDDO
2593
2594	disp = -1
2595	disp_f = -1
2596	disp_n = -1
2597	IF ( found ) THEN
2598
2599	IF ( cyclic_fill_mode ) THEN
2600
2601	CALL rrd_mpi_io_surface_cyclic_fill
2602
2603	ELSE
2604
2605	IF ( MAXVAL( m_global_start ) == -1 ) RETURN ! Nothing to do on this PE
2606	DO i = nxl, nxr
2607	DO j = nys, nyn
2608
2609	IF ( m_global_start(j,i) > 0 ) THEN
2610	disp = array_position+(m_global_start(j,i)-1) * wp
2611	nr_words = m_end_index(j,i)-m_start_index(j,i)+1
2612	nr_bytes = nr_words * wp
2613	ENDIF
2614	IF ( disp >= 0 .AND. disp_f == -1 ) THEN ! First entry
2615	disp_f = disp
2616	nr_bytes_f = 0
2617	i_f = i
2618	j_f = j
2619	ENDIF
2620	IF ( j == nyn .AND. i == nxr ) THEN ! Last entry
2621	disp_n = -1
2622	IF ( nr_bytes > 0 ) THEN
2623	nr_bytes_f = nr_bytes_f+nr_bytes
2624	ENDIF
2625	ELSEIF ( j == nyn ) THEN ! Next x
2626	IF ( m_global_start(nys,i+1) > 0 .AND. disp > 0 ) THEN
2627	disp_n = array_position + ( m_global_start(nys,i+1) - 1 ) * wp
2628	ELSE
2629	CYCLE
2630	ENDIF
2631	ELSE
2632	IF ( m_global_start(j+1,i) > 0 .AND. disp > 0 ) THEN
2633	disp_n = array_position + ( m_global_start(j+1,i) - 1 ) * wp
2634	ELSE
2635	CYCLE
2636	ENDIF
2637	ENDIF
2638
2639
2640	IF ( disp + nr_bytes == disp_n ) THEN ! Contiguous block
2641	nr_bytes_f = nr_bytes_f + nr_bytes
2642	ELSE ! Read
2643	#if defined( __parallel )
2644	CALL MPI_FILE_SEEK( fhs, disp_f, MPI_SEEK_SET, ierr )
2645	nr_words = nr_bytes_f / wp
2646	CALL MPI_FILE_READ( fhs, data(m_start_index(j_f,i_f)), nr_words, MPI_REAL, status, &
2647	ierr )
2648	#else
2649	CALL posix_lseek( fh, disp_f )
2650	CALL posix_read( fh, data(m_start_index(j_f,i_f):), nr_bytes_f )
2651	#endif
2652	disp_f = disp
2653	nr_bytes_f = nr_bytes
2654	i_f = i
2655	j_f = j
2656	ENDIF
2657
2658	ENDDO
2659	ENDDO
2660	ENDIF
2661
2662
2663	ELSE
2664
2665	message_string = 'surface array "' // TRIM( name ) // '" not found in restart file'
2666	CALL message( 'rrd_mpi_io_global_array_int_1d', 'PA0722', 3, 2, 0, 6, 0 )
2667
2668	ENDIF
2669
2670	! IF ( lo_first_index == 1 ) THEN
2671	! IF ( debug_level >= 2 .AND. nr_val > 0 ) WRITE(9,*) 'r_surf_1 ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
2672	! ELSE
2673	! IF ( debug_level >= 2 .AND. nr_val > 0 ) WRITE(9,*) 'r_surf_next ', TRIM( name ), ' ', &
2674	! lo_first_index,nr_val, SUM( data(1:nr_val) )
2675	! ENDIF
2676
2677
2678	CONTAINS
2679
2680	SUBROUTINE rrd_mpi_io_surface_cyclic_fill
2681
2682	IMPLICIT NONE
2683
2684	INTEGER(iwp) :: i !<
2685	INTEGER(iwp) :: ie !<
2686	INTEGER(iwp) :: ierr !<
2687	INTEGER(iwp) :: is !<
2688	INTEGER(iwp) :: i_remote !<
2689	INTEGER(iwp) :: j !<
2690	INTEGER(iwp) :: je !<
2691	INTEGER(iwp) :: js !<
2692	INTEGER(iwp) :: j_remote !<
2693	INTEGER(iwp) :: nval !<
2694	INTEGER(iwp) :: rem_pe !<
2695
2696	INTEGER(KIND=MPI_ADDRESS_KIND) :: rem_offs !<
2697
2698	LOGICAL :: write_done !<
2699
2700
2701	!
2702	!-- In the current version, there is only 1 value per grid cell allowed.
2703	!-- In this special case, the cyclical repetition can be done with the same method as for 2d-real
2704	!-- array.
2705	CALL prerun_grid%activate_grid_from_this_class()
2706
2707	IF ( pe_active_for_read ) THEN
2708	rmabuf_2d = -1.0
2709	DO i = nxl, nxr
2710	DO j = nys, nyn
2711
2712	IF ( m_global_start(j,i) > 0 ) THEN
2713	disp = array_position+(m_global_start(j,i)-1) * wp
2714	nr_words = m_end_index(j,i)-m_start_index(j,i)+1
2715	nr_bytes = nr_words * wp
2716	ENDIF
2717	IF ( disp >= 0 .AND. disp_f == -1 ) THEN ! First entry
2718	disp_f = disp
2719	nr_bytes_f = 0
2720	write_done = .TRUE.
2721	ENDIF
2722	IF( write_done ) THEN
2723	i_f = i
2724	j_f = j
2725	write_done = .FALSE.
2726	ENDIF
2727
2728	IF ( j == nyn .AND. i == nxr ) THEN ! Last entry
2729	disp_n = -1
2730	IF ( nr_bytes > 0 ) THEN
2731	nr_bytes_f = nr_bytes_f+nr_bytes
2732	ENDIF
2733	ELSEIF ( j == nyn ) THEN ! Next x
2734	IF ( m_global_start(nys,i+1) > 0 .AND. disp > 0 ) THEN
2735	disp_n = array_position + ( m_global_start(nys,i+1) - 1 ) * wp
2736	ELSE
2737	CYCLE
2738	ENDIF
2739	ELSE
2740	IF ( m_global_start(j+1,i) > 0 .AND. disp > 0 ) THEN
2741	disp_n = array_position + ( m_global_start(j+1,i) - 1 ) * wp
2742	ELSE
2743	CYCLE
2744	ENDIF
2745	ENDIF
2746
2747
2748	IF ( disp + nr_bytes == disp_n ) THEN ! Contiguous block
2749	nr_bytes_f = nr_bytes_f + nr_bytes
2750	ELSE ! Read
2751	#if defined( __parallel )
2752	CALL MPI_FILE_SEEK( fhs, disp_f, MPI_SEEK_SET, ierr )
2753	nr_words = nr_bytes_f / wp
2754	CALL MPI_FILE_READ( fhs, rmabuf_2d(j_f,i_f), nr_words, MPI_REAL, status, ierr )
2755	#else
2756	CALL posix_lseek( fh, disp_f )
2757	CALL posix_read( fh, rmabuf_2d(j_f,i_f), nr_bytes_f )
2758	#endif
2759
2760	disp_f = disp
2761	nr_bytes_f = nr_bytes
2762	write_done = .TRUE.
2763	ENDIF
2764
2765	ENDDO
2766	ENDDO
2767
2768	ENDIF
2769
2770	CALL mainrun_grid%activate_grid_from_this_class()
2771
2772	#if defined( __parallel )
2773	!
2774	!-- Close RMA window to allow remote access
2775	CALL MPI_WIN_FENCE( 0, rmawin_2d, ierr )
2776	#endif
2777
2778	IF ( .NOT. pe_active_for_read ) THEN
2779
2780	is = nxl
2781	ie = nxr
2782	js = nys
2783	je = nyn
2784
2785	ELSE
2786
2787	is = nxl
2788	ie = nxr
2789	js = prerun_grid%nys+1
2790	je = nyn
2791
2792	DO i = is, ie
2793	DO j = js, je
2794	i_remote = MOD(i,nx_on_file+1)
2795	j_remote = MOD(j,ny_on_file+1)
2796	rem_pe = remote_pe(i_remote,j_remote)
2797	rem_offs = rma_offset(i_remote,j_remote)
2798	nval = 1
2799
2800	#if defined( __parallel )
2801	IF ( rem_pe /= myid ) THEN
2802	CALL MPI_GET( data(m_start_index(j,i)), nval, MPI_REAL, rem_pe, rem_offs, nval, &
2803	MPI_REAL, rmawin_2d, ierr)
2804	ELSE
2805	data(m_start_index(j,i)) = rmabuf_2d(j_remote,i_remote)
2806	ENDIF
2807	#else
2808	data(m_start_index(j,i)) = array_2d(i_remote,j_remote)
2809	#endif
2810	ENDDO
2811	ENDDO
2812	is = prerun_grid%nxr+1
2813	ie = nxr
2814	js = nys
2815	je = nyn
2816
2817	ENDIF
2818
2819	DO i = is, ie
2820	DO j = js, je
2821	i_remote = MOD(i,nx_on_file+1)
2822	j_remote = MOD(j,ny_on_file+1)
2823	rem_pe = remote_pe(i_remote,j_remote)
2824	rem_offs = rma_offset(i_remote,j_remote)
2825	nval = 1
2826
2827	#if defined( __parallel )
2828	IF ( rem_pe /= myid ) THEN
2829	CALL MPI_GET( data(m_start_index(j,i)), nval, MPI_REAL, rem_pe, rem_offs, nval, &
2830	MPI_REAL, rmawin_2d, ierr)
2831	ELSE
2832	data(m_start_index(j,i)) = rmabuf_2d(j_remote,i_remote)
2833	ENDIF
2834	#else
2835	data(m_tart_index(j,i)) = array_2d(i_remote,j_remote)
2836	#endif
2837	ENDDO
2838	ENDDO
2839
2840	#if defined( __parallel )
2841	!
2842	!-- Reopen RMA window to allow filling
2843	CALL MPI_WIN_FENCE( 0, rmawin_2d, ierr )
2844	#endif
2845
2846	END SUBROUTINE rrd_mpi_io_surface_cyclic_fill
2847
2848	END SUBROUTINE rrd_mpi_io_surface
2849
2850
2851
2852	!--------------------------------------------------------------------------------------------------!
2853	! Description:
2854	! ------------
2855	!> Read 2d-REAL surface data array with MPI-IO.
2856	!> These consist of multiple 1d-REAL surface data arrays.
2857	!--------------------------------------------------------------------------------------------------!
2858	SUBROUTINE rrd_mpi_io_surface_2d( name, data )
2859
2860	IMPLICIT NONE
2861
2862	CHARACTER(LEN=*), INTENT(IN) :: name !<
2863
2864	INTEGER(iwp) :: i !<
2865
2866	REAL(wp), INTENT(OUT), DIMENSION(:,:) :: data !<
2867	REAL(wp), DIMENSION(SIZE( data,2)) :: tmp !<
2868
2869
2870	DO i = 1, SIZE( data, 1 )
2871	CALL rrd_mpi_io_surface( name, tmp, first_index = i )
2872	data(i,:) = tmp
2873	ENDDO
2874
2875	!
2876	!-- Comment from Klaus Ketelsen (September 2018)
2877	!-- The intention of the following loop was to let the compiler do the copying on return.
2878	!-- In my understanding it is standard conform to pass the second dimension to a 1d-array inside a
2879	!-- subroutine and the compiler is responsible to generate code for copying. Acually this works fine
2880	!-- for INTENT(IN) variables (wrd_mpi_io_surface_2d). For INTENT(OUT) like in this case the code
2881	!-- works on pgi compiler. But both, the Intel 16 and the Cray compiler show wrong answers using
2882	!-- this loop. That is the reason why the above auxiliary array tmp was introduced.
2883	! DO i = 1, SIZE( data,1)
2884	! CALL rrd_mpi_io_surface( name, data(i,:), first_index = i )
2885	! ENDDO
2886
2887	END SUBROUTINE rrd_mpi_io_surface_2d
2888
2889
2890
2891	!--------------------------------------------------------------------------------------------------!
2892	! Description:
2893	! ------------
2894	!> Write 1d-REAL surface data array with MPI-IO.
2895	!--------------------------------------------------------------------------------------------------!
2896	SUBROUTINE wrd_mpi_io_surface( name, data, first_index )
2897
2898	IMPLICIT NONE
2899
2900	CHARACTER(LEN=*), INTENT(IN) :: name !<
2901
2902	#if defined( __parallel )
2903	INTEGER(KIND=rd_offset_kind) :: disp !<
2904	#endif
2905	INTEGER(iwp), OPTIONAL :: first_index !<
2906	#if defined( __parallel )
2907	INTEGER(iwp) :: i !<
2908	#endif
2909	INTEGER(iwp) :: lo_first_index !<
2910	INTEGER(KIND=rd_offset_kind) :: offset !<
2911
2912	#if defined( __parallel )
2913	INTEGER, DIMENSION(rd_status_size) :: status !<
2914	#endif
2915
2916	REAL(wp), INTENT(IN), DIMENSION(:), TARGET :: data !<
2917
2918
2919	offset = 0
2920	lo_first_index = 1
2921
2922	IF ( PRESENT( first_index ) ) THEN
2923	lo_first_index = first_index
2924	ENDIF
2925	!
2926	!-- In case of 2d-data, name is written only once
2927	IF ( lo_first_index == 1 ) THEN
2928
2929	IF ( header_array_index == max_nr_arrays ) THEN
2930	STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
2931	ENDIF
2932
2933	array_names(header_array_index) = name
2934	array_offset(header_array_index) = array_position
2935	header_array_index = header_array_index + 1
2936
2937	ENDIF
2938
2939	#if defined( __parallel )
2940	IF ( sm_io%is_sm_active() ) THEN
2941	DO i = 1, nr_val
2942	array_1d(i+local_start) = data(i)
2943	ENDDO
2944	ELSE
2945	! array_1d => data !kk Did not work in all cases why???
2946	ALLOCATE( array_1d( SIZE( data ) ) )
2947	array_1d = data
2948	ENDIF
2949
2950	CALL sm_io%sm_node_barrier() ! Has no effect if I/O on limited number of cores is inactive
2951	IF ( sm_io%iam_io_pe ) THEN
2952	IF ( all_pes_write ) THEN
2953	CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_surf, 'native', MPI_INFO_NULL, &
2954	ierr )
2955	CALL MPI_FILE_WRITE_ALL( fh, array_1d, nr_iope, MPI_REAL, status, ierr )
2956	ELSE
2957	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
2958	IF ( nr_val > 0 ) THEN
2959	disp = array_position + 8 * ( glo_start - 1 )
2960	CALL MPI_FILE_SEEK( fh, disp, MPI_SEEK_SET, ierr )
2961	CALL MPI_FILE_WRITE( fh, array_1d, nr_iope, MPI_REAL, status, ierr )
2962	ENDIF
2963	ENDIF
2964	ENDIF
2965	CALL sm_io%sm_node_barrier()
2966	IF( .NOT. sm_io%is_sm_active() ) DEALLOCATE( array_1d )
2967	#else
2968	CALL posix_lseek( fh, array_position )
2969	CALL posix_write( fh, data, nr_val )
2970	#endif
2971	array_position = array_position + total_number_of_surface_values * wp
2972
2973	! IF ( lo_first_index == 1 ) THEN
2974	! IF ( debug_level >= 2 .AND. nr_val > 0 ) WRITE(9,*) 'w_surf_1 ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
2975	! ELSE
2976	! IF ( debug_level >= 2 .AND. nr_val > 0 ) WRITE(9,*) 'w_surf_n ', TRIM( name ), ' ', &
2977	! lo_first_index, nr_val, SUM( data(1:nr_val) )
2978	! ENDIF
2979
2980	END SUBROUTINE wrd_mpi_io_surface
2981
2982
2983
2984	!--------------------------------------------------------------------------------------------------!
2985	! Description:
2986	! ------------
2987	!> Read 2d-REAL surface data array with MPI-IO.
2988	!> This consist of multiple 1d-REAL surface data arrays.
2989	!--------------------------------------------------------------------------------------------------!
2990	SUBROUTINE wrd_mpi_io_surface_2d( name, data )
2991
2992	IMPLICIT NONE
2993
2994	CHARACTER(LEN=*), INTENT(IN) :: name !<
2995
2996	INTEGER(iwp) :: i !<
2997
2998	REAL(wp), INTENT(IN), DIMENSION(:,:) :: data !<
2999
3000
3001	DO i = 1, SIZE( data, 1 )
3002	CALL wrd_mpi_io_surface( name, data(i,:), first_index = i )
3003	ENDDO
3004
3005	END SUBROUTINE wrd_mpi_io_surface_2d
3006
3007
3008
3009	!--------------------------------------------------------------------------------------------------!
3010	! Description:
3011	! ------------
3012	!> Close restart file for MPI-IO
3013	!--------------------------------------------------------------------------------------------------!
3014	SUBROUTINE rd_mpi_io_close
3015
3016	IMPLICIT NONE
3017
3018	INTEGER(iwp) :: gh_size !<
3019	INTEGER(KIND=rd_offset_kind) :: offset !<
3020	#if defined( __parallel )
3021	INTEGER, DIMENSION(rd_status_size) :: status !<
3022	#endif
3023
3024	#if ! defined( __parallel )
3025	TYPE(C_PTR) :: buf_ptr !<
3026	#endif
3027
3028
3029	offset = 0
3030
3031	IF ( wr_flag .AND. sm_io%iam_io_pe ) THEN
3032
3033	tgh%nr_int = header_int_index - 1
3034	tgh%nr_char = header_char_index - 1
3035	tgh%nr_real = header_real_index - 1
3036	tgh%nr_arrays = header_array_index - 1
3037	tgh%total_nx = iog%nx + 1
3038	tgh%total_ny = iog%ny + 1
3039	IF ( include_total_domain_boundaries ) THEN ! Not sure, if LOGICAL interpretation is the same for all compilers,
3040	tgh%i_outer_bound = 1 ! therefore store as INTEGER in general header
3041	ELSE
3042	tgh%i_outer_bound = 0
3043	ENDIF
3044	!
3045	!-- Check for big/little endian format. This check is currently not used, and could be removed
3046	!-- if we can assume little endian as the default on all machines.
3047	CALL rd_mpi_io_check_endian( tgh%endian )
3048
3049	!
3050	!-- Set default view
3051	#if defined( __parallel )
3052	CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
3053	#endif
3054	!
3055	!-- Write header file
3056	gh_size = storage_size(tgh) / 8
3057	IF ( myid == 0 ) THEN ! myid = 0 always performs I/O, even if I/O is limited to some cores
3058	#if defined( __parallel )
3059	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
3060	CALL MPI_FILE_WRITE( fh, tgh, gh_size, MPI_INT, status, ierr )
3061	header_position = header_position + gh_size
3062	!
3063	!-- INTEGER values
3064	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
3065	CALL MPI_FILE_WRITE( fh, int_names, SIZE( int_names ) * 32, MPI_CHAR, status, ierr )
3066	header_position = header_position + SIZE( int_names ) * 32
3067
3068	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
3069	CALL MPI_FILE_WRITE( fh, int_values, SIZE( int_values ), MPI_INT, status, ierr )
3070	header_position = header_position + SIZE( int_values ) * iwp
3071	!
3072	!-- Character entries
3073	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
3074	CALL MPI_FILE_WRITE( fh, text_lines, SIZE( text_lines ) * 128, MPI_CHAR, status, ierr )
3075	header_position = header_position + SIZE( text_lines ) * 128
3076	!
3077	!--- REAL values
3078	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
3079	CALL MPI_FILE_WRITE( fh, real_names, SIZE( real_names ) * 32, MPI_CHAR, status, ierr )
3080	header_position = header_position + SIZE( real_names ) * 32
3081
3082	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
3083	CALL MPI_FILE_WRITE( fh, real_values, SIZE( real_values ), MPI_REAL, status, ierr )
3084	header_position = header_position + SIZE( real_values ) * wp
3085	!
3086	!-- 2d- and 3d- distributed array headers, all replicated array headers
3087	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
3088	CALL MPI_FILE_WRITE( fh, array_names, SIZE( array_names ) * 32, MPI_CHAR, status, ierr )
3089	header_position = header_position + SIZE( array_names ) * 32
3090
3091	CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
3092	CALL MPI_FILE_WRITE( fh, array_offset, SIZE( array_offset ) * MPI_OFFSET_KIND, MPI_BYTE, &
3093	status, ierr ) ! There is no I*8 datatype in Fortran
3094	header_position = header_position + SIZE( array_offset ) * rd_offset_kind
3095	#else
3096	CALL posix_lseek( fh, header_position )
3097	buf_ptr = C_LOC( tgh )
3098	CALL posix_write( fh, buf_ptr, gh_size )
3099	header_position = header_position + gh_size
3100	!
3101	!-- INTEGER values
3102	CALL posix_lseek( fh, header_position )
3103	CALL posix_write( fh, int_names )
3104	header_position = header_position + SIZE( int_names ) * 32
3105
3106	CALL posix_lseek( fh, header_position )
3107	CALL posix_write( fh, int_values, SIZE( int_values ) )
3108	header_position = header_position + SIZE( int_values ) * iwp
3109	!
3110	!-- Character entries
3111	CALL posix_lseek( fh, header_position )
3112	CALL posix_write( fh, text_lines )
3113	header_position = header_position + SIZE( text_lines ) * 128
3114	!
3115	!-- REAL values
3116	CALL posix_lseek( fh, header_position )
3117	CALL posix_write( fh, real_names )
3118	header_position = header_position + SIZE( real_names ) * 32
3119
3120	CALL posix_lseek( fh, header_position )
3121	CALL posix_write( fh, real_values, SIZE( real_values ) )
3122	header_position = header_position + SIZE( real_values ) * wp
3123	!
3124	!-- 2d- and 3d-distributed array headers, all replicated array headers
3125	CALL posix_lseek( fh, header_position )
3126	CALL posix_write( fh, array_names )
3127	header_position = header_position + SIZE( array_names ) * 32
3128
3129	CALL posix_lseek( fh, header_position )
3130	CALL posix_write( fh, array_offset, SIZE( array_offset ) )
3131	header_position = header_position + SIZE( array_offset ) * rd_offset_kind
3132	#endif
3133	IF ( debug_output ) CALL rd_mpi_io_print_header
3134	ENDIF
3135
3136	ENDIF
3137
3138	!
3139	!-- Free file types
3140	CALL rd_mpi_io_free_filetypes
3141
3142	!
3143	!-- Close MPI-IO files
3144	#if defined( __parallel )
3145	!
3146	!-- Restart file has been opened with comm2d
3147	IF ( fhs /= -1 ) THEN
3148	CALL MPI_FILE_CLOSE( fhs, ierr )
3149	ENDIF
3150	#endif
3151	!
3152	!-- Free RMA windows
3153	IF ( cyclic_fill_mode ) THEN
3154	CALL MPI_WIN_FREE( rmawin_2di, ierr )
3155	CALL MPI_WIN_FREE( rmawin_2d, ierr )
3156	CALL MPI_WIN_FREE( rmawin_3d, ierr )
3157	ENDIF
3158
3159	IF (.NOT. pe_active_for_read ) RETURN
3160	!
3161	!-- TODO: better explain the following message
3162	!-- In case on non cyclic read, pe_active_for_read is set .TRUE.
3163	IF ( sm_io%iam_io_pe ) THEN
3164
3165	#if defined( __parallel )
3166	CALL MPI_FILE_CLOSE( fh, ierr )
3167	#else
3168	CALL posix_close( fh )
3169	#endif
3170
3171	ENDIF
3172
3173	restart_file_size = array_position / ( 1024.0_dp * 1024.0_dp )
3174
3175	END SUBROUTINE rd_mpi_io_close
3176
3177
3178
3179	!--------------------------------------------------------------------------------------------------!
3180	! Description:
3181	! ------------
3182	!> This subroutine prepares a filetype and some variables for the I/O of surface data.
3183	!> Whenever a new set of start_index and end_index is used, rd_mpi_io_surface_filetypes has to be
3184	!> called. A main feature of this subroutine is computing the global start indices of the 1d- and
3185	!> 2d- surface arrays.
3186	!> Even if I/O is done by a limited number of cores only, the surface data are read by ALL cores!
3187	!> Reading them by some cores and then distributing the data would result in complicated code
3188	!> which is suspectable for errors and overloads the reading subroutine. Since reading of surface
3189	!> data is not time critical (data size is comparably small), it will be read by all cores.
3190	!--------------------------------------------------------------------------------------------------!
3191	SUBROUTINE rd_mpi_io_surface_filetypes( start_index, end_index, data_to_write, global_start )
3192
3193	IMPLICIT NONE
3194
3195	INTEGER(iwp) :: i !< loop index
3196	INTEGER(iwp) :: j !< loop index
3197	INTEGER(KIND=rd_offset_kind) :: offset !<
3198
3199	INTEGER(iwp), DIMENSION(1) :: dims1 !<
3200	INTEGER(iwp), DIMENSION(1) :: lize1 !<
3201	INTEGER(iwp), DIMENSION(1) :: start1 !<
3202
3203	INTEGER(iwp), DIMENSION(0:numprocs-1) :: all_nr_val !< number of values for all PEs
3204	INTEGER(iwp), DIMENSION(0:numprocs-1) :: lo_nr_val !< local number of values in x and y direction
3205
3206
3207	INTEGER, INTENT(INOUT), DIMENSION(nys:nyn,nxl:nxr) :: end_index !<
3208	INTEGER, INTENT(OUT), DIMENSION(nys:nyn,nxl:nxr) :: global_start !<
3209	INTEGER, INTENT(INOUT), DIMENSION(nys:nyn,nxl:nxr) :: start_index !<
3210
3211	LOGICAL, INTENT(OUT) :: data_to_write !< returns, if surface data have to be written
3212
3213	!
3214	!-- Actions during reading
3215	IF ( rd_flag ) THEN
3216	!
3217	!-- Set start index and end index for the mainrun grid.
3218	!-- ATTENTION: This works only for horizontal surfaces with one vale per grid cell!!!
3219	IF ( cyclic_fill_mode ) THEN
3220	DO i = nxl, nxr
3221	DO j = nys, nyn
3222	start_index (j,i) = (i-nxl) * nny + j - nys + 1
3223	end_index (j,i) = start_index(j,i)
3224	ENDDO
3225	ENDDO
3226	ENDIF
3227
3228	IF ( .NOT. ALLOCATED( m_start_index ) ) ALLOCATE( m_start_index(nys:nyn,nxl:nxr) )
3229	IF ( .NOT. ALLOCATED( m_end_index ) ) ALLOCATE( m_end_index(nys:nyn,nxl:nxr) )
3230	IF ( .NOT. ALLOCATED( m_global_start ) ) ALLOCATE( m_global_start(nys:nyn,nxl:nxr) )
3231	!
3232	!-- Save arrays for later reading
3233	m_start_index = start_index
3234	m_end_index = end_index
3235	m_global_start = global_start
3236	nr_val = MAXVAL( end_index )
3237
3238	ENDIF
3239
3240	IF ( .NOT. pe_active_for_read ) RETURN
3241
3242	IF ( cyclic_fill_mode ) CALL prerun_grid%activate_grid_from_this_class()
3243
3244	offset = 0
3245	lo_nr_val= 0
3246	lo_nr_val(myid) = MAXVAL( end_index )
3247	#if defined( __parallel )
3248	CALL MPI_ALLREDUCE( lo_nr_val, all_nr_val, numprocs, MPI_INTEGER, MPI_SUM, comm2d, ierr )
3249	IF ( ft_surf /= -1 .AND. sm_io%iam_io_pe ) THEN
3250	CALL MPI_TYPE_FREE( ft_surf, ierr ) ! If set, free last surface filetype
3251	ENDIF
3252
3253	IF ( win_surf /= -1 ) THEN
3254	IF ( sm_io%is_sm_active() ) THEN
3255	CALL MPI_WIN_FREE( win_surf, ierr )
3256	ENDIF
3257	win_surf = -1
3258	ENDIF
3259
3260	IF ( sm_io%is_sm_active() .AND. rd_flag ) THEN
3261	IF ( fhs == -1 ) THEN
3262	CALL MPI_FILE_OPEN( comm2d, TRIM( io_file_name ), MPI_MODE_RDONLY, MPI_INFO_NULL, fhs, &
3263	ierr )
3264	ENDIF
3265	ELSE
3266	fhs = fh
3267	ENDIF
3268	#else
3269	all_nr_val(myid) = lo_nr_val(myid)
3270	#endif
3271	nr_val = lo_nr_val(myid)
3272
3273	total_number_of_surface_values = 0
3274	DO i = 0, numprocs-1
3275	IF ( i == myid ) THEN
3276	glo_start = total_number_of_surface_values + 1
3277	ENDIF
3278	total_number_of_surface_values = total_number_of_surface_values + all_nr_val(i)
3279	ENDDO
3280
3281	!
3282	!-- Actions during reading
3283	IF ( rd_flag ) THEN
3284
3285	#if defined( __parallel )
3286	CALL MPI_FILE_SET_VIEW( fhs, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
3287	#endif
3288	ENDIF
3289
3290	IF ( cyclic_fill_mode ) CALL mainrun_grid%activate_grid_from_this_class()
3291
3292	!
3293	!-- Actions during writing
3294	IF ( wr_flag ) THEN
3295	!
3296	!-- Create surface filetype
3297	ft_surf = -1
3298	global_start = start_index + glo_start - 1
3299
3300	WHERE ( end_index < start_index )
3301	global_start = -1
3302	ENDWHERE
3303
3304	#if defined( __parallel )
3305	IF ( sm_io%is_sm_active() ) THEN
3306	IF ( sm_io%iam_io_pe ) THEN
3307	!
3308	!-- Calculate number of values of all PEs of an I/O group
3309	nr_iope = 0
3310	DO i = myid, myid+sm_io%sh_npes-1
3311	nr_iope = nr_iope + all_nr_val(i)
3312	ENDDO
3313	ELSE
3314	local_start = 0
3315	DO i = myid-sm_io%sh_rank, myid-1
3316	local_start = local_start + all_nr_val(i)
3317	ENDDO
3318	ENDIF
3319	!
3320	!-- Get the size of shared memory window on all PEs
3321	CALL MPI_BCAST( nr_iope, 1, MPI_INTEGER, 0, sm_io%comm_shared, ierr )
3322	CALL sm_io%sm_allocate_shared( array_1d, 1, MAX( 1, nr_iope ), win_surf )
3323	ELSE
3324	nr_iope = nr_val
3325	ENDIF
3326	#else
3327	nr_iope = nr_val
3328	#endif
3329
3330	!
3331	!-- Check, if surface data exist on this PE
3332	data_to_write = .TRUE.
3333	IF ( total_number_of_surface_values == 0 ) THEN
3334	data_to_write = .FALSE.
3335	RETURN
3336	ENDIF
3337
3338	IF ( sm_io%iam_io_pe ) THEN
3339
3340	all_pes_write = ( MINVAL( all_nr_val ) > 0 )
3341
3342	IF ( all_pes_write ) THEN
3343	dims1(1) = total_number_of_surface_values
3344	lize1(1) = nr_iope
3345	start1(1) = glo_start-1
3346
3347	#if defined( __parallel )
3348	IF ( total_number_of_surface_values > 0 ) THEN
3349	CALL MPI_TYPE_CREATE_SUBARRAY( 1, dims1, lize1, start1, MPI_ORDER_FORTRAN, &
3350	MPI_REAL, ft_surf, ierr )
3351	CALL MPI_TYPE_COMMIT( ft_surf, ierr )
3352	ENDIF
3353	#endif
3354	ENDIF
3355	ENDIF
3356
3357	ENDIF
3358
3359	END SUBROUTINE rd_mpi_io_surface_filetypes
3360
3361
3362
3363	!--------------------------------------------------------------------------------------------------!
3364	! Description:
3365	! ------------
3366	!> This subroutine creates file types to access 2d-/3d-REAL arrays and 2d-INTEGER arrays
3367	!> distributed in blocks among processes to a single file that contains the global arrays.
3368	!--------------------------------------------------------------------------------------------------!
3369	SUBROUTINE rd_mpi_io_create_filetypes
3370
3371	IMPLICIT NONE
3372
3373	INTEGER, DIMENSION(2) :: dims2 !<
3374	INTEGER, DIMENSION(2) :: lize2 !<
3375	INTEGER, DIMENSION(2) :: start2 !<
3376
3377	INTEGER, DIMENSION(3) :: dims3 !<
3378	INTEGER, DIMENSION(3) :: lize3 !<
3379	INTEGER, DIMENSION(3) :: start3 !<
3380
3381	TYPE(domain_decomposition_grid_features) :: save_io_grid !< temporary variable to store grid settings
3382
3383
3384	IF ( sm_io%is_sm_active() ) THEN
3385	save_io_grid = sm_io%io_grid
3386	ENDIF
3387
3388	IF( .NOT. pe_active_for_read ) RETURN
3389
3390	IF ( cyclic_fill_mode ) CALL prerun_grid%activate_grid_from_this_class()
3391
3392	!
3393	!-- Decision, if storage with or without ghost layers.
3394	!-- Please note that the indexing of the global array always starts at 0, even in Fortran.
3395	!-- Therefore the PE local indices have to be shifted by nbgp in the case with ghost layers.
3396	IF ( include_total_domain_boundaries ) THEN
3397
3398	iog%nxl = nxl + nbgp
3399	iog%nxr = nxr + nbgp
3400	iog%nys = nys + nbgp
3401	iog%nyn = nyn + nbgp
3402	iog%nnx = nnx
3403	iog%nny = nny
3404	iog%nx = nx + 2 * nbgp
3405	iog%ny = ny + 2 * nbgp
3406	IF ( myidx == 0 ) THEN
3407	iog%nxl = iog%nxl - nbgp
3408	iog%nnx = iog%nnx + nbgp
3409	ENDIF
3410	IF ( myidx == npex-1 .OR. npex == -1 ) THEN ! npex == 1 if -D__parallel not set
3411	iog%nxr = iog%nxr + nbgp
3412	iog%nnx = iog%nnx + nbgp
3413	ENDIF
3414	IF ( myidy == 0 ) THEN
3415	iog%nys = iog%nys - nbgp
3416	iog%nny = iog%nny + nbgp
3417	ENDIF
3418	IF ( myidy == npey-1 .OR. npey == -1 ) THEN ! npey == 1 if -D__parallel not set
3419	iog%nyn = iog%nyn + nbgp
3420	iog%nny = iog%nny + nbgp
3421	ENDIF
3422
3423	CALL sm_io%sm_adjust_outer_boundary()
3424
3425	ELSE
3426
3427	iog%nxl = nxl
3428	iog%nxr = nxr
3429	iog%nys = nys
3430	iog%nyn = nyn
3431	iog%nnx = nnx
3432	iog%nny = nny
3433	iog%nx = nx
3434	iog%ny = ny
3435
3436	ENDIF
3437
3438	IF ( sm_io%is_sm_active() ) THEN
3439	#if defined( __parallel )
3440	CALL sm_io%sm_allocate_shared( array_2d, sm_io%io_grid%nxl, sm_io%io_grid%nxr, &
3441	sm_io%io_grid%nys, sm_io%io_grid%nyn, win_2dr )
3442	CALL sm_io%sm_allocate_shared( array_2di, save_io_grid%nxl, save_io_grid%nxr, &
3443	save_io_grid%nys, save_io_grid%nyn, win_2di )
3444	CALL sm_io%sm_allocate_shared( array_3d, nzb, nzt+1, sm_io%io_grid%nxl, sm_io%io_grid%nxr, &
3445	sm_io%io_grid%nys, sm_io%io_grid%nyn, win_3dr )
3446	#endif
3447	ELSE
3448	ALLOCATE( array_2d(iog%nxl:iog%nxr,iog%nys:iog%nyn) )
3449	ALLOCATE( array_2di(nxl:nxr,nys:nyn) )
3450	ALLOCATE( array_3d(nzb:nzt+1,iog%nxl:iog%nxr,iog%nys:iog%nyn) )
3451	sm_io%io_grid = iog
3452	ENDIF
3453
3454	!
3455	!-- Create filetype for 2d-REAL array with ghost layers around the total domain
3456	dims2(1) = iog%nx + 1
3457	dims2(2) = iog%ny + 1
3458
3459	lize2(1) = sm_io%io_grid%nnx
3460	lize2(2) = sm_io%io_grid%nny
3461
3462	start2(1) = sm_io%io_grid%nxl
3463	start2(2) = sm_io%io_grid%nys
3464
3465	#if defined( __parallel )
3466	IF ( sm_io%iam_io_pe ) THEN
3467	CALL MPI_TYPE_CREATE_SUBARRAY( 2, dims2, lize2, start2, MPI_ORDER_FORTRAN, MPI_REAL, &
3468	ft_2d, ierr )
3469	CALL MPI_TYPE_COMMIT( ft_2d, ierr )
3470	ENDIF
3471	#endif
3472	!
3473	!-- Create filetype for 2d-INTEGER array without ghost layers around the total domain
3474	dims2(1) = nx + 1
3475	dims2(2) = ny + 1
3476
3477	IF ( sm_io%is_sm_active() ) THEN
3478
3479	lize2(1) = save_io_grid%nnx
3480	lize2(2) = save_io_grid%nny
3481
3482	start2(1) = save_io_grid%nxl
3483	start2(2) = save_io_grid%nys
3484
3485	ELSE
3486
3487	lize2(1) = nnx
3488	lize2(2) = nny
3489
3490	start2(1) = nxl
3491	start2(2) = nys
3492
3493	ENDIF
3494
3495	#if defined( __parallel )
3496	IF ( sm_io%iam_io_pe ) THEN
3497	CALL MPI_TYPE_CREATE_SUBARRAY( 2, dims2, lize2, start2, MPI_ORDER_FORTRAN, MPI_INTEGER, &
3498	ft_2di_nb, ierr )
3499	CALL MPI_TYPE_COMMIT( ft_2di_nb, ierr )
3500	ENDIF
3501	#endif
3502	!
3503	!-- Create filetype for 3d-REAL array
3504	dims3(1) = nz + 2
3505	dims3(2) = iog%nx + 1
3506	dims3(3) = iog%ny + 1
3507
3508	lize3(1) = dims3(1)
3509	lize3(2) = sm_io%io_grid%nnx
3510	lize3(3) = sm_io%io_grid%nny
3511
3512	start3(1) = nzb
3513	start3(2) = sm_io%io_grid%nxl
3514	start3(3) = sm_io%io_grid%nys
3515
3516	#if defined( __parallel )
3517	IF ( sm_io%iam_io_pe ) THEN
3518	CALL MPI_TYPE_CREATE_SUBARRAY( 3, dims3, lize3, start3, MPI_ORDER_FORTRAN, MPI_REAL, ft_3d, &
3519	ierr )
3520	CALL MPI_TYPE_COMMIT( ft_3d, ierr )
3521	ENDIF
3522	#endif
3523
3524	IF ( cyclic_fill_mode ) CALL mainrun_grid%activate_grid_from_this_class()
3525
3526	END SUBROUTINE rd_mpi_io_create_filetypes
3527
3528
3529
3530	!--------------------------------------------------------------------------------------------------!
3531	! Description:
3532	! ------------
3533	!> This subroutine creates file types to access 3d-soil arrays distributed in blocks among processes
3534	!> to a single file that contains the global arrays. It is not required for the serial mode.
3535	!--------------------------------------------------------------------------------------------------!
3536	#if defined( __parallel )
3537	SUBROUTINE rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
3538
3539	IMPLICIT NONE
3540
3541	INTEGER, INTENT(IN) :: nzb_soil !<
3542	INTEGER, INTENT(IN) :: nzt_soil !<
3543
3544	INTEGER, DIMENSION(3) :: dims3 !<
3545	INTEGER, DIMENSION(3) :: lize3 !<
3546	INTEGER, DIMENSION(3) :: start3 !<
3547
3548
3549	IF ( sm_io%is_sm_active() ) THEN
3550	CALL sm_io%sm_allocate_shared( array_3d_soil, nzb_soil, nzt_soil, sm_io%io_grid%nxl, &
3551	sm_io%io_grid%nxr, sm_io%io_grid%nys, sm_io%io_grid%nyn, &
3552	win_3ds )
3553	ELSE
3554	ALLOCATE( array_3d_soil(nzb_soil:nzt_soil,iog%nxl:iog%nxr,iog%nys:iog%nyn) )
3555	sm_io%io_grid = iog
3556	ENDIF
3557
3558	!
3559	!-- Create filetype for 3d-soil array
3560	dims3(1) = nzt_soil - nzb_soil + 1
3561	dims3(2) = iog%nx + 1
3562	dims3(3) = iog%ny + 1
3563
3564	lize3(1) = dims3(1)
3565	lize3(2) = sm_io%io_grid%nnx
3566	lize3(3) = sm_io%io_grid%nny
3567
3568	start3(1) = nzb_soil
3569	start3(2) = sm_io%io_grid%nxl
3570	start3(3) = sm_io%io_grid%nys
3571
3572	IF ( sm_io%iam_io_pe ) THEN
3573	CALL MPI_TYPE_CREATE_SUBARRAY( 3, dims3, lize3, start3, MPI_ORDER_FORTRAN, MPI_REAL, &
3574	ft_3dsoil, ierr )
3575	CALL MPI_TYPE_COMMIT( ft_3dsoil, ierr )
3576	ENDIF
3577
3578	END SUBROUTINE rd_mpi_io_create_filetypes_3dsoil
3579	#endif
3580
3581
3582
3583	!--------------------------------------------------------------------------------------------------!
3584	! Description:
3585	! ------------
3586	!> Free all file types that have been created for MPI-IO.
3587	!--------------------------------------------------------------------------------------------------!
3588	SUBROUTINE rd_mpi_io_free_filetypes
3589
3590	IMPLICIT NONE
3591
3592
3593	#if defined( __parallel )
3594	IF ( filetypes_created ) THEN
3595
3596	IF ( sm_io%iam_io_pe ) THEN
3597	CALL MPI_TYPE_FREE( ft_2d, ierr )
3598	CALL MPI_TYPE_FREE( ft_2di_nb, ierr )
3599	CALL MPI_TYPE_FREE( ft_3d, ierr )
3600	ENDIF
3601
3602	IF ( sm_io%is_sm_active() ) THEN
3603	CALL sm_io%sm_free_shared( win_2dr )
3604	CALL sm_io%sm_free_shared( win_2di )
3605	CALL sm_io%sm_free_shared( win_3dr )
3606	ELSE
3607	DEALLOCATE( array_2d, array_2di, array_3d )
3608	ENDIF
3609
3610	ENDIF
3611	!
3612	!-- Free last surface filetype
3613	IF ( sm_io%iam_io_pe .AND. ft_surf /= -1 ) THEN
3614	CALL MPI_TYPE_FREE( ft_surf, ierr )
3615	ENDIF
3616
3617	IF ( sm_io%is_sm_active() .AND. win_surf /= -1 ) THEN
3618	CALL sm_io%sm_free_shared( win_surf )
3619	ENDIF
3620
3621	ft_surf = -1
3622	win_surf = -1
3623	#else
3624	DEALLOCATE( array_2d, array_2di, array_3d )
3625	#endif
3626
3627	END SUBROUTINE rd_mpi_io_free_filetypes
3628
3629
3630
3631	!--------------------------------------------------------------------------------------------------!
3632	! Description:
3633	! ------------
3634	!> Print the restart data file header (MPI-IO format) for debugging.
3635	!--------------------------------------------------------------------------------------------------!
3636	SUBROUTINE rd_mpi_io_print_header
3637
3638	IMPLICIT NONE
3639
3640	INTEGER(iwp) :: i !<
3641
3642
3643	WRITE (9,*) 'header position after reading the restart file header: ', header_position
3644	WRITE (9,*) ' '
3645	WRITE (9,*) 'restart file header content:'
3646	WRITE (9,*) '----------------------------'
3647	WRITE (9,*) ' '
3648
3649	WRITE (9,*) ' CHARACTER header values Total number: ', tgh%nr_char
3650	WRITE (9,*) ' '
3651	DO i = 1, tgh%nr_char
3652	WRITE( 9, '(I3,A,1X,A)' ) i, ': ', text_lines(i)(1:80)
3653	ENDDO
3654	WRITE (9,*) ' '
3655
3656	WRITE (9,*) ' INTEGER header variables and values Total number: ', tgh%nr_int
3657	WRITE (9,*) ' '
3658	DO i = 1, tgh%nr_int
3659	WRITE(9,*) ' variable: ', int_names(i), ' value: ', int_values(i)
3660	ENDDO
3661	WRITE (9,*) ' '
3662
3663	WRITE (9,*) ' REAL header variables and values Total number: ', tgh%nr_real
3664	WRITE (9,*) ' '
3665	DO i = 1, tgh%nr_real
3666	WRITE(9,*) ' variable: ', real_names(i), ' value: ', real_values(i)
3667	ENDDO
3668	WRITE (9,*) ' '
3669
3670	WRITE (9,*) ' Header entries with offset (2d/3d arrays) Total number: ', tgh%nr_arrays
3671	WRITE (9,*) ' '
3672	DO i = 1, tgh%nr_arrays
3673	WRITE(9,*) ' variable: ', array_names(i), ' offset: ', array_offset(i)
3674	ENDDO
3675	WRITE (9,*) ' '
3676
3677	END SUBROUTINE rd_mpi_io_print_header
3678
3679
3680
3681	!--------------------------------------------------------------------------------------------------!
3682	! Description:
3683	! ------------
3684	!> Check if big/little endian data format is used.
3685	!> An int4 pointer is set to a int8 variable, the int*8 is set to 1, and then it is checked, if
3686	!> the first 4 bytes of the pointer are equal 1 (little endian) or not.
3687	!--------------------------------------------------------------------------------------------------!
3688	SUBROUTINE rd_mpi_io_check_endian( i_endian )
3689
3690	IMPLICIT NONE
3691
3692	INTEGER, INTENT(out) :: i_endian !<
3693	INTEGER(KIND=8), TARGET :: int8 !<
3694
3695	INTEGER, DIMENSION(1) :: bufshape !<
3696	INTEGER(KIND=4), POINTER, DIMENSION(:) :: int4 !<
3697
3698	TYPE(C_PTR) :: ptr !<
3699
3700
3701	ptr = C_LOC( int8 )
3702	bufshape(1) = 2
3703	CALL C_F_POINTER( ptr, int4, bufshape )
3704
3705	int8 = 1
3706
3707	IF ( int4(1) == 1 ) THEN
3708	i_endian = 1 ! Little endian
3709	ELSE
3710	i_endian = 2 ! Big endian
3711	ENDIF
3712
3713	END SUBROUTINE rd_mpi_io_check_endian
3714
3715	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 |