Index: /palm/trunk/SOURCE/Makefile
===================================================================
--- /palm/trunk/SOURCE/Makefile	(revision 4533)
+++ /palm/trunk/SOURCE/Makefile	(revision 4534)
@@ -25,4 +25,7 @@
 # -----------------
 # $Id$
+# shared_memory_io_mod included and respective dependencies added
+# 
+# 4525 2020-05-10 17:05:07Z raasch
 # dependency for salsa_mod updated
 # 
@@ -289,4 +292,5 @@
 	run_control.f90 \
 	salsa_mod.f90 \
+	shared_memory_io_mod.f90 \
 	singleton_mod.f90 \
 	sor.f90 \
@@ -1094,5 +1098,6 @@
 	mod_kinds.o \
 	exchange_horiz_mod.o \
-	posix_interface_mod.o
+	posix_interface_mod.o \
+	shared_memory_io_mod.o
 run_control.o: \
 	cpulog_mod.o \
@@ -1118,4 +1123,6 @@
 singleton_mod.o: \
 	mod_kinds.o
+shared_memory_io_mod.o: \
+	modules.o
 sor.o: \
 	exchange_horiz_mod.o \
Index: /palm/trunk/SOURCE/check_parameters.f90
===================================================================
--- /palm/trunk/SOURCE/check_parameters.f90	(revision 4533)
+++ /palm/trunk/SOURCE/check_parameters.f90	(revision 4534)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! adjustments for I/O on reduced number of cores using shared memory MPI
+! 
+! 4514 2020-04-30 16:29:59Z suehring
 ! Enable output of qsurf and ssurf
 ! 
@@ -256,5 +259,6 @@
 !-- Check and set the restart data format variables
     IF ( TRIM( restart_data_format ) /= 'fortran_binary'  .AND.                                    &
-         TRIM( restart_data_format ) /= 'mpi' )  THEN
+         TRIM( restart_data_format ) /= 'mpi'             .AND.                                    &
+         TRIM( restart_data_format ) /= 'mpi_shared_memory' )  THEN
        message_string = 'illegal restart data format "' // TRIM( restart_data_format ) // '"'
        CALL message( 'check_parameters', 'PA....', 1, 2, 0, 6, 0 )
@@ -269,5 +273,6 @@
 
     IF ( TRIM( restart_data_format_input ) /= 'fortran_binary'  .AND.                              &
-         TRIM( restart_data_format_input ) /= 'mpi' )  THEN
+         TRIM( restart_data_format_input ) /= 'mpi'             .AND.                              &
+         TRIM( restart_data_format_input ) /= 'mpi_shared_memory' )  THEN
        message_string = 'illegal restart data input format "' //                                   &
                         TRIM( restart_data_format_input ) // '"'
@@ -275,7 +280,16 @@
     ENDIF
     IF ( TRIM( restart_data_format_output ) /= 'fortran_binary'  .AND.                             &
-         TRIM( restart_data_format_output ) /= 'mpi' )  THEN
+         TRIM( restart_data_format_output ) /= 'mpi'             .AND.                             &
+         TRIM( restart_data_format_output ) /= 'mpi_shared_memory' )  THEN
        message_string = 'illegal restart data output format "' //                                  &
                         TRIM( restart_data_format_output ) // '"'
+       CALL message( 'check_parameters', 'PA....', 1, 2, 0, 6, 0 )
+    ENDIF
+
+    IF ( ( TRIM( restart_data_format_input ) == 'mpi_shared_memory'  .AND.                         &
+           TRIM( restart_data_format_output ) /= 'mpi_shared_memory' )  .OR.                       &
+         ( TRIM( restart_data_format_input ) /= 'mpi_shared_memory'  .AND.                         &
+           TRIM( restart_data_format_output ) == 'mpi_shared_memory' ) )  THEN
+       message_string = 'restart data formats both must be set "mpi_shared_memory"'
        CALL message( 'check_parameters', 'PA....', 1, 2, 0, 6, 0 )
     ENDIF
Index: /palm/trunk/SOURCE/land_surface_model_mod.f90
===================================================================
--- /palm/trunk/SOURCE/land_surface_model_mod.f90	(revision 4533)
+++ /palm/trunk/SOURCE/land_surface_model_mod.f90	(revision 4534)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! bugfix for switching on restart data output with MPI-IO
+! 
+! 4517 2020-05-03 14:29:30Z raasch
 ! added restart with MPI-IO for reading local arrays
 ! 
@@ -6660,5 +6663,5 @@
        ENDDO
 
-    ELSEIF ( TRIM( restart_data_format_output ) == 'fortran_binary' )  THEN
+    ELSEIF ( TRIM( restart_data_format_output ) == 'mpi' )  THEN
 
        IF ( ALLOCATED( c_liq_av ) )  CALL wrd_mpi_io( 'c_liq_av',  c_liq_av )
Index: /palm/trunk/SOURCE/read_restart_data_mod.f90
===================================================================
--- /palm/trunk/SOURCE/read_restart_data_mod.f90	(revision 4533)
+++ /palm/trunk/SOURCE/read_restart_data_mod.f90	(revision 4534)
@@ -25,4 +25,7 @@
 ! -----------------
 ! $Id$
+! adjustments for I/O on reduced number of cores using shared memory MPI
+! 
+! 4518 2020-05-04 15:44:28Z suehring
 ! Move input of diagnostic output quantities to doq_rrd_local
 ! 
@@ -833,5 +836,5 @@
        CALL close_file( 13 )
 
-    ELSEIF ( TRIM( restart_data_format_input ) == 'mpi' )  THEN
+    ELSEIF ( restart_data_format_input(1:3) == 'mpi' )  THEN
 !
 !--    Read global restart data using MPI-IO
@@ -843,5 +846,6 @@
 !
 !--    Open the MPI-IO restart file.
-       CALL rd_mpi_io_open( 'read', 'BININ' // TRIM( coupling_char ) )
+       CALL rd_mpi_io_open( 'read', 'BININ' // TRIM( coupling_char ),                              &
+                            open_for_global_io_only = .TRUE. )
 
 !
@@ -2158,5 +2162,5 @@
 
 
-    ELSEIF ( TRIM( restart_data_format_input ) == 'mpi' )  THEN
+    ELSEIF ( restart_data_format_input(1:3) == 'mpi' )  THEN
 
 !
Index: /palm/trunk/SOURCE/restart_data_mpi_io_mod.f90
===================================================================
--- /palm/trunk/SOURCE/restart_data_mpi_io_mod.f90	(revision 4533)
+++ /palm/trunk/SOURCE/restart_data_mpi_io_mod.f90	(revision 4534)
@@ -24,4 +24,7 @@
 ! -----------------
 ! $Id$
+! I/O on reduced number of cores added (using shared memory MPI)
+! 
+! 4500 2020-04-17 10:12:45Z suehring
 ! Fix too long lines
 ! 
@@ -60,5 +63,5 @@
 
     USE control_parameters,                                                                        &
-        ONLY:  include_total_domain_boundaries
+        ONLY:  include_total_domain_boundaries, restart_data_format_input
 
     USE exchange_horiz_mod,                                                                        &
@@ -73,58 +76,78 @@
         ONLY:  comm1dx, comm1dy, comm2d, myid, myidx, myidy, npex, npey, numprocs, pdims
 
-
-    IMPLICIT NONE
+    USE shared_memory_io_mod,                                                                      &
+        ONLY:  local_boundaries, sm_class
+
+
+    IMPLICIT NONE
+
+    CHARACTER(LEN=128) :: io_file_name  !> internal variable to communicate filename between
+                                        !> different subroutines
+
+#if defined( __parallel )
+    INTEGER(iwp)            ::  ierr                              !< error status of MPI-calls
+    INTEGER(iwp), PARAMETER ::  rd_offset_kind = MPI_OFFSET_KIND  !< Adress or Offset kind
+    INTEGER(iwp), PARAMETER ::  rd_status_size = MPI_STATUS_SIZE  !<
+#else
+    INTEGER(iwp), PARAMETER ::  rd_offset_kind = C_SIZE_T         !<
+    INTEGER(iwp), PARAMETER ::  rd_status_size = 1       !< Not required in sequential mode
+#endif
+
+    INTEGER(iwp)            ::  debug_level = 1 !< TODO: replace with standard debug output steering
+
+    INTEGER(iwp)            ::  comm_io       !< Communicator for MPI-IO
+    INTEGER(iwp)            ::  fh            !< MPI-IO file handle
+#if defined( __parallel )
+    INTEGER(iwp)            ::  fhs = -1      !< MPI-IO file handle to open file with comm2d always
+#endif
+    INTEGER(iwp)            ::  ft_surf = -1  !< MPI filetype surface data
+#if defined( __parallel )
+    INTEGER(iwp)            ::  ft_2di_nb     !< MPI filetype 2D array INTEGER no outer boundary
+    INTEGER(iwp)            ::  ft_2d         !< MPI filetype 2D array REAL with outer boundaries
+    INTEGER(iwp)            ::  ft_3d         !< MPI filetype 3D array REAL with outer boundaries
+    INTEGER(iwp)            ::  ft_3dsoil     !< MPI filetype for 3d-soil array
+#endif
+    INTEGER(iwp)            ::  glo_start     !< global start index on this PE
+#if defined( __parallel )
+    INTEGER(iwp)            ::  local_start   !<
+#endif
+    INTEGER(iwp)            ::  nr_iope       !<
+    INTEGER(iwp)            ::  nr_val        !< local number of values in x and y direction
+#if defined( __parallel )
+    INTEGER(iwp)            ::  win_2di
+    INTEGER(iwp)            ::  win_2dr
+    INTEGER(iwp)            ::  win_3dr
+    INTEGER(iwp)            ::  win_3ds
+    INTEGER(iwp)            ::  win_surf = -1
+#endif
+    INTEGER(iwp)            ::  total_number_of_surface_values    !< total number of values for one variable
+
+    INTEGER(KIND=rd_offset_kind) ::  array_position   !<
+    INTEGER(KIND=rd_offset_kind) ::  header_position  !<
+
+    INTEGER(iwp), DIMENSION(:,:), POINTER, CONTIGUOUS ::  array_2di  !<
+
+    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_end_index     !<
+    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_start_index   !<
+    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_global_start  !<
 
     LOGICAL ::  all_pes_write                     !< all PEs have data to write
     LOGICAL ::  filetypes_created                 !<
+    LOGICAL ::  io_on_limited_cores_per_node      !< switch to shared memory MPI-IO
     LOGICAL ::  print_header_now = .TRUE.         !<
     LOGICAL ::  rd_flag                           !< file is opened for read
     LOGICAL ::  wr_flag                           !< file is opened for write
 
-#if defined( __parallel )
-    INTEGER(iwp)                     :: ierr                              !< error status of MPI-calls
-    INTEGER(iwp), PARAMETER          :: rd_offset_kind = MPI_OFFSET_KIND  !< Adress or Offset kind
-    INTEGER(iwp), PARAMETER          :: rd_status_size = MPI_STATUS_SIZE  !<
-#else
-    INTEGER(iwp), PARAMETER          :: rd_offset_kind = C_SIZE_T         !<
-    INTEGER(iwp), PARAMETER          :: rd_status_size = 1       !< Not required in sequential mode
-#endif
-
-    INTEGER(iwp)                     :: debug_level = 1 !< TODO: replace with standard debug output steering
-
-    INTEGER(iwp)                     :: fh            !< MPI-IO file handle
-    INTEGER(iwp)                     :: ft_surf = -1  !< MPI filetype surface data
-#if defined( __parallel )
-    INTEGER(iwp)                     :: ft_2di_nb     !< MPI filetype 2D array INTEGER no outer boundary
-    INTEGER(iwp)                     :: ft_2d         !< MPI filetype 2D array REAL with outer boundaries
-    INTEGER(iwp)                     :: ft_3d         !< MPI filetype 3D array REAL with outer boundaries
-    INTEGER(iwp)                     :: ft_3dsoil     !< MPI filetype for 3d-soil array
-#endif
-    INTEGER(iwp)                     :: glo_start     !< global start index on this PE
-    INTEGER(iwp)                     :: nr_val        !< local number of values in x and y direction
-    INTEGER(iwp)                     :: total_number_of_surface_values    !< total number of values for one variable
-
-    INTEGER(KIND=rd_offset_kind) ::  array_position   !<
-    INTEGER(KIND=rd_offset_kind) ::  header_position  !<
-
-    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_end_index     !<
-    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_start_index   !<
-    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_global_start  !<
-
     REAL(KIND=wp) ::  mb_processed  !<
 
+#if defined( __parallel )
+    REAL(wp), DIMENSION(:), POINTER, CONTIGUOUS     ::  array_1d       !<
+#endif
+    REAL(wp), DIMENSION(:,:), POINTER, CONTIGUOUS   ::  array_2d       !<
+    REAL(wp), DIMENSION(:,:,:), POINTER, CONTIGUOUS ::  array_3d       !<
+    REAL(wp), DIMENSION(:,:,:), POINTER, CONTIGUOUS ::  array_3d_soil  !<
+
 !
 !-- Handling of outer boundaries
-    TYPE local_boundaries
-       INTEGER(iwp) ::  nnx
-       INTEGER(iwp) ::  nny
-       INTEGER(iwp) ::  nx
-       INTEGER(iwp) ::  nxl
-       INTEGER(iwp) ::  nxr
-       INTEGER(iwp) ::  ny
-       INTEGER(iwp) ::  nyn
-       INTEGER(iwp) ::  nys
-    END TYPE local_boundaries
-
     TYPE(local_boundaries) ::  lb  !<
 
@@ -144,4 +167,6 @@
     TYPE(general_header), TARGET ::  tgh
 
+    TYPE(sm_class)               ::  sm_io
+
 !
 !-- Declaration of varibales for file header section
@@ -250,5 +275,5 @@
 !> Open restart file for read or write with MPI-IO
 !--------------------------------------------------------------------------------------------------!
- SUBROUTINE rd_mpi_io_open( action, file_name )
+ SUBROUTINE rd_mpi_io_open( action, file_name, open_for_global_io_only )
 
     IMPLICIT NONE
@@ -266,4 +291,7 @@
 #endif
 
+    LOGICAL, INTENT(IN), OPTIONAL ::  open_for_global_io_only          !<
+    LOGICAL                       ::  set_filetype                     !<
+
 #if ! defined( __parallel )
     TYPE(C_PTR)                   ::  buf_ptr                          !<
@@ -275,4 +303,32 @@
     rd_flag = ( TRIM( action ) == 'READ'  .OR. TRIM( action ) == 'read'  )
     wr_flag = ( TRIM( action ) == 'WRITE' .OR. TRIM( action ) == 'write' )
+
+!
+!-- Store name of I/O file to communicate it internally within this module.
+    io_file_name = file_name
+!
+!-- Setup for IO on a limited number of threads per node (using shared memory MPI)
+    IF ( TRIM( restart_data_format_input ) == 'mpi_shared_memory' )  THEN
+       io_on_limited_cores_per_node = .TRUE.
+       set_filetype                 = .TRUE.
+    ENDIF
+!
+!-- Shared memory MPI is not used for reading of global data
+    IF ( PRESENT( open_for_global_io_only )  .AND.  rd_flag )  THEN
+       IF ( open_for_global_io_only )  THEN
+          io_on_limited_cores_per_node = .FALSE.
+          set_filetype                 = .FALSE.
+       ENDIF
+    ENDIF
+
+    CALL sm_io%sm_init_comm( io_on_limited_cores_per_node )
+
+!
+!-- Set communicator to be used. If all cores are doing I/O, comm2d is used as usual.
+    IF( sm_io%is_sm_active() )  THEN
+       comm_io = sm_io%comm_io
+    ELSE
+       comm_io = comm2d
+    ENDIF
 
 !
@@ -284,5 +340,5 @@
 !-- further below.
     IF ( wr_flag)  THEN
-       CALL rs_mpi_io_create_filetypes
+       CALL rd_mpi_io_create_filetypes
        filetypes_created = .TRUE.
     ENDIF
@@ -291,26 +347,29 @@
 !-- Open file for MPI-IO
 #if defined( __parallel )
-    IF ( rd_flag )  THEN
-       CALL MPI_FILE_OPEN( comm2d, TRIM( file_name ), MPI_MODE_RDONLY, MPI_INFO_NULL, fh, ierr )
-       WRITE (9,*) 'Open MPI-IO restart file for read  ==> ', TRIM( file_name )
-    ELSEIF ( wr_flag )  THEN
-       CALL MPI_FILE_OPEN( comm2d, TRIM( file_name ), MPI_MODE_CREATE+MPI_MODE_WRONLY,             &
-                           MPI_INFO_NULL, fh, ierr )
-       WRITE (9,*) 'Open MPI-IO restart file for write ==> ', TRIM( file_name )
-    ELSE
-       CALL rs_mpi_io_error( 1 )
+    IF ( sm_io%iam_io_pe )  THEN
+       IF ( rd_flag )  THEN
+          CALL MPI_FILE_OPEN( comm_io, TRIM( io_file_name ), MPI_MODE_RDONLY, MPI_INFO_NULL, fh,   &
+                              ierr )
+          WRITE (9,*) 'Open MPI-IO restart file for read  ==> ', TRIM( io_file_name )
+       ELSEIF ( wr_flag )  THEN
+          CALL MPI_FILE_OPEN( comm_io, TRIM( io_file_name ), MPI_MODE_CREATE+MPI_MODE_WRONLY,      &
+                              MPI_INFO_NULL, fh, ierr )
+          WRITE (9,*) 'Open MPI-IO restart file for write ==> ', TRIM( io_file_name )
+       ELSE
+          CALL rd_mpi_io_error( 1 )
+       ENDIF
     ENDIF
 #else
     IF ( rd_flag )  THEN
-       fh = posix_open( TRIM( file_name ), .TRUE. )
-       WRITE (9,*) 'Open sequential restart file for read  ==> ', TRIM( file_name ), ' ', fh
+       fh = posix_open( TRIM( io_file_name ), .TRUE. )
+       WRITE (9,*) 'Open sequential restart file for read  ==> ', TRIM( io_file_name ), ' ', fh
     ELSEIF ( wr_flag )  THEN
-       fh = posix_open( TRIM( file_name ), .FALSE. )
-       WRITE (9,*) 'Open sequential restart file for write ==> ', TRIM(file_name), ' ', fh
+       fh = posix_open( TRIM( io_file_name ), .FALSE. )
+       WRITE (9,*) 'Open sequential restart file for write ==> ', TRIM( io_file_name ), ' ', fh
     ELSE
-       CALL rs_mpi_io_error( 1 )
-    ENDIF
-
-    IF ( fh < 0 )  CALL rs_mpi_io_error( 6 )
+       CALL rd_mpi_io_error( 1 )
+    ENDIF
+
+    IF ( fh < 0 )  CALL rd_mpi_io_error( 6 )
 #endif
 
@@ -347,17 +406,24 @@
 
     IF ( rd_flag )  THEN
-!
-!--    File is open for read.
-#if defined( __parallel )
-!--    Set the default view
-       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
-!
-!--    Read the file header size
-       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ( fh, tgh, gh_size, MPI_BYTE, status, ierr )
+       IF ( sm_io%iam_io_pe )  THEN
+!
+!--       File is open for read.
+#if defined( __parallel )
+!--       Set the default view
+          CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+!
+!--       Read the file header size
+          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ( fh, tgh, gh_size, MPI_BYTE, status, ierr )
 #else
-       CALL posix_lseek( fh, header_position )
-       buf_ptr = C_LOC( tgh )
-       CALL posix_read( fh, buf_ptr, gh_size )
+          CALL posix_lseek( fh, header_position )
+          buf_ptr = C_LOC( tgh )
+          CALL posix_read( fh, buf_ptr, gh_size )
+#endif
+       ENDIF
+#if defined( __parallel )
+       IF ( sm_io%is_sm_active() )  THEN
+          CALL MPI_BCAST( tgh, gh_size, MPI_BYTE, 0, sm_io%comm_shared, ierr )
+       ENDIF
 #endif
        header_position = header_position + gh_size
@@ -368,78 +434,114 @@
 !--    File types depend on if boundaries of the total domain is included in data. This has been
 !--    checked with the previous statement.
-       CALL rs_mpi_io_create_filetypes
-       filetypes_created = .TRUE.
-
-#if defined( __parallel )
-!
-!--    Read INTEGER values
-       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ( fh, int_names, SIZE( int_names ) * 32, MPI_CHAR, status, ierr )
-       header_position = header_position + SIZE( int_names ) * 32
-
-       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ (fh, int_values, SIZE( int_values ), MPI_INT, status, ierr )
-       header_position = header_position + SIZE( int_values ) * iwp
-!
-!--    Character entries
-       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ( fh, text_lines, SIZE( text_lines ) * 128, MPI_CHAR, status, ierr )
-       header_position = header_position+size(text_lines) * 128
-!
-!--    REAL values
-       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ( fh, real_names, SIZE( real_names ) * 32, MPI_CHAR, status, ierr )
-       header_position = header_position + SIZE( real_names ) * 32
-
-       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ( fh, real_values, SIZE( real_values ), MPI_REAL, status, ierr )
-       header_position = header_position + SIZE( real_values ) * wp
-!
-!--    2d- and 3d-array headers
-       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ( fh, array_names, SIZE( array_names ) * 32, MPI_CHAR, status, ierr )
-       header_position = header_position + SIZE( array_names ) * 32
-
-       CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ( fh, array_offset, SIZE( array_offset ) * MPI_OFFSET_KIND, MPI_BYTE,     &
-                           status,ierr )   ! there is no I*8 datatype in Fortran
-       header_position = header_position + SIZE( array_offset ) * rd_offset_kind
+       IF ( set_filetype )  THEN
+          CALL rd_mpi_io_create_filetypes
+          filetypes_created = .TRUE.
+       ENDIF
+
+       IF ( sm_io%iam_io_pe )  THEN
+#if defined( __parallel )
+!
+!--       Read INTEGER values
+          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ( fh, int_names, SIZE( int_names ) * 32, MPI_CHAR, status, ierr )
+          header_position = header_position + SIZE( int_names ) * 32
+
+          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ (fh, int_values, SIZE( int_values ), MPI_INT, status, ierr )
+          header_position = header_position + SIZE( int_values ) * iwp
+!
+!--       Character entries
+          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ( fh, text_lines, SIZE( text_lines ) * 128, MPI_CHAR, status, ierr )
+          header_position = header_position+size(text_lines) * 128
+!
+!--       REAL values
+          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ( fh, real_names, SIZE( real_names ) * 32, MPI_CHAR, status, ierr )
+          header_position = header_position + SIZE( real_names ) * 32
+
+          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ( fh, real_values, SIZE( real_values ), MPI_REAL, status, ierr )
+          header_position = header_position + SIZE( real_values ) * wp
+!
+!--       2d- and 3d-array headers
+          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ( fh, array_names, SIZE( array_names ) * 32, MPI_CHAR, status, ierr )
+          header_position = header_position + SIZE( array_names ) * 32
+
+          CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ( fh, array_offset, SIZE( array_offset ) * MPI_OFFSET_KIND, MPI_BYTE,  &
+                              status,ierr )   ! there is no I*8 datatype in Fortran
+          header_position = header_position + SIZE( array_offset ) * rd_offset_kind
 #else
-       CALL posix_lseek( fh, header_position )
-       CALL posix_read( fh, int_names )
-       header_position = header_position + SIZE( int_names ) * 32
-
-       CALL posix_lseek( fh, header_position )
-       CALL posix_read( fh, int_values, SIZE( int_values ) )
-       header_position = header_position + SIZE( int_values ) * iwp
-!
-!--    Character entries
-       CALL posix_lseek( fh, header_position )
-       CALL posix_read( fh, text_lines )
-       header_position = header_position + SIZE( text_lines ) * 128
-!
-!--    REAL values
-       CALL posix_lseek( fh, header_position )
-       CALL posix_read( fh, real_names )
-       header_position = header_position + SIZE( real_names ) * 32
-
-       CALL posix_lseek( fh, header_position )
-       CALL posix_read( fh, real_values, SIZE( real_values ) )
-       header_position = header_position + SIZE( real_values ) * wp
-!
-!--    2d- and 3d-array headers
-       CALL posix_lseek( fh, header_position )
-       CALL posix_read( fh, array_names )
-       header_position = header_position + SIZE( array_names ) * 32
-
-       CALL posix_lseek( fh, header_position )
-       CALL posix_read( fh, array_offset, SIZE( array_offset ) ) ! there is no I*8 datatype in Fortran
-       header_position = header_position + SIZE( array_offset ) * rd_offset_kind
-#endif
-       IF ( debug_level >= 2 )  THEN
-          WRITE (9,*) 'header positio after array metadata  ', header_position
-       ENDIF
-
-       IF ( print_header_now )  CALL rs_mpi_io_print_header
+          CALL posix_lseek( fh, header_position )
+          CALL posix_read( fh, int_names )
+          header_position = header_position + SIZE( int_names ) * 32
+
+          CALL posix_lseek( fh, header_position )
+          CALL posix_read( fh, int_values, SIZE( int_values ) )
+          header_position = header_position + SIZE( int_values ) * iwp
+!
+!--       Character entries
+          CALL posix_lseek( fh, header_position )
+          CALL posix_read( fh, text_lines )
+          header_position = header_position + SIZE( text_lines ) * 128
+!
+!--       REAL values
+          CALL posix_lseek( fh, header_position )
+          CALL posix_read( fh, real_names )
+          header_position = header_position + SIZE( real_names ) * 32
+
+          CALL posix_lseek( fh, header_position )
+          CALL posix_read( fh, real_values, SIZE( real_values ) )
+          header_position = header_position + SIZE( real_values ) * wp
+!
+!--       2d- and 3d-array headers
+          CALL posix_lseek( fh, header_position )
+          CALL posix_read( fh, array_names )
+          header_position = header_position + SIZE( array_names ) * 32
+
+          CALL posix_lseek( fh, header_position )
+          CALL posix_read( fh, array_offset, SIZE( array_offset ) ) ! there is no I*8 datatype in Fortran
+          header_position = header_position + SIZE( array_offset ) * rd_offset_kind
+#endif
+          IF ( debug_level >= 2 )  THEN
+             WRITE (9,*) 'header positio after array metadata  ', header_position
+          ENDIF
+
+          IF ( print_header_now )  CALL rd_mpi_io_print_header
+
+       ENDIF
+
+#if defined( __parallel )
+!
+!--    Broadcast header to all remaining cores that are not involved in I/O
+       IF ( sm_io%is_sm_active() )  THEN
+!
+!--        Not sure, that it is possible to broadcast CHARACTER array in one MPI_Bcast call
+           DO  i = 1, SIZE( int_names )
+              CALL MPI_BCAST( int_names(i), 32, MPI_CHARACTER, 0, sm_io%comm_shared, ierr )
+           ENDDO
+           CALL MPI_BCAST( int_values, SIZE( int_values ), MPI_INTEGER, 0, sm_io%comm_shared, ierr )
+
+           DO  i = 1, SIZE( text_lines )
+              CALL MPI_BCAST( text_lines(i), 128, MPI_CHARACTER, 0, sm_io%comm_shared, ierr )
+           ENDDO
+
+           DO  i = 1, SIZE( real_names )
+              CALL MPI_BCAST( real_names(i), 32, MPI_CHARACTER, 0, sm_io%comm_shared, ierr )
+           ENDDO
+           CALL MPI_BCAST( real_values, SIZE( real_values ), MPI_REAL, 0, sm_io%comm_shared, ierr )
+
+           DO  i = 1, SIZE( array_names )
+              CALL MPI_BCAST( array_names(i), 32, MPI_CHARACTER, 0, sm_io%comm_shared, ierr )
+           ENDDO
+           CALL MPI_BCAST( array_offset, SIZE( array_offset )*8, MPI_BYTE, 0, sm_io%comm_shared,   &
+                           ierr )  ! there is no I*8 datatype in Fortran (array_offset is I*8!)
+
+           CALL MPI_BCAST( header_position, rd_offset_kind, MPI_BYTE, 0, sm_io%comm_shared, ierr )
+
+       ENDIF
+#endif
 
     ENDIF
@@ -515,5 +617,5 @@
     IF ( .NOT. lo_found )  THEN
        WRITE(9,*)  'INTEGER not found ', name
-       CALL rs_mpi_io_error( 3 )
+       CALL rd_mpi_io_error( 3 )
     ENDIF
 
@@ -560,5 +662,5 @@
     IF ( .NOT. lo_found )  THEN
        WRITE(9,*) 'REAL value not found ', name
-       CALL rs_mpi_io_error(3)
+       CALL rd_mpi_io_error(3)
     ENDIF
 
@@ -586,6 +688,4 @@
 
     REAL(wp), INTENT(INOUT), DIMENSION(nysg:nyng,nxlg:nxrg) ::  data
-
-    REAL(KIND=wp), DIMENSION(lb%nxl:lb%nxr,lb%nys:lb%nyn)   ::  array_2d
 
 
@@ -600,31 +700,36 @@
     ENDDO
 
-     IF ( found )  THEN
-#if defined( __parallel )
-        CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL, ierr )
-        CALL MPI_FILE_READ_ALL( fh, array_2d, SIZE( array_2d), MPI_REAL, status, ierr )
+    IF ( found )  THEN
+#if defined( __parallel )
+       CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
+       IF ( sm_io%iam_io_pe )  THEN
+          CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL,   &
+                                  ierr )
+          CALL MPI_FILE_READ_ALL( fh, array_2d, SIZE( array_2d ), MPI_REAL, status, ierr )
+       ENDIF
+       CALL sm_io%sm_node_barrier()
 #else
-        CALL posix_lseek( fh, array_position )
-        CALL posix_read( fh, array_2d, SIZE( array_2d ) )
-#endif
-
-        IF ( include_total_domain_boundaries)  THEN
-           DO  i = lb%nxl, lb%nxr
-              data(lb%nys-nbgp:lb%nyn-nbgp,i-nbgp) = array_2d(i,lb%nys:lb%nyn)
-           ENDDO
-           IF ( debug_level >= 2)  WRITE(9,*) 'r2f_ob ', TRIM(name),' ', SUM( data(nys:nyn,nxl:nxr) )
-        ELSE
-           DO  i = nxl, nxr
-              data(nys:nyn,i) = array_2d(i,nys:nyn)
-           ENDDO
-           IF ( debug_level >= 2) WRITE(9,*) 'r2f ', TRIM( name ),' ', SUM( data)
-        ENDIF
-
-        CALL exchange_horiz_2d( data )
-
-     ELSE
-        WRITE(9,*) 'array_2D not found ', name
-        CALL rs_mpi_io_error( 2 )
-     ENDIF
+       CALL posix_lseek( fh, array_position )
+       CALL posix_read( fh, array_2d, SIZE( array_2d ) )
+#endif
+
+       IF ( include_total_domain_boundaries)  THEN
+          DO  i = lb%nxl, lb%nxr
+             data(lb%nys-nbgp:lb%nyn-nbgp,i-nbgp) = array_2d(i,lb%nys:lb%nyn)
+          ENDDO
+          IF ( debug_level >= 2)  WRITE(9,*) 'r2f_ob ', TRIM(name),' ', SUM( data(nys:nyn,nxl:nxr) )
+       ELSE
+          DO  i = nxl, nxr
+             data(nys:nyn,i) = array_2d(i,nys:nyn)
+          ENDDO
+          IF ( debug_level >= 2) WRITE(9,*) 'r2f ', TRIM( name ),' ', SUM( data(nys:nyn,nxl:nxr) )
+       ENDIF
+
+       CALL exchange_horiz_2d( data )
+
+    ELSE
+       WRITE(9,*) 'array_2D not found ', name
+       CALL rd_mpi_io_error( 2 )
+    ENDIF
 
  END SUBROUTINE rrd_mpi_io_real_2d
@@ -649,6 +754,4 @@
     INTEGER, DIMENSION(rd_status_size)  ::  status
 #endif
-
-    INTEGER, DIMENSION(nxl:nxr,nys:nyn) ::  array_2d
 
     INTEGER(KIND=iwp), INTENT(INOUT), DIMENSION(:,:) ::  data
@@ -675,5 +778,5 @@
 !--                  would be dimensioned in the caller subroutine like this:
 !--                  INTEGER, DIMENSION(nysg:nyng,nxlg:nxrg)::  data
-          CALL rs_mpi_io_error( 2 )
+          CALL rd_mpi_io_error( 2 )
 
        ELSEIF ( (nxr-nxl+1) == SIZE( data, 2 ) )  THEN
@@ -684,23 +787,27 @@
 
 #if defined( __parallel )
-          CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native',            &
-                                  MPI_INFO_NULL, ierr )
-          CALL MPI_FILE_READ_ALL( fh, array_2d, SIZE( array_2d), MPI_INTEGER, status, ierr )
+          CALL sm_io%sm_node_barrier() ! has no effect if I/O on limited number of cores is inactive
+          IF ( sm_io%iam_io_pe )  THEN
+             CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native',         &
+                                     MPI_INFO_NULL, ierr )
+             CALL MPI_FILE_READ_ALL( fh, array_2d, SIZE( array_2d ), MPI_INTEGER, status, ierr )
+          ENDIF
+          CALL sm_io%sm_node_barrier()
 #else
           CALL posix_lseek( fh, array_position )
-          CALL posix_read( fh, array_2d, SIZE( array_2d ) )
+          CALL posix_read( fh, array_2di, SIZE( array_2di ) )
 #endif
 
           DO  j = nys, nyn
              DO  i = nxl, nxr
-                data(j-nys+1,i-nxl+1) = array_2d(i,j)
+                data(j-nys+1,i-nxl+1) = array_2di(i,j)
              ENDDO
           ENDDO
 
-          IF ( debug_level >= 2 )  WRITE(9,*) 'r2i ', TRIM( name ),' ', SUM( array_2d )
+          IF ( debug_level >= 2 )  WRITE(9,*) 'r2i ', TRIM( name ),' ', SUM( array_2di )
 
        ELSE
           WRITE (9,*) '### rrd_mpi_io_int_2d  array: ', TRIM( name )
-          CALL rs_mpi_io_error( 4 )
+          CALL rd_mpi_io_error( 4 )
        ENDIF
 
@@ -708,5 +815,5 @@
 
        WRITE(9,*) 'array_2D not found ', name
-       CALL rs_mpi_io_error( 2 )
+       CALL rd_mpi_io_error( 2 )
 
     ENDIF
@@ -734,6 +841,4 @@
 
     LOGICAL                            ::  found
-
-    REAL(KIND=wp), DIMENSION(nzb:nzt+1,lb%nxl:lb%nxr,lb%nys:lb%nyn)   ::  array_3d
 
     REAL(wp), INTENT(INOUT), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  data
@@ -752,6 +857,11 @@
     IF ( found )  THEN
 #if defined( __parallel )
-       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL, ierr )
-       CALL MPI_FILE_READ_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
+       CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
+       IF( sm_io%iam_io_pe )  THEN
+          CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL,    &
+                                  ierr )
+          CALL MPI_FILE_READ_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
+       ENDIF
+       CALL sm_io%sm_node_barrier()
 #else
        CALL posix_lseek( fh, array_position )
@@ -767,5 +877,5 @@
              data(:,nys:nyn,i) = array_3d(:,i,nys:nyn)
           ENDDO
-          IF ( debug_level >= 2 )  WRITE(9,*) 'r3f ', TRIM( name ),' ', SUM( data )
+          IF ( debug_level >= 2 )  WRITE(9,*) 'r3f ', TRIM( name ),' ', SUM( data(:,nys:nyn,nxl:nxr) )
        ENDIF
 
@@ -774,5 +884,5 @@
     ELSE
        WRITE(9,*)  'array_3D not found ', name
-       CALL rs_mpi_io_error(2)
+       CALL rd_mpi_io_error(2)
     ENDIF
 
@@ -803,6 +913,4 @@
 
     LOGICAL                            ::  found
-
-    REAL(KIND=wp), DIMENSION(nzb_soil:nzt_soil,lb%nxl:lb%nxr,lb%nys:lb%nyn)   ::  array_3d
 
     REAL(wp), INTENT(INOUT), DIMENSION(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ::  data
@@ -822,11 +930,15 @@
 #if defined( __parallel )
        CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
-       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3dsoil, 'native', MPI_INFO_NULL,   &
-                               ierr )
-       CALL MPI_FILE_READ_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
-       CALL MPI_TYPE_FREE( ft_3dsoil, ierr )
+       CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
+       IF ( sm_io%iam_io_pe )  THEN
+          CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3dsoil, 'native', MPI_INFO_NULL,&
+                                  ierr )
+          CALL MPI_FILE_READ_ALL( fh, array_3d_soil, SIZE( array_3d_soil ), MPI_REAL, status, ierr )
+          CALL MPI_TYPE_FREE( ft_3dsoil, ierr )
+       ENDIF
+       CALL sm_io%sm_node_barrier()
 #else
        CALL posix_lseek( fh, array_position )
-       CALL posix_read( fh, array_3d, SIZE( array_3d ) )
+       CALL posix_read( fh, array_3d_soil, SIZE( array_3d_soil ) )
 #endif
        IF ( include_total_domain_boundaries )  THEN
@@ -844,5 +956,5 @@
     ELSE
        WRITE(9,*)  'array_3D not found ', name
-       CALL rs_mpi_io_error( 2 )
+       CALL rd_mpi_io_error( 2 )
     ENDIF
 
@@ -890,5 +1002,5 @@
     IF ( .NOT. lo_found )  THEN
        WRITE(9,*)  'Character variable not found ', name
-         CALL rs_mpi_io_error( 3 )
+       CALL rd_mpi_io_error( 3 )
     ENDIF
 
@@ -975,6 +1087,4 @@
     INTEGER, DIMENSION(rd_status_size) ::  status
 #endif
-
-    REAL(KIND=wp), DIMENSION(lb%nxl:lb%nxr,lb%nys:lb%nyn)  :: array_2d
 
     REAL(wp), INTENT(IN), DIMENSION(nysg:nyng,nxlg:nxrg)    :: data
@@ -998,10 +1108,15 @@
           array_2d(i,lb%nys:lb%nyn) = data(nys:nyn,i)
        ENDDO
-       IF ( debug_level >= 2 )  WRITE(9,*)  'w2f ', TRIM( name ),' ', SUM( array_2d )
-    ENDIF
-
-#if defined( __parallel )
-    CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL, ierr )
-    CALL MPI_FILE_WRITE_ALL( fh, array_2d, SIZE( array_2d), MPI_REAL, status, ierr )
+       IF ( debug_level >= 2 )  WRITE(9,*)  'w2f ', TRIM( name ),' ',                              &
+                                            SUM( array_2d(nxl:nxr, lb%nys:lb%nyn) )
+    ENDIF
+
+#if defined( __parallel )
+    CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
+    IF ( sm_io%iam_io_pe )  THEN
+       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_2d, 'native', MPI_INFO_NULL, ierr )
+       CALL MPI_FILE_WRITE_ALL( fh, array_2d, SIZE( array_2d), MPI_REAL, status, ierr )
+    ENDIF
+    CALL sm_io%sm_node_barrier()
 #else
     CALL posix_lseek( fh, array_position )
@@ -1036,5 +1151,4 @@
 #endif
     INTEGER(KIND=iwp), INTENT(IN), DIMENSION(:,:) ::  data
-    INTEGER, DIMENSION(nxl:nxr,nys:nyn)           ::  array_2d
 
     LOGICAl, OPTIONAL                             ::  ar_found
@@ -1051,5 +1165,5 @@
 !--    INTEGER, DIMENSION(nysg:nyng,nxlg:nxrg) ::  data
        WRITE (9,*) '### wrd_mpi_io_int_2d  IF  array: ', TRIM( name )
-       CALL rs_mpi_io_error( 4 )
+       CALL rd_mpi_io_error( 4 )
 
     ELSEIF ( ( nxr-nxl+1 ) == SIZE( data, 2 ) )  THEN
@@ -1060,15 +1174,19 @@
        DO  j = nys, nyn
           DO  i = nxl, nxr
-             array_2d(i,j) = data(j-nys+1,i-nxl+1)
+             array_2di(i,j) = data(j-nys+1,i-nxl+1)
           ENDDO
        ENDDO
-       IF ( debug_level >= 2 )  WRITE(9,*) 'w2i ', TRIM( name ), ' ', SUM( array_2d )
-#if defined( __parallel )
-       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native', MPI_INFO_NULL,&
-                               ierr )
-       CALL MPI_FILE_WRITE_ALL( fh, array_2d, SIZE( array_2d ), MPI_INTEGER, status, ierr )
+       IF ( debug_level >= 2 )  WRITE(9,*) 'w2i ', TRIM( name ), ' ', SUM( array_2di(nxl:nxr,nys:nyn) ), SUM( data )
+#if defined( __parallel )
+       CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
+       IF ( sm_io%iam_io_pe )  THEN
+          CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native',            &
+                                  MPI_INFO_NULL, ierr )
+          CALL MPI_FILE_WRITE_ALL( fh, array_2di, SIZE( array_2di ), MPI_INTEGER, status, ierr )
+       ENDIF
+       CALL sm_io%sm_node_barrier()
 #else
        CALL posix_lseek( fh, array_position )
-       CALL posix_write( fh, array_2d, SIZE( array_2d ) )
+       CALL posix_write( fh, array_2di, SIZE( array_2di ) )
 #endif
 !
@@ -1080,5 +1198,5 @@
     ELSE
        WRITE (9,*) '### wrd_mpi_io_int_2d  array: ', TRIM( name )
-       CALL rs_mpi_io_error( 4 )
+       CALL rd_mpi_io_error( 4 )
     ENDIF
 
@@ -1105,6 +1223,4 @@
 #endif
     REAL(wp), INTENT(IN), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg)  ::  data
-
-    REAL(KIND=wp), DIMENSION(nzb:nzt+1,lb%nxl:lb%nxr,lb%nys:lb%nyn) ::  array_3d
 
 
@@ -1130,9 +1246,13 @@
            array_3d(:,i,lb%nys:lb%nyn) = data(:,nys:nyn,i)
        ENDDO
-       IF ( debug_level >= 2 )  WRITE(9,*)  'w3f ', TRIM( name ),' ', SUM( array_3d )
-    ENDIF
-#if defined( __parallel )
-    CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL, ierr )
-    CALL MPI_FILE_WRITE_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
+       IF ( debug_level >= 2 )  WRITE(9,*)  'w3f ', TRIM( name ),' ', SUM( data(:,nys:nyn,nxl:nxr) )
+    ENDIF
+#if defined( __parallel )
+    CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
+    IF ( sm_io%iam_io_pe )  THEN
+       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3d, 'native', MPI_INFO_NULL, ierr )
+       CALL MPI_FILE_WRITE_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
+    ENDIF
+    CALL sm_io%sm_node_barrier()
 #else
     CALL posix_lseek( fh, array_position )
@@ -1173,10 +1293,12 @@
     REAL(wp), INTENT(IN), DIMENSION(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg)  ::  data
 
-    REAL(KIND=wp), DIMENSION(nzb_soil:nzt_soil,lb%nxl:lb%nxr,lb%nys:lb%nyn) ::  array_3d
-
 
     array_names(header_arr_index)  = name
     array_offset(header_arr_index) = array_position
     header_arr_index = header_arr_index + 1
+
+#if defined( __parallel )
+    CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
+#endif
 
     IF ( include_total_domain_boundaries)  THEN
@@ -1188,5 +1310,5 @@
 !--    For this reason, the original PALM data need to be swaped.
        DO  i = lb%nxl, lb%nxr
-          array_3d(:,i,lb%nys:lb%nyn) = data(:,lb%nys-nbgp:lb%nyn-nbgp,i-nbgp)
+          array_3d_soil(:,i,lb%nys:lb%nyn) = data(:,lb%nys-nbgp:lb%nyn-nbgp,i-nbgp)
        ENDDO
        IF ( debug_level >= 2 )  WRITE(9,*) 'w3f_ob_soil ', TRIM( name ), ' ', SUM( data(:,nys:nyn,nxl:nxr) )
@@ -1195,16 +1317,19 @@
 !--    Prepare output of 3d-REAL-array without ghost layers
        DO  i = nxl, nxr
-          array_3d(:,i,lb%nys:lb%nyn) = data(:,nys:nyn,i)
+          array_3d_soil(:,i,lb%nys:lb%nyn) = data(:,nys:nyn,i)
        ENDDO
        IF ( debug_level >= 2 )  WRITE(9,*) 'w3f_soil ', TRIM( name ), ' ', SUM( array_3d )
     ENDIF
 #if defined( __parallel )
-    CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
-    CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3dsoil, 'native', MPI_INFO_NULL, ierr )
-    CALL MPI_FILE_WRITE_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
-    CALL MPI_TYPE_FREE( ft_3dsoil, ierr )
+    CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
+    IF ( sm_io%iam_io_pe )  THEN
+       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_3dsoil, 'native', MPI_INFO_NULL,   &
+                               ierr )
+       CALL MPI_FILE_WRITE_ALL( fh, array_3d_soil, SIZE( array_3d_soil ), MPI_REAL, status, ierr )
+    ENDIF
+    CALL sm_io%sm_node_barrier()
 #else
     CALL posix_lseek( fh, array_position )
-    CALL posix_write( fh, array_3d, SIZE( array_3d ) )
+    CALL posix_write( fh, array_3d_soil, SIZE( array_3d_soil ) )
 #endif
 !
@@ -1309,7 +1434,12 @@
 !--    Set default view
 #if defined( __parallel )
-       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
-       CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ_ALL( fh, data, SIZE( data ), MPI_REAL, status, ierr )
+       IF ( sm_io%iam_io_pe )  THEN
+          CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+          CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ_ALL( fh, data, SIZE( data ), MPI_REAL, status, ierr )
+       ENDIF
+       IF ( sm_io%is_sm_active() )  THEN
+          CALL MPI_BCAST( data, SIZE(data), MPI_REAL, 0, sm_io%comm_shared, ierr )
+       ENDIF
 #else
        CALL posix_lseek( fh, array_position )
@@ -1319,5 +1449,5 @@
     ELSE
        WRITE(9,*)  'replicated array_1D not found ', name
-       CALL rs_mpi_io_error( 2 )
+       CALL rd_mpi_io_error( 2 )
     ENDIF
 
@@ -1459,7 +1589,12 @@
 !--    Set default view
 #if defined( __parallel )
-       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
-       CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
-       CALL MPI_FILE_READ_ALL( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
+       IF ( sm_io%iam_io_pe )  THEN
+          CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+          CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
+          CALL MPI_FILE_READ_ALL( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
+       ENDIF
+       IF ( sm_io%is_sm_active() )  THEN
+          CALL MPI_BCAST( data, SIZE(data), MPI_INTEGER, 0, sm_io%comm_shared, ierr )
+       ENDIF
 #else
        CALL posix_lseek( fh, array_position )
@@ -1472,7 +1607,7 @@
        ELSE
           WRITE (9,*) '### rrd_mpi_io_global_array_int_1d ', TRIM( name )
-          CALL rs_mpi_io_error( 4 )
+          CALL rd_mpi_io_error( 4 )
           WRITE(9,*)  'replicated array_1D not found ', name
-          CALL rs_mpi_io_error( 2 )
+          CALL rd_mpi_io_error( 2 )
        ENDIF
     ENDIF
@@ -1515,5 +1650,7 @@
 !--    Set default view
 #if defined( __parallel )
-       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+       IF ( sm_io%iam_io_pe )  THEN
+          CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+       ENDIF
 !
 !--    Only PE 0 writes replicated data
@@ -1656,5 +1793,7 @@
 !-- Set default view
 #if defined( __parallel )
-    CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+    IF ( sm_io%iam_io_pe )  THEN
+       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+    ENDIF
 !
 !-- Only PE 0 writes replicated data
@@ -1767,12 +1906,12 @@
              ELSE                                          ! read
 #if defined( __parallel )
-                IF ( debug_level >= 2 )  WRITE(9,'(a,8i8)') 'read block ', j, i, j_f, i_f, &
-                                                            m_start_index(j_f,i_f), nr_bytes_f, disp_f
-                CALL MPI_FILE_SEEK( fh, disp_f, MPI_SEEK_SET, ierr )
+                IF ( debug_level >= 2 )  WRITE(9,'(a,4i4,4i10)') 'read block ', j, i, j_f, i_f,    &
+                                                         m_start_index(j_f,i_f), nr_bytes_f, disp_f
+                CALL MPI_FILE_SEEK( fhs, disp_f, MPI_SEEK_SET, ierr )
                 nr_words = nr_bytes_f / wp
-                CALL MPI_FILE_READ( fh, data(m_start_index(j_f,i_f)), nr_words, MPI_REAL, status, ierr )
+                CALL MPI_FILE_READ( fhs, data(m_start_index(j_f,i_f)), nr_words, MPI_REAL, status, ierr )
 #else
                 CALL posix_lseek( fh, disp_f )
-!                CALL posix_read( fh, data(m_start_index(j_f:,i_f:)), nr_bytes_f )
+                CALL posix_read( fh, data(m_start_index(j_f,i_f):), nr_bytes_f )
 #endif
                 disp_f     = disp
@@ -1787,9 +1926,9 @@
     ELSE
        WRITE(9,*) 'surface array not found ', name
-       CALL rs_mpi_io_error( 2 )
+       CALL rd_mpi_io_error( 2 )
     ENDIF
 
       IF ( lo_first_index == 1 )  THEN
-         IF ( debug_level >= 2 .AND. nr_val > 0 )  WRITE(9,*)  'r_surf ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
+         IF ( debug_level >= 2 .AND. nr_val > 0 )  WRITE(9,*)  'r_surf_1 ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
       ELSE
          IF ( debug_level >= 2 .AND. nr_val > 0 )  WRITE(9,*)  'r_surf_next ', TRIM( name ), ' ', &
@@ -1856,4 +1995,7 @@
 #endif
     INTEGER(iwp), OPTIONAL             ::  first_index
+#if defined( __parallel )
+    INTEGER(iwp)                       ::  i
+#endif
     INTEGER(iwp)                       ::  lo_first_index
     INTEGER(KIND=rd_offset_kind)       ::  offset
@@ -1863,5 +2005,5 @@
 #endif
 
-    REAL(wp), INTENT(IN), DIMENSION(:) ::  data
+    REAL(wp), INTENT(IN), DIMENSION(:), TARGET ::  data
 
 
@@ -1880,15 +2022,32 @@
     ENDIF
 #if defined( __parallel )
-    IF ( all_pes_write )  THEN
-       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_surf, 'native', MPI_INFO_NULL, ierr )
-       CALL MPI_FILE_WRITE_ALL( fh, data, nr_val, MPI_REAL, status, ierr )
+    IF ( sm_io%is_sm_active() )  THEN
+       DO  i = 1, nr_val
+          array_1d(i+local_start) = data(i)
+       ENDDO
+       IF ( debug_level >= 2 )  WRITE(9,*) 'w_surf ', TRIM( name ), ' ', SUM( array_1d(local_start+1:local_start+nr_val)),sum(data)
     ELSE
-       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
-       IF ( nr_val > 0 )  THEN
-          disp = array_position + 8 * ( glo_start - 1 )
-          CALL MPI_FILE_SEEK( fh, disp, MPI_SEEK_SET, ierr )
-          CALL MPI_FILE_WRITE( fh, data, nr_val, MPI_REAL, status, ierr )
-       ENDIF
-    ENDIF
+!       array_1d => data                           !kk Did not work in all cases    why???
+       ALLOCATE( array_1d( SIZE( data ) ) )
+       array_1d = data
+    ENDIF
+
+    CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
+    IF ( sm_io%iam_io_pe )  THEN
+       IF ( all_pes_write )  THEN
+          CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_surf, 'native', MPI_INFO_NULL,  &
+                                  ierr )
+          CALL MPI_FILE_WRITE_ALL( fh, array_1d, nr_iope, MPI_REAL, status, ierr )
+       ELSE
+          CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+          IF ( nr_val > 0 )  THEN
+             disp = array_position + 8 * ( glo_start - 1 )
+             CALL MPI_FILE_SEEK( fh, disp, MPI_SEEK_SET, ierr )
+             CALL MPI_FILE_WRITE( fh, array_1d, nr_iope, MPI_REAL, status, ierr )
+          ENDIF
+       ENDIF
+    ENDIF
+    CALL sm_io%sm_node_barrier()
+    IF( .NOT. sm_io%is_sm_active() )  DEALLOCATE( array_1d )
 #else
     CALL posix_lseek( fh, array_position )
@@ -1898,5 +2057,5 @@
 
     IF ( lo_first_index == 1 )  THEN
-       IF ( debug_level >= 2 .AND. nr_val  > 0 )  WRITE(9,*) 'w_surf ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
+       IF ( debug_level >= 2 .AND. nr_val  > 0 )  WRITE(9,*) 'w_surf_1 ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
     ELSE
        IF ( debug_level >= 2 .AND. nr_val  > 0 ) WRITE(9,*) 'w_surf_n ', TRIM( name ), ' ', &
@@ -1955,5 +2114,5 @@
     offset = 0
 
-    IF ( wr_flag )  THEN
+    IF ( wr_flag  .AND.  sm_io%iam_io_pe )  THEN
 
        tgh%nr_int    = header_int_index - 1
@@ -1971,5 +2130,5 @@
 !--    Check for big/little endian format. This check is currently not used, and could be removed
 !--    if we can assume little endian as the default on all machines.
-       CALL rs_mpi_io_check_endian( tgh%endian )
+       CALL rd_mpi_io_check_endian( tgh%endian )
 
 !
@@ -1981,5 +2140,5 @@
 !--    Write header file
        gh_size = storage_size(tgh) / 8
-       IF ( myid == 0 )  THEN
+       IF ( myid == 0 )  THEN   ! myid = 0 always performs I/O, even if I/O is limited to some cores
 #if defined( __parallel )
           CALL MPI_FILE_SEEK( fh, header_position, MPI_SEEK_SET, ierr )
@@ -2060,5 +2219,5 @@
           ENDIF
 
-          IF ( print_header_now )  CALL rs_mpi_io_print_header
+          IF ( print_header_now )  CALL rd_mpi_io_print_header
        ENDIF
 
@@ -2067,13 +2226,25 @@
 !
 !-- Free file types
-    CALL rs_mpi_io_free_filetypes
-
-!
-!-- Close MPI-IO file
-#if defined( __parallel )
-    CALL MPI_FILE_CLOSE( fh, ierr )
+    CALL rd_mpi_io_free_filetypes
+
+!
+!-- Close MPI-IO files
+#if defined( __parallel )
+!
+!-- Restart file has been opened with comm2d
+    IF ( fhs /= -1 )  THEN
+       CALL MPI_FILE_CLOSE( fhs, ierr )
+    ENDIF
+#endif
+
+    IF ( sm_io%iam_io_pe )  THEN
+
+#if defined( __parallel )
+       CALL MPI_FILE_CLOSE( fh, ierr )
 #else
-    CALL posix_close( fh )
-#endif
+       CALL posix_close( fh )
+#endif
+
+    ENDIF
 
     mb_processed = array_position / ( 1024.0_dp * 1024.0_dp )
@@ -2090,4 +2261,8 @@
 !> called. A main feature of this subroutine is computing the global start indices of the 1d- and
 !> 2d- surface arrays.
+!> Even if I/O is done by a limited number of cores only, the surface data are read by ALL cores!
+!> Reading them by some cores and then distributing the data would result in complicated code
+!> which is suspectable for errors and overloads the reading subroutine. Since reading of surface
+!> data is not time critical (data size is comparably small), it will be read by all cores.
 !--------------------------------------------------------------------------------------------------!
  SUBROUTINE rd_mpi_io_surface_filetypes( start_index, end_index, data_to_write, global_start )
@@ -2116,6 +2291,22 @@
 #if defined( __parallel )
     CALL MPI_ALLREDUCE( lo_nr_val, all_nr_val, numprocs, MPI_INTEGER, MPI_SUM, comm2d, ierr )
-    IF ( ft_surf /= -1 )  THEN
+    IF ( ft_surf /= -1  .AND.  sm_io%iam_io_pe )  THEN
        CALL MPI_TYPE_FREE( ft_surf, ierr )    ! if set, free last surface filetype
+    ENDIF
+
+    IF ( win_surf /= -1 )  THEN
+       IF ( sm_io%is_sm_active() )  THEN
+          CALL MPI_WIN_FREE( win_surf, ierr )
+       ENDIF
+       win_surf = -1
+    ENDIF
+
+    IF ( sm_io%is_sm_active() .AND. rd_flag )  THEN
+       IF ( fhs == -1 )  THEN
+          CALL MPI_FILE_OPEN( comm2d, TRIM( io_file_name ), MPI_MODE_RDONLY, MPI_INFO_NULL, fhs,   &
+                              ierr )
+       ENDIF
+    ELSE
+       fhs = fh
     ENDIF
 #else
@@ -2146,5 +2337,5 @@
 
 #if defined( __parallel )
-       CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
+       CALL MPI_FILE_SET_VIEW( fhs, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
 #endif
     ENDIF
@@ -2161,4 +2352,30 @@
           global_start = -1
        ENDWHERE
+
+#if defined( __parallel )
+       IF ( sm_io%is_sm_active() )  THEN
+          IF ( sm_io%iam_io_pe )  THEN
+!
+!--          Calculate number of values of all PEs of an I/O group
+             nr_iope = 0
+             DO  i = myid, myid+sm_io%sh_npes-1
+                nr_iope = nr_iope + all_nr_val(i)
+             ENDDO
+          ELSE
+             local_start = 0
+             DO  i = myid-sm_io%sh_rank, myid-1
+                local_start = local_start + all_nr_val(i)
+             ENDDO
+          ENDIF
+!
+!--       Get the size of shared memory window on all PEs
+          CALL MPI_BCAST( nr_iope, 1, MPI_INTEGER, 0, sm_io%comm_shared, ierr )
+          CALL sm_io%sm_allocate_shared( array_1d, 1, MAX( 1, nr_iope ), win_surf )
+       ELSE
+          nr_iope = nr_val
+       ENDIF
+#else
+       nr_iope = nr_val
+#endif
 
 !
@@ -2170,18 +2387,23 @@
        ENDIF
 
-       all_pes_write = ( MINVAL( all_nr_val ) > 0 )
-
-       IF ( all_pes_write )  THEN
-          dims1(1)  = total_number_of_surface_values
-          lize1(1)  = nr_val
-          start1(1) = glo_start-1
-
-#if defined( __parallel )
-          IF ( total_number_of_surface_values > 0 )  THEN
-              CALL MPI_TYPE_CREATE_SUBARRAY( 1, dims1, lize1, start1, MPI_ORDER_FORTRAN, MPI_REAL, ft_surf, ierr )
-              CALL MPI_TYPE_COMMIT( ft_surf, ierr )
+       IF ( sm_io%iam_io_pe )  THEN
+
+          all_pes_write = ( MINVAL( all_nr_val ) > 0 )
+
+          IF ( all_pes_write )  THEN
+             dims1(1)  = total_number_of_surface_values
+             lize1(1)  = nr_iope
+             start1(1) = glo_start-1
+
+#if defined( __parallel )
+             IF ( total_number_of_surface_values > 0 )  THEN
+                 CALL MPI_TYPE_CREATE_SUBARRAY( 1, dims1, lize1, start1, MPI_ORDER_FORTRAN,        &
+                                                MPI_REAL, ft_surf, ierr )
+                 CALL MPI_TYPE_COMMIT( ft_surf, ierr )
+             ENDIF
+#endif
           ENDIF
-#endif
-       ENDIF
+       ENDIF
+
     ENDIF
 
@@ -2196,5 +2418,5 @@
 !> distributed in blocks among processes to a single file that contains the global arrays.
 !--------------------------------------------------------------------------------------------------!
-  SUBROUTINE rs_mpi_io_create_filetypes
+  SUBROUTINE rd_mpi_io_create_filetypes
 
     IMPLICIT NONE
@@ -2208,4 +2430,10 @@
     INTEGER, DIMENSION(3) ::  start3
 
+    TYPE(local_boundaries) ::  save_io_grid  !< temporary variable to store grid settings
+
+
+    IF ( sm_io%is_sm_active() )  THEN
+       save_io_grid = sm_io%io_grid
+    ENDIF
 
 !
@@ -2240,4 +2468,6 @@
        ENDIF
 
+       CALL sm_io%sm_adjust_outer_boundary()
+
     ELSE
 
@@ -2253,4 +2483,20 @@
     ENDIF
 
+    IF ( sm_io%is_sm_active() )  THEN
+#if defined( __parallel )
+       CALL sm_io%sm_allocate_shared( array_2d,  sm_io%io_grid%nxl, sm_io%io_grid%nxr,             &
+                                      sm_io%io_grid%nys, sm_io%io_grid%nyn, win_2dr )
+       CALL sm_io%sm_allocate_shared( array_2di, save_io_grid%nxl, save_io_grid%nxr,               &
+                                      save_io_grid%nys, save_io_grid%nyn, win_2di )
+       CALL sm_io%sm_allocate_shared( array_3d, nzb, nzt+1, sm_io%io_grid%nxl, sm_io%io_grid%nxr,  &
+                                      sm_io%io_grid%nys, sm_io%io_grid%nyn, win_3dr )
+#endif
+    ELSE
+       ALLOCATE( array_2d(lb%nxl:lb%nxr,lb%nys:lb%nyn) )
+       ALLOCATE( array_2di(nxl:nxr,nys:nyn) )
+       ALLOCATE( array_3d(nzb:nzt+1,lb%nxl:lb%nxr,lb%nys:lb%nyn) )
+       sm_io%io_grid = lb
+    ENDIF
+
 !
 !-- Create filetype for 2d-REAL array with ghost layers around the total domain
@@ -2258,13 +2504,16 @@
     dims2(2)  = lb%ny + 1
 
-    lize2(1)  = lb%nnx
-    lize2(2)  = lb%nny
-
-    start2(1) = lb%nxl
-    start2(2) = lb%nys
-
-#if defined( __parallel )
-    CALL MPI_TYPE_CREATE_SUBARRAY( 2, dims2, lize2, start2, MPI_ORDER_FORTRAN, MPI_REAL, ft_2d, ierr )
-    CALL MPI_TYPE_COMMIT( ft_2d, ierr )
+    lize2(1)  = sm_io%io_grid%nnx
+    lize2(2)  = sm_io%io_grid%nny
+
+    start2(1) = sm_io%io_grid%nxl
+    start2(2) = sm_io%io_grid%nys
+
+#if defined( __parallel )
+    IF ( sm_io%iam_io_pe )  THEN
+       CALL MPI_TYPE_CREATE_SUBARRAY( 2, dims2, lize2, start2, MPI_ORDER_FORTRAN, MPI_REAL,        &
+                                      ft_2d, ierr )
+       CALL MPI_TYPE_COMMIT( ft_2d, ierr )
+    ENDIF
 #endif
 !
@@ -2273,13 +2522,28 @@
     dims2(2)  = ny + 1
 
-    lize2(1)  = nnx
-    lize2(2)  = nny
-
-    start2(1) = nxl
-    start2(2) = nys
-
-#if defined( __parallel )
-    CALL MPI_TYPE_CREATE_SUBARRAY( 2, dims2, lize2, start2, MPI_ORDER_FORTRAN, MPI_INTEGER, ft_2di_nb, ierr )
-    CALL MPI_TYPE_COMMIT( ft_2di_nb, ierr )
+    IF ( sm_io%is_sm_active() )  THEN
+
+       lize2(1)  = save_io_grid%nnx
+       lize2(2)  = save_io_grid%nny
+
+       start2(1) = save_io_grid%nxl
+       start2(2) = save_io_grid%nys
+
+    ELSE
+
+       lize2(1)  = nnx
+       lize2(2)  = nny
+
+       start2(1) = nxl
+       start2(2) = nys
+
+    ENDIF
+
+#if defined( __parallel )
+    IF ( sm_io%iam_io_pe )  THEN
+       CALL MPI_TYPE_CREATE_SUBARRAY( 2, dims2, lize2, start2, MPI_ORDER_FORTRAN, MPI_INTEGER,     &
+                                      ft_2di_nb, ierr )
+       CALL MPI_TYPE_COMMIT( ft_2di_nb, ierr )
+    ENDIF
 #endif
 !
@@ -2290,17 +2554,20 @@
 
     lize3(1)  = dims3(1)
-    lize3(2)  = lb%nnx
-    lize3(3)  = lb%nny
+    lize3(2)  = sm_io%io_grid%nnx
+    lize3(3)  = sm_io%io_grid%nny
 
     start3(1) = nzb
-    start3(2) = lb%nxl
-    start3(3) = lb%nys
-
-#if defined( __parallel )
-    CALL MPI_TYPE_CREATE_SUBARRAY( 3, dims3, lize3, start3, MPI_ORDER_FORTRAN, MPI_REAL, ft_3d, ierr )
-    CALL MPI_TYPE_COMMIT( ft_3d, ierr )
-#endif
-
- END SUBROUTINE rs_mpi_io_create_filetypes
+    start3(2) = sm_io%io_grid%nxl
+    start3(3) = sm_io%io_grid%nys
+
+#if defined( __parallel )
+    IF ( sm_io%iam_io_pe )  THEN
+       CALL MPI_TYPE_CREATE_SUBARRAY( 3, dims3, lize3, start3, MPI_ORDER_FORTRAN, MPI_REAL, ft_3d, &
+                                      ierr )
+       CALL MPI_TYPE_COMMIT( ft_3d, ierr )
+    ENDIF
+#endif
+
+ END SUBROUTINE rd_mpi_io_create_filetypes
 
 
@@ -2325,4 +2592,13 @@
     INTEGER, DIMENSION(3) ::  start3
 
+
+    IF ( sm_io%is_sm_active() )  THEN
+       CALL sm_io%sm_allocate_shared( array_3d_soil, nzb_soil, nzt_soil, sm_io%io_grid%nxl,        &
+                                      sm_io%io_grid%nxr, sm_io%io_grid%nys, sm_io%io_grid%nyn,     &
+                                      win_3ds )
+    ELSE
+       ALLOCATE( array_3d_soil(nzb_soil:nzt_soil,lb%nxl:lb%nxr,lb%nys:lb%nyn) )
+       sm_io%io_grid = lb
+    ENDIF
 
 !
@@ -2333,14 +2609,16 @@
 
     lize3(1)  = dims3(1)
-    lize3(2)  = lb%nnx
-    lize3(3)  = lb%nny
+    lize3(2)  = sm_io%io_grid%nnx
+    lize3(3)  = sm_io%io_grid%nny
 
     start3(1) = nzb_soil
-    start3(2) = lb%nxl
-    start3(3) = lb%nys
-
-    CALL MPI_TYPE_CREATE_SUBARRAY( 3, dims3, lize3, start3, MPI_ORDER_FORTRAN, MPI_REAL,           &
-                                   ft_3dsoil, ierr )
-    CALL MPI_TYPE_COMMIT( ft_3dsoil, ierr )
+    start3(2) = sm_io%io_grid%nxl
+    start3(3) = sm_io%io_grid%nys
+
+    IF ( sm_io%iam_io_pe )  THEN
+       CALL MPI_TYPE_CREATE_SUBARRAY( 3, dims3, lize3, start3, MPI_ORDER_FORTRAN, MPI_REAL,        &
+                                      ft_3dsoil, ierr )
+       CALL MPI_TYPE_COMMIT( ft_3dsoil, ierr )
+    ENDIF
 
  END SUBROUTINE rd_mpi_io_create_filetypes_3dsoil
@@ -2354,5 +2632,5 @@
 !> Free all file types that have been created for MPI-IO.
 !--------------------------------------------------------------------------------------------------!
- SUBROUTINE rs_mpi_io_free_filetypes
+ SUBROUTINE rd_mpi_io_free_filetypes
 
     IMPLICIT NONE
@@ -2361,16 +2639,37 @@
 #if defined( __parallel )
     IF ( filetypes_created )  THEN
-       CALL MPI_TYPE_FREE( ft_2d, ierr )
-       CALL MPI_TYPE_FREE( ft_2di_nb, ierr )
-       CALL MPI_TYPE_FREE( ft_3d, ierr )
+
+       IF ( sm_io%iam_io_pe )  THEN
+          CALL MPI_TYPE_FREE( ft_2d, ierr )
+          CALL MPI_TYPE_FREE( ft_2di_nb, ierr )
+          CALL MPI_TYPE_FREE( ft_3d, ierr )
+       ENDIF
+
+       IF ( sm_io%is_sm_active() )  THEN
+          CALL sm_io%sm_free_shared( win_2dr )
+          CALL sm_io%sm_free_shared( win_2di )
+          CALL sm_io%sm_free_shared( win_3dr )
+       ELSE
+          DEALLOCATE( array_2d, array_2di, array_3d )
+       ENDIF
+
     ENDIF
 !
 !-- Free last surface filetype
-    IF ( ft_surf /= -1 )  THEN  
+    IF ( sm_io%iam_io_pe .AND. ft_surf /= -1 )  THEN
        CALL MPI_TYPE_FREE( ft_surf, ierr )
     ENDIF
-#endif
-
- END SUBROUTINE rs_mpi_io_free_filetypes
+
+    IF ( sm_io%is_sm_active() .AND.  win_surf /= -1 )  THEN
+       CALL sm_io%sm_free_shared( win_surf )
+    ENDIF
+
+    ft_surf  = -1
+    win_surf = -1
+#else
+    DEALLOCATE( array_2d, array_2di, array_3d )
+#endif
+
+ END SUBROUTINE rd_mpi_io_free_filetypes
 
 
@@ -2381,5 +2680,5 @@
 !> Print the restart data file header (MPI-IO format) for debugging.
 !--------------------------------------------------------------------------------------------------!
- SUBROUTINE rs_mpi_io_print_header
+ SUBROUTINE rd_mpi_io_print_header
 
     IMPLICIT NONE
@@ -2422,5 +2721,5 @@
     print_header_now = .FALSE.
 
- END SUBROUTINE rs_mpi_io_print_header
+ END SUBROUTINE rd_mpi_io_print_header
 
 
@@ -2431,5 +2730,5 @@
 !> Print error messages for reading/writing restart data with MPI-IO
 !--------------------------------------------------------------------------------------------------!
- SUBROUTINE rs_mpi_io_error( error_number )
+ SUBROUTINE rd_mpi_io_error( error_number )
 
     IMPLICIT NONE
@@ -2454,5 +2753,5 @@
              WRITE(6,*)  'posix IO: ERROR Opening Restart File'
           CASE DEFAULT
-             WRITE(6,*)  'rs_mpi_io_error: illegal error number: ',error_number
+             WRITE(6,*)  'rd_mpi_io_error: illegal error number: ',error_number
 
        END SELECT
@@ -2460,11 +2759,13 @@
     ENDIF
 #if defined( __parallel )
-    CALL MPI_BARRIER( comm2d, ierr )
-    CALL MPI_ABORT( comm2d, 1, ierr )
+    IF ( .NOT. sm_io%iam_io_pe )  RETURN
+
+    CALL MPI_BARRIER( comm_io, ierr )
+    CALL MPI_ABORT( comm_io, 1, ierr )
 #else
     CALL ABORT
 #endif
 
- END SUBROUTINE rs_mpi_io_error
+ END SUBROUTINE rd_mpi_io_error
 
 
@@ -2477,5 +2778,5 @@
 !> the first 4 bytes of the pointer are equal 1 (little endian) or not.
 !--------------------------------------------------------------------------------------------------!
- SUBROUTINE rs_mpi_io_check_endian( i_endian )
+ SUBROUTINE rd_mpi_io_check_endian( i_endian )
 
     IMPLICIT NONE
@@ -2502,5 +2803,5 @@
     ENDIF
 
- END SUBROUTINE rs_mpi_io_check_endian
+ END SUBROUTINE rd_mpi_io_check_endian
 
  END MODULE restart_data_mpi_io_mod
Index: /palm/trunk/SOURCE/shared_memory_io_mod.f90
===================================================================
--- /palm/trunk/SOURCE/shared_memory_io_mod.f90	(revision 4534)
+++ /palm/trunk/SOURCE/shared_memory_io_mod.f90	(revision 4534)
@@ -0,0 +1,671 @@
+!> @file shared_memory_io_mod.f90
+!------------------------------------------------------------------------------!
+! This file is part of PALM.
+!
+! PALM is free software: you can redistribute it and/or modify it under the
+! terms of the GNU General Public License as published by the Free Software
+! Foundation, either version 3 of the License, or (at your option) any later
+! version.
+!
+! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
+! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
+! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
+!
+! You should have received a copy of the GNU General Public License along with
+! PALM. If not, see <http://www.gnu.org/licenses/>.
+!
+! Copyright 1997-2020 Leibniz Universitaet Hannover
+!------------------------------------------------------------------------------!
+!
+! Current revisions:
+! -----------------
+! 
+! 
+! Former revisions:
+! $Id$
+!
+! Initial version (Klaus Ketelsen)
+!
+! 
+!
+! Description:
+! ------------
+!> handle MPI-IO or NetCDF-IO shared memory arrays.
+!> This module performs the organization of new communicators, adapted PE-grids
+!> and allocation of shared memory arrays. The IO itself is not done here.
+!------------------------------------------------------------------------------!
+ MODULE shared_memory_io_mod
+
+#if defined( __parallel )
+#if defined( __mpifh )
+    INCLUDE "mpif.h"
+#else
+    USE MPI
+#endif
+#endif
+
+    USE, INTRINSIC ::  ISO_C_BINDING
+
+    USE control_parameters,                                                                        &
+        ONLY:  maximum_grid_level, mg_switch_to_pe0_level, message_string
+
+    USE indices,                                                                                   &
+        ONLY: nbgp, nnx, nny, nnz, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg, nzb, nzt
+
+    USE kinds,                                                                                     &
+        ONLY: wp, iwp
+
+    USE transpose_indices,                                                                         &
+        ONLY:  nys_x, nyn_x, nys_z, nyn_z, nxl_z, nxr_z
+
+    USE pegrid,                                                                                    &
+        ONLY: comm1dx, comm1dy, comm2d, ierr, myid, myidx, myidy, npex, npey, numprocs, pdims,     &
+              pleft, pnorth, pright, psouth, sendrecvcount_xy
+#if defined( __parallel )
+    USE pegrid,                                                                                    &
+        ONLY: pcoord, reorder
+#endif
+
+    IMPLICIT NONE
+
+    PRIVATE
+
+    SAVE
+
+!
+!-- Type to store grid information
+    TYPE, PUBLIC ::  local_boundaries
+
+       INTEGER(iwp) ::  nxl
+       INTEGER(iwp) ::  nxr
+       INTEGER(iwp) ::  nys
+       INTEGER(iwp) ::  nyn
+       INTEGER(iwp) ::  nnx
+       INTEGER(iwp) ::  nny
+       INTEGER(iwp) ::  nx
+       INTEGER(iwp) ::  ny
+
+    END TYPE local_boundaries
+
+!
+!-- Class definition for shared memory instances.
+!-- For every use of shared memory IO, one instance of this class is created.
+    TYPE, PUBLIC ::  sm_class
+
+       INTEGER(iwp)         :: nr_io_pe_per_node = 2   !< typical configuration, 2 sockets per node
+       LOGICAL              :: no_shared_Memory_in_this_run
+!
+!--    Variables for the shared memory communicator
+       INTEGER(iwp), PUBLIC ::  comm_shared    !< Communicator for processes with shared array
+       INTEGER(iwp), PUBLIC ::  sh_npes
+       INTEGER(iwp), PUBLIC ::  sh_rank
+
+       LOGICAL, PUBLIC ::  iam_io_pe = .TRUE.  !< This PE is an IO-PE
+!
+!--    Variables for the I/O virtual grid
+       INTEGER(iwp), PUBLIC ::  comm_io        !< Communicator for all IO processes
+       INTEGER(iwp), PUBLIC ::  io_npes
+       INTEGER(iwp), PUBLIC ::  io_rank
+
+       TYPE( local_boundaries ), PUBLIC ::  io_grid
+
+!
+!--    Variables for the node local communicator
+       INTEGER(iwp)         ::  comm_node      !< Communicator for all processes of current node
+       INTEGER(iwp)         ::  io_pe_global_rank
+       INTEGER(iwp)         ::  n_npes
+       INTEGER(iwp)         ::  n_rank
+
+       CONTAINS
+
+          PRIVATE
+
+          PROCEDURE, PASS(this), PUBLIC ::  is_sm_active
+          PROCEDURE, PASS(this), PUBLIC ::  sm_adjust_outer_boundary
+          PROCEDURE, PASS(this), PUBLIC ::  sm_free_shared
+          PROCEDURE, PASS(this), PUBLIC ::  sm_init_comm
+          PROCEDURE, PASS(this), PUBLIC ::  sm_node_barrier
+#if defined( __parallel )
+          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_1d
+          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_2d
+          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_2di
+          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_3d
+
+          GENERIC, PUBLIC ::  sm_allocate_shared =>  sm_allocate_shared_1d, sm_allocate_shared_2d, &
+                                                     sm_allocate_shared_2di, sm_allocate_shared_3d
+#endif
+    END TYPE sm_class
+
+
+CONTAINS
+
+
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Create the communicator for shared memory groups and IO-PEs.
+!> Setup the grid for shared memory IO.
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE sm_init_comm( this, sm_active )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(INOUT) ::  this  !< pointer to access internal variables of this call
+
+#if defined( __parallel )
+    INTEGER             :: color
+    INTEGER             :: max_n_npes  !< Maximum number of PEs/node
+#endif
+
+    LOGICAL, INTENT(IN) ::  sm_active  !< Flag to activate shared-memory IO
+
+
+    this%no_shared_memory_in_this_run = .NOT. sm_active
+
+    IF ( this%no_shared_memory_in_this_run )  THEN
+       this%iam_io_pe = .TRUE.
+       RETURN
+    ENDIF
+
+#if defined( __parallel )
+!
+!-- Determine, how many MPI threads are running on a node
+    this%iam_io_pe = .FALSE.
+    CALL MPI_COMM_SPLIT_TYPE( comm2d, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL, this%comm_node, ierr )
+    CALL MPI_COMM_SIZE( this%comm_node, this%n_npes, ierr )
+    CALL MPI_COMM_RANK( this%comm_node, this%n_rank, ierr )
+
+    CALL MPI_ALLREDUCE( this%n_npes, max_n_npes, 1, MPI_INTEGER, MPI_MAX, comm2d, ierr )
+!
+!-- Decide, if the configuration can run with shared-memory IO
+    IF ( max_n_npes > 64 )  THEN
+!
+!--    Special configuration on the HLRN-IV system with 4 shared memory blocks/node
+       this%nr_io_pe_per_node = 4
+
+    ELSEIF ( max_n_npes <= 32 )  THEN
+!
+!--    No shared memory IO with less than 32 threads/node
+       this%no_shared_memory_in_this_run = .TRUE.
+       this%iam_io_pe = .TRUE.
+       RETURN
+    ENDIF
+
+!
+!-- No shared memory IO with small setups
+    IF ( nx < 24  .OR.  ny < 24 )  THEN
+       this%no_shared_memory_in_this_run = .TRUE.
+       this%iam_io_pe = .TRUE.
+       RETURN
+    ENDIF
+
+!
+!-- Divide a node into shared memory groups, depending on the virtual x-y grid
+    CALL compute_color( color )
+!
+!-- If no shared memory IO possible, nothing is left to be done here.
+    IF ( this%no_shared_memory_in_this_run )  RETURN
+
+!
+!-- Setup the shared memory area
+    CALL MPI_COMM_SPLIT( this%comm_node, color, 0, this%comm_shared, ierr )
+    CALL MPI_COMM_SIZE( this%comm_shared, this%sh_npes, ierr )
+    CALL MPI_COMM_RANK( this%comm_shared, this%sh_rank, ierr )
+
+!
+!-- Setup the communicator across the nodes depending on the shared memory rank.
+!-- All threads with shared memory rank 0 will be I/O threads.
+    color = this%sh_rank
+    CALL MPI_COMM_SPLIT( comm2d, color, 0, this%comm_io, ierr )
+
+    IF ( this%comm_io /= MPI_COMM_NULL )  THEN
+       CALL MPI_COMM_SIZE( this%comm_io, this%io_npes, ierr )
+       CALL MPI_COMM_RANK( this%comm_io, this%io_rank, ierr )
+    ELSE
+       this%io_npes = -1
+       this%io_rank = -1
+    ENDIF
+
+    IF ( this%sh_rank == 0 )  THEN
+       this%iam_io_pe = .TRUE.
+       this%io_pe_global_rank = myid
+    ENDIF
+    CALL MPI_BCAST( this%io_pe_global_rank, 1, MPI_INTEGER, 0, this%comm_shared, ierr )
+
+#else
+    this%iam_io_pe = .TRUE.
+#endif
+
+#if defined( __parallel )
+ CONTAINS
+
+ SUBROUTINE compute_color( color )
+
+    IMPLICIT NONE
+
+    INTEGER(iwp), INTENT(OUT) :: color
+
+    INTEGER(iwp) ::  group_start
+    INTEGER(iwp) ::  n
+    INTEGER(iwp) ::  my_color
+    INTEGER(iwp) ::  pe
+    INTEGER(iwp) ::  sh_group_size
+
+    INTEGER(iwp), DIMENSION(4,0:this%n_npes-1) ::  local_dim_s
+    INTEGER(iwp), DIMENSION(4,0:this%n_npes-1) ::  local_dim_r
+
+    TYPE(local_boundaries), DIMENSION(32) ::  node_grid
+
+#if defined( __intel_compiler )
+    IF ( numprocs < this%n_npes )  THEN  !  NO shared memory IO on one node jobs
+       this%no_shared_Memory_in_this_run = .TRUE.
+       RETURN
+    ENDIF
+#endif
+
+    local_dim_s = 0
+    local_dim_s(1,this%n_rank) = nxl
+    local_dim_s(2,this%n_rank) = nxr
+    local_dim_s(3,this%n_rank) = nys
+    local_dim_s(4,this%n_rank) = nyn
+
+    node_grid%nyn = -1
+!
+!-- Distribute the x-y layout of all cores of a node to all node processes
+    CALL MPI_ALLREDUCE( local_dim_s, local_dim_r, SIZE( local_dim_s ), MPI_INTEGER, MPI_SUM,       &
+                        this%comm_node, ierr )
+    sh_group_size = ( max_n_npes + this%nr_io_pe_per_node - 1 ) / this%nr_io_pe_per_node
+
+    pe       = 0
+    my_color = 1  ! color is used to split the shared memory communicator into a communicator for
+                  ! io groups
+    group_start = pe
+    node_grid(my_color)%nxl = local_dim_r(1,group_start)
+    node_grid(my_color)%nxr = local_dim_r(2,group_start)
+    node_grid(my_color)%nys = local_dim_r(3,group_start)
+
+    DO  n = 1, this%n_npes-1
+
+       pe =  n
+       IF ( n > 0  .AND.  MOD( n,sh_group_size ) == 0 )  THEN
+!
+!--       If group boundary, start new IO group
+          node_grid(my_color)%nyn = local_dim_r(4,pe-1)
+          my_color = my_color + 1
+          group_start = pe
+          node_grid(my_color)%nxl = local_dim_r(1,group_start)
+          node_grid(my_color)%nxr = local_dim_r(2,group_start)
+          node_grid(my_color)%nys = local_dim_r(3,group_start)
+
+       ELSEIF ( local_dim_r(1,pe) /= node_grid(my_color)%nxl )  THEN
+!
+!--       If nxl changes, start new IO group
+          node_grid(my_color)%nyn = local_dim_r(4,pe-1)
+          my_color = my_color+1
+          group_start = pe
+          node_grid(my_color)%nxl = local_dim_r(1,group_start)
+          node_grid(my_color)%nxr = local_dim_r(2,group_start)
+          node_grid(my_color)%nys = local_dim_r(3,group_start)
+       ENDIF
+!
+!--    Save values for local PE
+       IF ( this%n_rank == pe )  THEN                                 !
+          color = my_color
+       ENDIF
+       IF ( n == this%n_npes-1 )  node_grid(my_color)%nyn = local_dim_r(4,pe)
+
+    ENDDO
+
+    IF ( this%n_rank == 0 )  THEN
+       color = 1
+    ENDIF
+
+ END SUBROUTINE compute_color
+#endif
+
+ END SUBROUTINE sm_init_comm
+
+
+
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Function to return if shared Memory IO is active.
+!--------------------------------------------------------------------------------------------------!
+ FUNCTION is_sm_active( this) RESULT( ac )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(inout) ::  this
+
+    LOGICAL ::  ac
+
+    ac = .NOT. this%no_shared_memory_in_this_run
+
+ END FUNCTION is_sm_active
+
+
+#if defined( __parallel )
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Allocate shared 1d-REAL array on ALL threads
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE sm_allocate_shared_1d( this, p1, d1, d2, win )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(inout)  ::  this
+
+    INTEGER(iwp)                    ::  disp_unit
+    INTEGER(iwp), INTENT(IN)        ::  d1
+    INTEGER(iwp), INTENT(IN)        ::  d2
+    INTEGER(iwp), SAVE              ::  pe_from = 0
+    INTEGER(KIND=MPI_ADDRESS_KIND)  ::  rem_size
+    INTEGER(iwp), INTENT(OUT)       ::  win
+    INTEGER(KIND=MPI_ADDRESS_KIND)  ::  wsize
+
+    INTEGER, DIMENSION(1)           ::  buf_shape
+
+    REAL(wp), DIMENSION(:), POINTER ::  buf
+    REAL(wp), DIMENSION(:), POINTER ::  p1
+
+    TYPE(C_PTR), SAVE               ::  base_ptr
+    TYPE(C_PTR), SAVE               ::  rem_ptr
+
+
+    IF ( this%no_shared_memory_in_this_run )  RETURN
+!
+!-- Allocate shared memory on node rank 0 threads.
+    IF ( this%sh_rank == pe_from )  THEN
+       wsize = d2 - d1 + 1
+    ELSE
+       wsize = 1
+    ENDIF
+    wsize = wsize * 8  ! Please note, size is always in bytes, independently of the displacement
+                       ! unit
+
+    CALL MPI_WIN_ALLOCATE_SHARED( wsize, 8, MPI_INFO_NULL, this%comm_shared,base_ptr, win, ierr )
+!
+!-- Get C-pointer of the memory located on node-rank pe_from (sh_rank == pe_from)
+    CALL MPI_WIN_SHARED_QUERY( win, pe_from, rem_size, disp_unit, rem_ptr, ierr )
+!
+!-- Convert C- to Fortran-pointer
+    buf_shape(1) = d2 - d1 + 1
+    CALL C_F_POINTER( rem_ptr, buf, buf_shape )
+    p1(d1:) => buf
+!
+!-- Allocate shared memory in round robin on all PEs of a node.
+    pe_from = MOD( pe_from, this%sh_npes )
+
+ END SUBROUTINE sm_allocate_shared_1d
+
+
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Allocate shared 2d-REAL array on ALL threads
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE sm_allocate_shared_2d( this, p2, n_nxlg, n_nxrg, n_nysg, n_nyng, win )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(INOUT)    ::  this
+
+    INTEGER(iwp)                      ::  disp_unit
+    INTEGER(iwp), INTENT(IN)          ::  n_nxlg
+    INTEGER(iwp), INTENT(IN)          ::  n_nxrg
+    INTEGER(iwp), INTENT(IN)          ::  n_nyng
+    INTEGER(iwp), INTENT(IN)          ::  n_nysg
+    INTEGER(iwp), SAVE                ::  pe_from = 0
+    INTEGER(KIND=MPI_ADDRESS_KIND)    ::  rem_size
+    INTEGER(iwp), INTENT(OUT)         ::  win
+    INTEGER(KIND=MPI_ADDRESS_KIND)    ::  wsize
+
+    INTEGER(iwp), DIMENSION(2)        ::  buf_shape
+
+    REAL(wp), DIMENSION(:,:), POINTER ::  buf
+    REAL(wp), DIMENSION(:,:), POINTER ::  p2
+
+    TYPE(C_PTR),SAVE                  ::  base_ptr
+    TYPE(C_PTR),SAVE                  ::  rem_ptr
+
+
+    IF ( this%no_shared_memory_in_this_run )  RETURN
+!
+!-- Allocate shared memory on node rank 0 threads.
+    IF ( this%sh_rank == pe_from )  THEN
+       wsize = ( n_nyng - n_nysg + 1 ) * ( n_nxrg - n_nxlg + 1 )
+    ELSE
+       wsize = 1
+    ENDIF
+
+    wsize = wsize * 8  ! Please note, size is always in bytes, independently of the displacement
+                       ! unit
+
+    CALL MPI_WIN_ALLOCATE_SHARED( wsize, 8, MPI_INFO_NULL, this%comm_shared, base_ptr, win, ierr )
+!
+!-- Get C-pointer of the memory located on node-rank pe_from (sh_rank == pe_from)
+    CALL MPI_WIN_SHARED_QUERY( win, pe_from, rem_size, disp_unit, rem_ptr, ierr )
+!
+!-- Convert C- to Fortran-pointer
+    buf_shape(2) = n_nyng - n_nysg + 1
+    buf_shape(1) = n_nxrg - n_nxlg + 1
+    CALL C_F_POINTER( rem_ptr, buf, buf_shape )
+    p2(n_nxlg:, n_nysg:) => buf
+!
+!-- Allocate shared memory in round robin on all PEs of a node.
+    pe_from = MOD( pe_from, this%sh_npes )
+
+ END SUBROUTINE sm_allocate_shared_2d
+
+
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Allocate shared 2d-INTEGER array on ALL threads
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE sm_allocate_shared_2di( this, p2i, n_nxlg, n_nxrg, n_nysg, n_nyng, win )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(inout)      :: this
+
+    INTEGER(iwp)                          ::  disp_unit
+    INTEGER(iwp), INTENT(IN)              ::  n_nxlg
+    INTEGER(iwp), INTENT(IN)              ::  n_nxrg
+    INTEGER(iwp), INTENT(IN)              ::  n_nyng
+    INTEGER(iwp), INTENT(IN)              ::  n_nysg
+    INTEGER(iwp), SAVE                    ::  pe_from = 0
+    INTEGER(kind=MPI_ADDRESS_KIND)        ::  rem_size
+    INTEGER(iwp), INTENT(OUT)             ::  win
+    INTEGER(kind=MPI_ADDRESS_KIND)        ::  wsize
+
+    INTEGER(iwp), DIMENSION(2)            ::  buf_shape
+
+    INTEGER(iwp), DIMENSION(:,:), POINTER ::  buf
+    INTEGER(iwp), DIMENSION(:,:), POINTER ::  p2i
+
+    TYPE(C_PTR),SAVE                      ::  base_ptr
+    TYPE(C_PTR),SAVE                      ::  rem_ptr
+
+
+    IF ( this%no_shared_memory_in_this_run )  RETURN
+!
+!-- Allocate shared memory on node rank 0 threads.
+    IF ( this%sh_rank == pe_from )  THEN
+       wsize = ( n_nyng - n_nysg + 1 ) * ( n_nxrg - n_nxlg + 1 )
+    ELSE
+       wsize = 1
+    ENDIF
+
+    wsize = wsize * 4  ! Please note, size is always in bytes, independently of the displacement
+                       ! unit
+
+    CALL MPI_WIN_ALLOCATE_SHARED( wsize, 4, MPI_INFO_NULL, this%comm_shared, base_ptr, win, ierr )
+!
+!-- Get C-pointer of the memory located on node-rank pe_from (sh_rank == pe_from)
+    CALL MPI_WIN_SHARED_QUERY( win, pe_from, rem_size, disp_unit, rem_ptr, ierr )
+!
+!-- Convert C- to Fortran-pointer
+    buf_shape(2) = n_nyng - n_nysg + 1
+    buf_shape(1) = n_nxrg - n_nxlg + 1
+    CALL C_F_POINTER( rem_ptr, buf, buf_shape )
+    p2i(n_nxlg:, n_nysg:) => buf
+!
+!-- Allocate shared memory in round robin on all PEs of a node.
+    pe_from = MOD( pe_from, this%sh_npes )
+
+ END SUBROUTINE sm_allocate_shared_2di
+
+
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Allocate shared 3d-REAL array on ALL threads
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE sm_allocate_shared_3d( this, p3, d1s, d1e, d2s, d2e, d3s, d3e, win )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(inout)      ::  this
+
+    INTEGER                             ::  disp_unit
+    INTEGER, INTENT(IN)                 ::  d1e
+    INTEGER, INTENT(IN)                 ::  d1s
+    INTEGER, INTENT(IN)                 ::  d2e
+    INTEGER, INTENT(IN)                 ::  d2s
+    INTEGER, INTENT(IN)                 ::  d3e
+    INTEGER, INTENT(IN)                 ::  d3s
+    INTEGER, SAVE                       ::  pe_from = 0
+    INTEGER(KIND=MPI_ADDRESS_KIND)      ::  rem_size
+    INTEGER, INTENT(OUT)                ::  win
+    INTEGER(KIND=MPI_ADDRESS_KIND)      ::  wsize
+
+    INTEGER, DIMENSION(3)               ::  buf_shape
+
+    REAL(wp), DIMENSION(:,:,:), POINTER ::  buf
+    REAL(wp), DIMENSION(:,:,:), POINTER ::  p3
+
+    TYPE(C_PTR), SAVE                   ::  base_ptr
+    TYPE(C_PTR), SAVE                   ::  rem_ptr
+
+
+    IF ( this%no_shared_memory_in_this_run )  RETURN
+!
+!-- Allocate shared memory on node rank 0 threads.
+    IF ( this%sh_rank == pe_from )  THEN
+       wsize = ( d3e - d3s + 1 ) * ( d2e - d2s + 1 ) * ( d1e - d1s + 1 )
+    ELSE
+       wsize = 1
+    ENDIF
+
+    wsize = wsize * 8  ! Please note, size is always in bytes, independently of the displacement
+                       ! unit
+
+    CALL MPI_WIN_ALLOCATE_SHARED( wsize, 8, MPI_INFO_NULL, this%comm_shared, base_ptr, win, ierr )
+!
+!-- Get C-pointer of the memory located on node-rank pe_from (sh_rank == pe_from)
+    CALL MPI_WIN_SHARED_QUERY( win, pe_from, rem_size, disp_unit, rem_ptr, ierr )
+!
+!-- Convert C- to Fortran-pointer
+    buf_shape(3) = d3e - d3s + 1
+    buf_shape(2) = d2e - d2s + 1
+    buf_shape(1) = d1e - d1s + 1
+    CALL C_F_POINTER( rem_ptr, buf, buf_shape )
+    p3(d1s:,d2s:,d3s:) => buf
+!
+!-- Allocate shared memory in round robin on all PEs of a node.
+    pe_from = MOD( pe_from, this%sh_npes )
+
+ END SUBROUTINE sm_allocate_shared_3d
+#endif
+
+
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> ???
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE sm_adjust_outer_boundary( this )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(inout) ::  this
+
+
+    IF ( this%no_shared_memory_in_this_run )  RETURN
+
+    IF ( this%io_grid%nxl == 0 )  THEN
+       this%io_grid%nxl = this%io_grid%nxl - nbgp
+       this%io_grid%nnx = this%io_grid%nnx + nbgp
+    ENDIF
+
+    IF ( this%io_grid%nxr == nx  .OR.  npex == -1 )  THEN   ! npex == -1 if -D__parallel not set
+       this%io_grid%nxr = this%io_grid%nxr + nbgp
+       this%io_grid%nnx = this%io_grid%nnx + nbgp
+    ENDIF
+
+    IF ( this%io_grid%nys == 0 )  THEN
+       this%io_grid%nys = this%io_grid%nys - nbgp
+       this%io_grid%nny = this%io_grid%nny + nbgp
+    ENDIF
+
+    IF ( this%io_grid%nyn == ny .OR.  npey == -1 )  THEN   ! npey == -1 if -D__parallel not set
+       this%io_grid%nyn = this%io_grid%nyn + nbgp
+       this%io_grid%nny = this%io_grid%nny + nbgp
+    ENDIF
+
+    this%io_grid%nxl = this%io_grid%nxl + nbgp
+    this%io_grid%nxr = this%io_grid%nxr + nbgp
+    this%io_grid%nys = this%io_grid%nys + nbgp
+    this%io_grid%nyn = this%io_grid%nyn + nbgp
+    this%io_grid%nnx = this%io_grid%nnx
+    this%io_grid%nny = this%io_grid%nny
+
+ END SUBROUTINE sm_adjust_outer_boundary
+
+
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> Deallocate shared aray and free related window.
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE sm_free_shared( this, win )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(inout) ::  this
+
+    INTEGER(iwp), INTENT(INOUT)    ::  win
+
+    IF ( this%no_shared_memory_in_this_run  .OR.  win == -1234567890 )  RETURN
+                     ! win is used just to avoid compile errors because of unused arguments
+#if defined( __parallel )
+    CALL MPI_WIN_FREE( win, ierr )
+#endif
+
+ END SUBROUTINE sm_free_shared
+
+
+!--------------------------------------------------------------------------------------------------!
+! Description:
+! ------------
+!> ...
+!--------------------------------------------------------------------------------------------------!
+ SUBROUTINE sm_node_barrier( this )
+
+    IMPLICIT NONE
+
+    CLASS(sm_class), INTENT(inout)     :: this
+
+
+    IF ( this%no_shared_memory_in_this_run )  RETURN
+
+#if defined( __parallel )
+    CALL MPI_BARRIER( this%comm_shared, ierr )
+#endif
+
+ END SUBROUTINE sm_node_barrier
+
+ END MODULE shared_memory_io_mod
