Changeset 4893


Ignore:
Timestamp:
Mar 2, 2021 4:39:14 PM (3 years ago)
Author:
raasch
Message:

revised output of surface data via MPI-IO for better performance

Location:
palm/trunk/SOURCE
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/init_pegrid.f90

    r4848 r4893  
    358358    CALL MPI_COMM_RANK( comm1dy, myidy, ierr )
    359359
    360 
    361360!
    362361!-- Calculate array bounds along x-direction for every PE.
  • palm/trunk/SOURCE/land_surface_model_mod.f90

    r4876 r4893  
    2424! -----------------
    2525! $Id$
     26! revised output of surface data via MPI-IO for better performance
     27!
     28! 4876 2021-02-17 12:27:36Z raasch
    2629! bugfix for instantaneous c_liq output
    2730!
     
    62496252    INTEGER(iwp)      ::  l      !< index variable for surface orientation
    62506253
    6251     INTEGER(iwp),DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index  !< index for surface data (MPI-IO)
     6254    INTEGER(iwp),DIMENSION(nys:nyn,nxl:nxr) ::  global_end_index    !< end index for surface data (MPI-IO)
     6255    INTEGER(iwp),DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index  !< start index for surface data (MPI-IO)
    62526256
    62536257    LOGICAL ::  surface_data_to_write  !< switch for MPI-I/O if PE has surface data to write
     
    64356439
    64366440          CALL rd_mpi_io_surface_filetypes( surf_lsm_h(l)%start_index, surf_lsm_h(l)%end_index,    &
    6437                                             surface_data_to_write, global_start_index )
    6438 
    6439           CALL wrd_mpi_io( 'lsm_start_index_h_' // dum,  surf_lsm_h(l)%start_index )
    6440           CALL wrd_mpi_io( 'lsm_end_index_h_' // dum,  surf_lsm_h(l)%end_index )
     6441                                            surface_data_to_write, global_start_index,             &
     6442                                            global_end_index )
     6443
    64416444          CALL wrd_mpi_io( 'lsm_global_start_index_h_' // dum, global_start_index )
     6445          CALL wrd_mpi_io( 'lsm_global_end_index_h_' // dum, global_end_index )
    64426446
    64436447          IF ( .NOT. surface_data_to_write )  CYCLE
    64446448
    64456449          CALL wrd_mpi_io_surface( 't_soil_h(' // dum // ')', t_soil_h(l)%var_2d )
    6446           CALL wrd_mpi_io_surface( 'm_soil_h(' // dum // ')',  m_soil_h(l)%var_2d )
     6450          CALL wrd_mpi_io_surface( 'm_soil_h(' // dum // ')', m_soil_h(l)%var_2d )
    64476451          CALL wrd_mpi_io_surface( 'm_liq_h(' // dum // ')', m_liq_h(l)%var_1d )
    64486452          CALL wrd_mpi_io_surface( 't_surface_h(' // dum // ')', t_surface_h(l)%var_1d )
     
    64546458
    64556459          CALL rd_mpi_io_surface_filetypes ( surf_lsm_v(l)%start_index, surf_lsm_v(l)%end_index,   &
    6456                                              surface_data_to_write, global_start_index )
    6457 
    6458           CALL wrd_mpi_io( 'lsm_start_index_v_' // dum,  surf_lsm_v(l)%start_index )
    6459           CALL wrd_mpi_io( 'lsm_end_index_v_' // dum,  surf_lsm_v(l)%end_index )
    6460           CALL wrd_mpi_io( 'lsm_global_start_index_v_' // dum , global_start_index )
     6460                                            surface_data_to_write, global_start_index,             &
     6461                                            global_end_index )
     6462
     6463          CALL wrd_mpi_io( 'lsm_global_start_index_v_' // dum, global_start_index )
     6464          CALL wrd_mpi_io( 'lsm_global_end_index_v_' // dum, global_end_index )
    64616465
    64626466          IF ( .NOT. surface_data_to_write )  CYCLE
     
    71337137    INTEGER(iwp) ::  l   !< running index surface orientation
    71347138
    7135     !INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  end_index
    7136     INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_start
    7137     !INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  start_index
    7138 
    7139     LOGICAL      :: array_found
    7140     LOGICAL      :: ldum
     7139    INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_end_index
     7140    INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index
     7141
     7142    LOGICAL ::  array_found
     7143    LOGICAL ::  data_to_read    !< switch to steer reading of data
    71417144
    71427145
     
    72057208       WRITE( dum, '(I1)')  l
    72067209
    7207        CALL rrd_mpi_io( 'lsm_start_index_h_' // dum,  surf_lsm_h(l)%start_index )
    7208        CALL rrd_mpi_io( 'lsm_end_index_h_' // dum,  surf_lsm_h(l)%end_index )
    7209        CALL rrd_mpi_io( 'lsm_global_start_index_h_' // dum, global_start )
    7210 
    7211        CALL rd_mpi_io_surface_filetypes ( surf_lsm_h(l)%start_index, surf_lsm_h(l)%end_index, ldum,&
    7212                                           global_start )
    7213 
    7214        IF ( MAXVAL( surf_lsm_h(l)%end_index ) <= 0 )  CYCLE
     7210!
     7211!--    surf_lsm_h(l)%start_index and surf_lsm_h(l)%end_index are already set and should not be read
     7212!--    from restart file.
     7213       CALL rrd_mpi_io( 'lsm_global_start_index_h_' // dum, global_start_index )
     7214       CALL rrd_mpi_io( 'lsm_global_end_index_h_' // dum, global_end_index )
     7215
     7216       CALL rd_mpi_io_surface_filetypes ( surf_lsm_h(l)%start_index, surf_lsm_h(l)%end_index,      &
     7217                                          data_to_read, global_start_index, global_end_index )
     7218       IF ( .NOT. data_to_read )  CYCLE
    72157219
    72167220       CALL rrd_mpi_io_surface( 't_soil_h(' // dum // ')', t_soil_h(l)%var_2d )
     
    72187222       CALL rrd_mpi_io_surface( 'm_liq_h(' // dum // ')', m_liq_h(l)%var_1d )
    72197223       CALL rrd_mpi_io_surface( 't_surface_h(' // dum // ')', t_surface_h(l)%var_1d )
     7224
    72207225    ENDDO
    72217226
     
    72247229       WRITE( dum, '(I1)')  l
    72257230
    7226 !kk    In case of nothing to do, the settings of start_index and end_index differ
    7227 !kk    between writing and reading restart file
    7228 !kk
    7229 !kk    Has to be discussed with the developers
    7230 
    7231        CALL rrd_mpi_io( 'lsm_start_index_v_' // dum,  surf_lsm_v(l)%start_index )
    7232        CALL rrd_mpi_io( 'lsm_end_index_v_' // dum,  surf_lsm_v(l)%end_index )
    7233        CALL rrd_mpi_io( 'lsm_global_start_index_v_' // dum , global_start )
    7234 
    7235        CALL rd_mpi_io_surface_filetypes( surf_lsm_v(l)%start_index, surf_lsm_v(l)%end_index, ldum, &
    7236                                          global_start )
    7237 
    7238        IF ( MAXVAL( surf_lsm_v(l)%end_index ) <= 0 )  CYCLE
     7231!
     7232!--    surf_lsm_v(l)%start_index and surf_lsm_v(l)%end_index are already set and should not be read
     7233!--    from restart file.
     7234       CALL rrd_mpi_io( 'lsm_global_start_index_v_' // dum , global_start_index )
     7235       CALL rrd_mpi_io( 'lsm_global_end_index_v_' // dum , global_end_index )
     7236
     7237       CALL rd_mpi_io_surface_filetypes( surf_lsm_v(l)%start_index, surf_lsm_v(l)%end_index,       &
     7238                                         data_to_read, global_start_index, global_end_index )
     7239       IF ( .NOT. data_to_read )  CYCLE
    72397240
    72407241       CALL rrd_mpi_io_surface( 't_soil_v(' // dum // ')', t_soil_v(l)%var_2d )
  • palm/trunk/SOURCE/read_restart_data_mod.f90

    r4848 r4893  
    2525! -----------------
    2626! $Id$
     27! revised output of surface data via MPI-IO for better performance,
     28! therefore binary version number has changed
     29!
     30! 4848 2021-01-21 15:51:51Z gronemeier
    2731! bugfix: removed syn_turb_gen from restart files
    2832!
     
    276280    CALL location_message( 'read global restart data', 'start' )
    277281
     282!
     283!-- Caution: When any of the read instructions have been changed, the
     284!-- -------  version number stored in the variable binary_version_global has
     285!--          to be increased. The same changes must also be done in wrd_write_global.
     286    binary_version_global = '5.3'
     287
    278288    IF ( TRIM( restart_data_format_input ) == 'fortran_binary' )  THEN
    279289!
     
    286296       READ ( 13 )  version_on_file
    287297
    288        binary_version_global = '5.2'
    289298       IF ( TRIM( version_on_file ) /= TRIM( binary_version_global ) )  THEN
    290299          WRITE( message_string, * ) 'version mismatch concerning ',           &
     
    373382!
    374383!--    Now read all control parameters:
    375 !--    Caution: When the following read instructions have been changed, the
    376 !--    -------  version number stored in the variable binary_version_global has
    377 !--             to be increased. The same changes must also be done in
    378 !--             wrd_write_global.
     384
    379385       READ ( 13 )  length
    380386       READ ( 13 )  restart_string(1:length)
     
    873879!--    Read global restart data using MPI-IO
    874880!--    ATTENTION: Arrays need to be read with routine rrd_mpi_io_global_array!
    875 !--    Caution: When any of the following read instructions have been changed, the
    876 !--    -------  version number stored in the variable binary_version_global has
    877 !--             to be increased. The same changes must also be done in
    878 !--             wrd_write_global.
     881
    879882!
    880883!--    Open the MPI-IO restart file.
     
    886889       CALL rrd_mpi_io( 'binary_version_global',  version_on_file )
    887890
    888        binary_version_global = '5.1'
    889891       IF ( TRIM( version_on_file ) /= TRIM( binary_version_global ) )  THEN
    890892          WRITE( message_string, * ) 'version mismatch concerning binary_version_global:',         &
  • palm/trunk/SOURCE/restart_data_mpi_io_mod.f90

    r4857 r4893  
    2525! -----------------
    2626! $Id$
     27! revised output of surface data via MPI-IO for better performance
     28!
     29! 4857 2021-01-26 07:24:41Z raasch
    2730! bugfix: allocation of 3d-int4 array moved from particle output to standard output
    2831!
     
    155158               myidx,                                                                              &
    156159               myidy,                                                                              &
    157                npex,                                                                               &
    158                npey,                                                                               &
    159160               numprocs,                                                                           &
    160161               pdims
     
    183184    INTEGER(iwp)            ::  fh = -1          !< MPI-IO file handle
    184185#if defined( __parallel )
    185     INTEGER(iwp)            ::  fhs = -1         !< MPI-IO file handle to open file with comm2d always
    186 #endif
    187186    INTEGER(iwp)            ::  ft_surf = -1     !< MPI filetype surface data
    188 #if defined( __parallel )
    189187    INTEGER(iwp)            ::  ft_2di_nb        !< MPI filetype 2D array INTEGER no outer boundary
    190188    INTEGER(iwp)            ::  ft_2d            !< MPI filetype 2D array REAL with outer boundaries
     
    196194    INTEGER(iwp)            ::  glo_start        !< global start index on this PE
    197195#if defined( __parallel )
    198     INTEGER(iwp)            ::  local_start      !<
    199 #endif
    200     INTEGER(iwp)            ::  nr_iope          !<
    201     INTEGER(iwp)            ::  nr_val           !< local number of values in x and y direction
    202 #if defined( __parallel )
    203196    INTEGER(iwp)            ::  win_2di          !<
    204197    INTEGER(iwp)            ::  win_2dr          !<
     
    207200    INTEGER(iwp)            ::  win_3dr          !<
    208201    INTEGER(iwp)            ::  win_3ds          !<
     202    INTEGER(iwp)            ::  win_end   = -1   !<
     203    INTEGER(iwp)            ::  win_glost = -1   !<
     204    INTEGER(iwp)            ::  win_out   = -1   !<
     205    INTEGER(iwp)            ::  win_start = -1   !<
    209206    INTEGER(iwp)            ::  win_surf = -1    !<
    210207#endif
     
    216213    INTEGER(iwp), DIMENSION(:,:), POINTER, CONTIGUOUS   ::  array_2di   !<
    217214
    218     INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_end_index     !<
    219     INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_global_start  !<
     215    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  e_end_index     !< extended end index, every grid cell has at least one value
     216    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  e_start_index   !<
     217    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_end_index     !< module copy of end_index
    220218    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  m_start_index   !<
     219    INTEGER(iwp), DIMENSION(:),   ALLOCATABLE ::  thread_index    !<
     220    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  target_thread   !<
     221    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  transfer_index  !<
     222    INTEGER(iwp), DIMENSION(:),   ALLOCATABLE ::  thread_values   !<
     223!
     224!-- Indices for cyclic fill
     225    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  o_start_index   !<
     226    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  c_start_index   !<
     227!#if defined( __parallel )
     228    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  o_end_index     !< extended end index, every grid cell has at least one value
     229    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  c_end_index     !< extended end index, every grid cell has at least one value
     230    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  c_global_start  !<
     231    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  c_global_end    !<
     232!#endif
    221233
    222234    INTEGER(isp), DIMENSION(:,:,:), POINTER, CONTIGUOUS ::  array_3di4  !<
    223235    INTEGER(idp), DIMENSION(:,:,:), POINTER, CONTIGUOUS ::  array_3di8  !<
    224236
    225     LOGICAL ::  all_pes_write                 !< all PEs have data to write
    226237    LOGICAL ::  filetypes_created             !<
    227238    LOGICAL ::  io_on_limited_cores_per_node  !< switch to shared memory MPI-IO
     
    229240    LOGICAL ::  wr_flag                       !< file is opened for write
    230241
     242#if defined( __parallel )
     243    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE       :: local_indices
     244#endif
     245
     246    REAL(wp), DIMENSION(:), POINTER, CONTIGUOUS     ::  array_out      !<
    231247#if defined( __parallel )
    232248    REAL(wp), DIMENSION(:), POINTER, CONTIGUOUS     ::  array_1d       !<
     
    251267       INTEGER(iwp) :: nr_int         !< number of INTEGER entries in header
    252268       INTEGER(iwp) :: nr_real        !< number of REAL entries in header
     269       INTEGER(iwp) :: pes_along_x    !< number of PEs along x-direction during writing restart file
     270       INTEGER(iwp) :: pes_along_y    !< number of PEs along y-direction during writing restart file
    253271       INTEGER(iwp) :: total_nx       !< total number of points in x-direction
    254272       INTEGER(iwp) :: total_ny       !< total number of points in y-direction
    255273    END TYPE general_header
    256274
    257     TYPE(general_header), TARGET ::  tgh    !<
     275    TYPE(general_header), TARGET, PUBLIC ::  tgh    !<
    258276
    259277    TYPE(sm_class)               ::  sm_io  !<
     
    421439            wrd_mpi_io_surface
    422440
    423 
    424441 CONTAINS
    425442
     
    452469    TYPE(C_PTR)                   ::  buf_ptr  !<
    453470#endif
     471
    454472
    455473    offset = 0
     
    467485    io_file_name = file_name
    468486!
    469 !-- Setup for IO on a limited number of threads per node (using shared memory MPI)
     487!-- Setup for IO on a limited number of PEs per node (using shared memory MPI)
    470488    IF ( rd_flag )  THEN
    471489       set_filetype = .TRUE.
     
    822840
    823841!
    824 !--    TODO: describe in more detail what is done here and why it is done
    825 !--    save grid of main run
     842!--    Save grid information of the mainrun, i.e. grid variables like nxl, nxr, nys, nyn and other
     843!--    values are stored within the mainrun_grid structure
    826844       CALL mainrun_grid%save_grid_into_this_class()
    827845
     
    834852       rma_offset_s = 0
    835853!
    836 !--    Determine, if gridpoints of the prerun are located on this thread.
     854!--    Determine, if gridpoints of the prerun are located on this PE.
    837855!--    Set the (cyclic) prerun grid.
    838856       nxr = MIN( nxr, nx_on_file )
     
    857875       ny = ny_on_file
    858876!
    859 !--    Determine, if this thread is doing IO
     877!--    Determine, if this PE is doing IO
    860878       IF ( nnx > 0  .AND.  nny > 0 )  THEN
    861879          color = 1
     
    892910#endif
    893911!
    894 !--    Allocate 2d buffers as RMA window, accessible on all threads
     912!--    Allocate 2d buffers as RMA window, accessible on all PEs
    895913       IF ( pe_active_for_read )  THEN
    896914          ALLOCATE( rmabuf_2di(nys:nyn,nxl:nxr) )
     
    918936
    919937!
    920 !--    Allocate 3d buffer as RMA window, accessable on all threads
     938!--    Allocate 3d buffer as RMA window, accessable on all PEs
    921939       IF ( pe_active_for_read )  THEN
    922940          ALLOCATE( rmabuf_3d(nzb:nzt+1,nys:nyn,nxl:nxr) )
     
    932950
    933951!
    934 !--    TODO: comment in more detail, what is done here, and why
    935 !--    save small grid
     952!--    Save grid of the prerun, i.e. grid variables like nxl, nxr, nys, nyn and other values
     953!--    are stored within the prerun_grid structure.
     954!--    The prerun grid can later be activated by calling prerun_grid%activate_grid_from_this_class()
    936955       CALL prerun_grid%save_grid_into_this_class()
    937956       prerun_grid%comm2d = comm_cyclic_fill
     
    11521171
    11531172
    1154 !kk       write(9,*) 'Here is rma_cylic_fill_real_2d ',nxl,nxr,nys,nyn; FLUSH(9)
    1155 
    11561173!
    11571174!--    Reading 2d real array on prerun grid
     
    12971314!--                  array would be dimensioned in the caller subroutine like this:
    12981315!--                  INTEGER, DIMENSION(nysg:nyng,nxlg:nxrg)::  data
    1299           message_string = '2d-INTEGER array "' // TRIM( name ) // '" to be read from restart ' // &
    1300                            'file is defined with illegal dimensions in the PALM code'
     1316          message_string = '2d-INTEGER array with nbgp "' // TRIM( name ) // '" to be read ' //    &
     1317                           'from restart file is defined with illegal dimensions in the PALM code'
    13011318          CALL message( 'rrd_mpi_io_int_2d', 'PA0723', 3, 2, 0, 6, 0 )
    13021319
     
    13741391
    13751392
    1376        CALL prerun_grid%activate_grid_from_this_class()
    1377 
    13781393       IF ( pe_active_for_read )  THEN
     1394          CALL prerun_grid%activate_grid_from_this_class()
     1395
    13791396#if defined( __parallel )
    13801397          CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_INTEGER, ft_2di_nb, 'native',            &
    13811398                                  MPI_INFO_NULL, ierr )
    13821399          CALL MPI_FILE_READ_ALL( fh, array_2di, SIZE( array_2di ), MPI_INTEGER, status, ierr )
     1400#else
     1401          CALL posix_lseek( fh, array_position )
     1402          CALL posix_read( fh, array_2di, SIZE( array_2di ) )
    13831403#endif
    13841404          DO  i = nxl, nxr
     
    13861406          ENDDO
    13871407          data(1:nny,1:nnx) = rmabuf_2di
    1388        ENDIF
    1389 
    1390        CALL mainrun_grid%activate_grid_from_this_class()
     1408
     1409          CALL mainrun_grid%activate_grid_from_this_class()
     1410       ENDIF
    13911411
    13921412#if defined( __parallel )
     
    13961416#endif
    13971417
    1398        IF ( .NOT. pe_active_for_read )  THEN
    1399 
    1400           is = nxl
    1401           ie = nxr
    1402           js = nys
    1403           je = nyn
    1404 
    1405        ELSE
    1406 
    1407           is = nxl
    1408           ie = nxr
    1409           js = prerun_grid%nys+1
    1410           je = nyn
    1411           DO  i = is, ie
    1412              DO  j = js, je
    1413                 i_remote = MOD(i,nx_on_file+1)
    1414                 j_remote = MOD(j,ny_on_file+1)
    1415                 rem_pe   = remote_pe(i_remote,j_remote)
    1416                 rem_offs = rma_offset(i_remote,j_remote)
    1417                 nval     = 1
    1418 
    1419 #if defined( __parallel )
    1420                 IF ( rem_pe /= myid )  THEN
    1421                    CALL MPI_GET( data(j-nys+1,i-nxl+1), nval, MPI_INTEGER, rem_pe, rem_offs, nval, &
    1422                                  MPI_INTEGER, rmawin_2di, ierr )
    1423                 ELSE
    1424                    data(j-nys+1,i-nxl+1) = rmabuf_2di(j_remote,i_remote)
    1425                 ENDIF
    1426 #else
    1427                 data(j-nys+1,i-nxl+1) = array_2di(i_remote,j_remote)
    1428 #endif
    1429              ENDDO
    1430           ENDDO
    1431           is = prerun_grid%nxr+1
    1432           ie = nxr
    1433           js = nys
    1434           je = nyn
    1435 
    1436        ENDIF
     1418       is = nxl
     1419       ie = nxr
     1420       js = nys
     1421       je = nyn
    14371422
    14381423       DO  i = is, ie
     
    17171702                                  ierr )
    17181703          CALL MPI_FILE_READ_ALL( fh, array_3d, SIZE( array_3d ), MPI_REAL, status, ierr )
     1704#else
     1705          CALL posix_lseek( fh, array_position )
     1706          CALL posix_read( fh, array_3d, SIZE( array_3d ) )
    17191707#endif
    17201708          DO  i = nxl, nxr
     
    17271715#if defined( __parallel )
    17281716!
    1729 !--     Close RMA window to allow remote access
    1730         CALL MPI_WIN_FENCE( 0, rmawin_3d, ierr )
    1731 #endif
    1732 
    1733        IF ( .NOT. pe_active_for_read )  THEN
    1734 
    1735           is = nxl
    1736           ie = nxr
    1737           js = nys
    1738           je = nyn
    1739 
    1740        ELSE
    1741 
    1742           is = nxl
    1743           ie = nxr
    1744           js = prerun_grid%nys+1
    1745           je = nyn
    1746 
    1747           DO  i = is, ie
    1748              DO  j = js, je
    1749                 i_remote = MOD(i,nx_on_file+1)
    1750                 j_remote = MOD(j,ny_on_file+1)
    1751                 rem_pe   = remote_pe(i_remote,j_remote)
    1752                 rem_offs = rma_offset(i_remote,j_remote)*(nzt-nzb+2)
    1753                 nval     = nzt-nzb+2
    1754 
    1755 #if defined( __parallel )
    1756                 IF(rem_pe /= myid)   THEN
    1757                    CALL MPI_GET( data(nzb,j,i), nval, MPI_REAL, rem_pe, rem_offs, nval, MPI_REAL,  &
    1758                                  rmawin_3d, ierr)
    1759                 ELSE
    1760                    data(:,j,i) = rmabuf_3d(:,j_remote,i_remote)
    1761                 ENDIF
    1762 #else
    1763                 data(:,j,i) = array_3d(:,i_remote,j_remote)
    1764 #endif
    1765              ENDDO
    1766           ENDDO
    1767           is = prerun_grid%nxr+1
    1768           ie = nxr
    1769           js = nys
    1770           je = nyn
    1771 
    1772        ENDIF
     1717!--    Close RMA window to allow remote access
     1718       CALL MPI_WIN_FENCE( 0, rmawin_3d, ierr )
     1719#endif
     1720
     1721       is = nxl
     1722       ie = nxr
     1723       js = nys
     1724       je = nyn
    17731725
    17741726       DO  i = is, ie
    17751727          DO  j = js, je
    1776              i_remote = MOD(i,nx_on_file+1)
    1777              j_remote = MOD(j,ny_on_file+1)
     1728             i_remote = MOD( i, nx_on_file+1 )
     1729             j_remote = MOD( j, ny_on_file+1 )
    17781730             rem_pe   = remote_pe(i_remote,j_remote)
    17791731             rem_offs = rma_offset(i_remote,j_remote) * ( nzt-nzb+2 )
     
    18501802
    18511803    IF ( found )  THEN
    1852 #if defined( __parallel )
    18531804       CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
     1805#if defined( __parallel )
    18541806       CALL sm_io%sm_node_barrier()  ! Has no effect if I/O on limited number of cores is inactive
    18551807       IF ( sm_io%iam_io_pe )  THEN
     
    18741826       ENDIF
    18751827
     1828#if defined( __parallel )
     1829      IF ( sm_io%is_sm_active() )  THEN
     1830         CALL MPI_WIN_FREE( win_3ds, ierr )
     1831      ELSE
     1832         DEALLOCATE( array_3d_soil )
     1833      ENDIF
     1834#else
     1835      DEALLOCATE( array_3d_soil )
     1836#endif
     1837
    18761838    ELSE
    18771839
     
    20422004
    20432005    IF ( header_array_index == max_nr_arrays )  THEN
    2044        STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     2006       message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     2007       CALL message( 'wrd_mpi_io_real_2d', 'PA0585', 1, 2, 0, 6, 0 )
    20452008    ENDIF
    20462009
     
    21072070
    21082071    IF ( header_array_index == max_nr_arrays )  THEN
    2109        STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     2072       message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     2073       CALL message( 'wrd_mpi_io_int_2d', 'PA0585', 1, 2, 0, 6, 0 )
    21102074    ENDIF
    21112075
     
    21822146
    21832147    IF ( header_array_index == max_nr_arrays )  THEN
    2184        STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     2148       message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     2149       CALL message( 'wrd_mpi_io_int4_3d', 'PA0585', 1, 2, 0, 6, 0 )
    21852150    ENDIF
    21862151
     
    22252190                                      INT( (iog%nx+1), KIND = rd_offset_kind ) * isp
    22262191
    2227     write(9,*) 'array_position int4_3d ',trim(name),' ',array_position
    2228 
    22292192 END SUBROUTINE wrd_mpi_io_int4_3d
    22302193
     
    22502213
    22512214    IF ( header_array_index == max_nr_arrays )  THEN
    2252        STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     2215       message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     2216       CALL message( 'wrd_mpi_io_int8_3d', 'PA0585', 1, 2, 0, 6, 0 )
    22532217    ENDIF
    22542218
     
    22932257                                      INT( (iog%nx+1), KIND = rd_offset_kind ) * dp
    22942258
    2295     write(9,*) 'array_position int8_3d ',trim(name),' ',array_position
    2296 
    22972259 END SUBROUTINE wrd_mpi_io_int8_3d
    22982260
     
    23182280
    23192281    IF ( header_array_index == max_nr_arrays )  THEN
    2320        STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     2282       message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     2283       CALL message( 'wrd_mpi_io_real_3d', 'PA0585', 1, 2, 0, 6, 0 )
    23212284    ENDIF
    23222285
     
    23952358
    23962359    IF ( header_array_index == max_nr_arrays )  THEN
    2397        STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     2360       message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     2361       CALL message( 'wrd_mpi_io_real_3d_soil', 'PA0585', 1, 2, 0, 6, 0 )
    23982362    ENDIF
    23992363
     
    24022366    header_array_index = header_array_index + 1
    24032367
    2404 #if defined( __parallel )
    24052368    CALL rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
    2406 #endif
    24072369
    24082370    IF ( include_total_domain_boundaries)  THEN
     
    24322394    ENDIF
    24332395    CALL sm_io%sm_node_barrier()
     2396
     2397    IF ( sm_io%is_sm_active() )  THEN
     2398       CALL MPI_WIN_FREE( win_3ds, ierr )
     2399    ELSE
     2400       DEALLOCATE( array_3d_soil )
     2401    ENDIF
     2402    IF ( sm_io%iam_io_pe )  THEN
     2403       CALL MPI_TYPE_FREE( ft_3dsoil, ierr )
     2404    ENDIF
    24342405#else
    24352406    CALL posix_lseek( fh, array_position )
    24362407    CALL posix_write( fh, array_3d_soil, SIZE( array_3d_soil ) )
     2408    DEALLOCATE( array_3d_soil )
    24372409#endif
    24382410!
     
    25892561         CALL MPI_BCAST( data, SIZE( data ), MPI_REAL, 0, comm2d, ierr )
    25902562       ELSE
    2591           IF ( sm_io%iam_io_pe )  THEN
     2563          IF( sm_io%iam_io_pe )  THEN
    25922564             CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
     2565          ENDIF
     2566          IF ( myid == 0 )  THEN
    25932567             CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
    2594              CALL MPI_FILE_READ_ALL( fh, data, SIZE( data ), MPI_REAL, status, ierr )
     2568             CALL MPI_FILE_READ( fh, data, SIZE( data ), MPI_REAL, status, ierr )
    25952569          ENDIF
    2596           IF ( sm_io%is_sm_active() )  THEN
    2597              CALL MPI_BCAST( data, SIZE( data ), MPI_REAL, 0, sm_io%comm_shared, ierr )
    2598           ENDIF
     2570          CALL MPI_BCAST( data, SIZE( data ), MPI_REAL, 0, comm2d, ierr )
    25992571       ENDIF
    26002572#else
     
    27492721             CALL MPI_FILE_READ_ALL( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
    27502722          ENDIF
    2751           CALL MPI_BCAST( data, SIZE( data ), MPI_REAL, 0, comm2d, ierr )
     2723          CALL MPI_BCAST( data, SIZE( data ), MPI_INTEGER, 0, comm2d, ierr )
    27522724       ELSE
    2753           IF ( sm_io%iam_io_pe )  THEN
     2725          IF( sm_io%iam_io_pe )  THEN
    27542726             CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
     2727          ENDIF
     2728          IF ( myid == 0 )  THEN
    27552729             CALL MPI_FILE_SEEK( fh, array_position, MPI_SEEK_SET, ierr )
    2756              CALL MPI_FILE_READ_ALL( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
     2730             CALL MPI_FILE_READ( fh, data, SIZE( data), MPI_INTEGER, status, ierr )
    27572731          ENDIF
    2758           IF ( sm_io%is_sm_active() )  THEN
    2759              CALL MPI_BCAST( data, SIZE( data ), MPI_INTEGER, 0, sm_io%comm_shared, ierr )
    2760           ENDIF
    2761        ENDIF
     2732          CALL MPI_BCAST( data, SIZE( data ), MPI_INTEGER, 0, comm2d, ierr )
     2733        ENDIF
    27622734#else
    27632735       CALL posix_lseek( fh, array_position )
     
    28002772
    28012773    IF ( header_array_index == max_nr_arrays )  THEN
    2802        STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     2774       message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     2775       CALL message( 'wrd_mpi_io_global_array_real_1d', 'PA0585', 1, 2, 0, 6, 0 )
    28032776    ENDIF
    28042777
     
    29392912
    29402913    IF ( header_array_index == max_nr_arrays )  THEN
    2941        STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     2914       message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     2915       CALL message( 'wrd_mpi_io_global_array_int_1d', 'PA0585', 1, 2, 0, 6, 0 )
    29422916    ENDIF
    29432917
     
    30303004       ENDDO
    30313005
    3032        write(9,*) 'particle_size_read ',particle_size,array_size,array_position,sum(prt_global_index)
    3033 
    30343006       ALLOCATE( prt_data(MAX(array_size,1)) )
    30353007
     
    30783050       array_position = prt_nr_bytes
    30793051
    3080        write(9,*) 'array_position after particle read ',array_position,prt_nr_bytes,rs
    3081 
    30823052       DEALLOCATE( prt_data )
    30833053
     
    30923062! ------------
    30933063!> Read 1d-REAL surface data array with MPI-IO.
    3094 !--------------------------------------------------------------------------------------------------!
    3095  SUBROUTINE rrd_mpi_io_surface( name, data, first_index )
     3064!> This is a recursive subroutine. In case of cyclic fill mode it may call itself for reading parts
     3065!> of the prerun grid.
     3066!--------------------------------------------------------------------------------------------------!
     3067 RECURSIVE SUBROUTINE rrd_mpi_io_surface( name, data, first_index )
    30963068
    30973069    IMPLICIT NONE
     
    30993071    CHARACTER(LEN=*), INTENT(IN) ::  name            !<
    31003072
     3073    INTEGER(iwp), OPTIONAL       ::  first_index     !<
     3074    INTEGER(iwp)                 ::  i               !<
     3075    INTEGER(iwp)                 ::  j               !<
     3076    INTEGER(iwp)                 ::  lo_first_index  !<
     3077
     3078#if defined( __parallel )
     3079    INTEGER(iwp)                 ::  buf_start       !<
    31013080    INTEGER(KIND=rd_offset_kind) ::  disp            !< displacement of actual indices
    3102     INTEGER(KIND=rd_offset_kind) ::  disp_f          !< displacement in file
    3103     INTEGER(KIND=rd_offset_kind) ::  disp_n          !< displacement of next column
    3104     INTEGER(iwp), OPTIONAL       ::  first_index     !<
    3105 
    3106     INTEGER(iwp)                 ::  i               !<
    3107     INTEGER(iwp)                 ::  i_f             !<
    3108     INTEGER(iwp)                 ::  j               !<
    3109     INTEGER(iwp)                 ::  j_f             !<
    3110     INTEGER(iwp)                 ::  lo_first_index  !<
    3111     INTEGER(iwp)                 ::  nr_bytes        !<
    3112     INTEGER(iwp)                 ::  nr_bytes_f      !<
    3113     INTEGER(iwp)                 ::  nr_words        !<
    3114 #if defined( __parallel )
    3115     INTEGER, DIMENSION(rd_status_size)  ::  status   !<
    3116 #else
    3117     TYPE(C_PTR)                         ::  buf      !<
    3118 #endif
    3119 
    3120     LOGICAL                             ::  found    !<
     3081    INTEGER(iwp)                 ::  ie              !<
     3082    INTEGER(iwp)                 ::  ind_gb          !<
     3083    INTEGER(iwp)                 ::  ind_out         !<
     3084    INTEGER(iwp)                 ::  is              !<
     3085    INTEGER(iwp)                 ::  n               !<
     3086    INTEGER(iwp)                 ::  n_trans         !<
     3087
     3088    INTEGER(iwp),DIMENSION(0:numprocs-1) ::  lo_index  !<
     3089    INTEGER, DIMENSION(rd_status_size)   ::  status    !<
     3090#endif
     3091    LOGICAL                      ::  found  !<
    31213092
    31223093    REAL(wp), INTENT(OUT), DIMENSION(:), TARGET ::  data  !<
     3094#if defined( __parallel )
     3095    REAL(wp),DIMENSION(:),ALLOCATABLE    ::  put_buffer  !<
     3096#endif
    31233097
    31243098
     
    31323106    DO  i = 1, tgh%nr_arrays
    31333107        IF ( TRIM( array_names(i) ) == TRIM( name ) )  THEN
     3108!
     3109!--        ATTENTION: The total_number_of_surface_values and wp MUST be INTERGER(8).
     3110!--        The compiler (at least Intel) first computes total_number_of_surface_values*wp
     3111!--        and then does the conversion to INTEGER(8).
     3112!--        This may lead to wrong results when total_number_of_surface_values*wp is > 2*10**6
    31343113           array_position = array_offset(i) + ( lo_first_index - 1 ) *                             &
    3135                             total_number_of_surface_values * wp
     3114                            INT( total_number_of_surface_values, idp ) * INT( wp, idp )
    31363115           found = .TRUE.
    31373116           EXIT
     
    31393118    ENDDO
    31403119
    3141     disp   = -1
    3142     disp_f = -1
    3143     disp_n = -1
     3120!
     3121!-- In case of 2d-data, name is written only once
     3122    IF ( lo_first_index == 1 )  THEN
     3123
     3124       IF ( header_array_index == max_nr_arrays )  THEN
     3125          message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     3126          CALL message( 'rrd_mpi_io_surface', 'PA0585', 1, 2, 0, 6, 0 )
     3127       ENDIF
     3128
     3129       array_names(header_array_index)  = name
     3130       array_offset(header_array_index) = array_position
     3131       header_array_index = header_array_index + 1
     3132
     3133    ENDIF
     3134
    31443135    IF ( found )  THEN
    3145 
    31463136       IF ( cyclic_fill_mode )  THEN
    31473137
    31483138          CALL rrd_mpi_io_surface_cyclic_fill
     3139          RETURN
    31493140
    31503141       ELSE
    3151 
    3152           IF ( MAXVAL( m_global_start ) == -1 )   RETURN   ! Nothing to do on this PE
     3142#if defined( __parallel )
     3143!
     3144!--       Read data from restart file
     3145          CALL sm_io%sm_node_barrier() ! has no effect if I/O on limited number of cores is inactive
     3146          IF ( sm_io%iam_io_pe )  THEN
     3147             CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_surf, 'native',              &
     3148                                     MPI_INFO_NULL, ierr )
     3149             CALL MPI_FILE_READ_ALL ( fh, array_out, SIZE(array_out), MPI_REAL, status, ierr )
     3150          ENDIF
     3151          CALL sm_io%sm_node_barrier()
     3152
     3153!
     3154!--       Copy data into transfer buffer. Data is organized in a way that only one MPI_PUT to the
     3155!--       respective PE ist required.
     3156          ALLOCATE( put_buffer(SUM( transfer_index(4,:) )) )
     3157
     3158          ind_gb = 1
     3159          DO  i = 1, SIZE( local_indices, 2 )
     3160             ind_out = local_indices(1,i)
     3161             DO  j = 1, local_indices(2,i)
     3162                put_buffer(ind_gb) = array_out(ind_out)
     3163                ind_out = ind_out + 1
     3164                ind_gb  = ind_gb  + 1
     3165             ENDDO
     3166          ENDDO
     3167!
     3168!--       Transfer data from I/O PEs to the respective PEs to which they belong.
     3169          CALL MPI_WIN_FENCE( 0, win_surf, ierr )
     3170
     3171          buf_start = 1
     3172          DO  n = 0, numprocs-1
     3173             n_trans = transfer_index(4,n)
     3174             IF ( n_trans > 0 )  THEN
     3175                disp = transfer_index(3,n) - 1
     3176                CALL MPI_PUT( put_buffer(buf_start), n_trans, MPI_REAL, n, disp, n_trans, MPI_REAL,&
     3177                              win_surf, ierr)
     3178                buf_start = buf_start + n_trans
     3179             ENDIF
     3180          ENDDO
     3181
     3182          CALL MPI_WIN_FENCE( 0, win_surf, ierr )
     3183          DEALLOCATE( put_buffer )
     3184!
     3185!--       Copy from RMA window into output array (data) to allow transfering data to target PEs.
     3186!--       Check, if the number of surface values per grid cell match the index setup.
     3187          lo_index = thread_values
    31533188          DO  i = nxl, nxr
    31543189             DO  j = nys, nyn
    3155 
    3156                 IF ( m_global_start(j,i) > 0 )  THEN
    3157                    disp     = array_position+(m_global_start(j,i)-1) * wp
    3158                    nr_words = m_end_index(j,i)-m_start_index(j,i)+1
    3159                    nr_bytes = nr_words * wp
     3190                is = lo_index(target_thread(j,i)) + 1
     3191                ie = is + m_end_index(j,i) - m_start_index(j,i)
     3192                data(m_start_index(j,i):m_end_index(j,i)) = array_1d(is:ie)
     3193                lo_index(target_thread(j,i)) = lo_index(target_thread(j,i)) +                      &
     3194                                               e_end_index(j,i) - e_start_index(j,i) + 1
     3195!
     3196!--             TODO: Test can be removed later.
     3197                IF ( e_end_index(j,i)-e_start_index(j,i)+1 /= NINT( array_1d(is-1) ) )  THEN
     3198                   WRITE( 9, '(A,6I8)' ) 'Nr surface values does not match ', j, i,                &
     3199                                         e_start_index(j,i), e_end_index(j,i),                     &
     3200                                         e_end_index(j,i)-e_start_index(j,i)+1 ,                   &
     3201                                         NINT( array_1d(is-1) )
     3202                   FLUSH( 9 )
     3203                   CALL MPI_ABORT( comm2d, 1, ierr )
    31603204                ENDIF
    3161                 IF ( disp >= 0  .AND.  disp_f == -1 )  THEN   ! First entry
    3162                    disp_f     = disp
    3163                    nr_bytes_f = 0
    3164                    i_f = i
    3165                    j_f = j
    3166                 ENDIF
    3167                 IF ( j == nyn  .AND.  i == nxr )  THEN        ! Last entry
    3168                    disp_n = -1
    3169                    IF (  nr_bytes > 0 )  THEN
    3170                       nr_bytes_f = nr_bytes_f+nr_bytes
    3171                    ENDIF
    3172                 ELSEIF ( j == nyn )  THEN                     ! Next x
    3173                    IF ( m_global_start(nys,i+1) > 0  .AND.  disp > 0 )  THEN
    3174                       disp_n = array_position + ( m_global_start(nys,i+1) - 1 ) * wp
    3175                    ELSE
    3176                       CYCLE
    3177                    ENDIF
    3178                 ELSE
    3179                    IF ( m_global_start(j+1,i) > 0  .AND.  disp > 0 )  THEN
    3180                       disp_n = array_position + ( m_global_start(j+1,i) - 1 ) * wp
    3181                    ELSE
    3182                       CYCLE
    3183                    ENDIF
    3184                 ENDIF
    3185 
    3186 
    3187                 IF ( disp + nr_bytes == disp_n )  THEN        ! Contiguous block
    3188                    nr_bytes_f = nr_bytes_f + nr_bytes
    3189                 ELSE                                          ! Read
    3190 #if defined( __parallel )
    3191                    CALL MPI_FILE_SEEK( fhs, disp_f, MPI_SEEK_SET, ierr )
    3192                    nr_words = nr_bytes_f / wp
    3193                    CALL MPI_FILE_READ( fhs, data(m_start_index(j_f,i_f)), nr_words, MPI_REAL, status, &
    3194                       ierr )
    3195 #else
    3196 !
    3197 !--                Use C_PTR here, because posix read does not work with indexed array
    3198                    buf = C_LOC( data(m_start_index(j_f,i_f)) )
    3199                    CALL posix_lseek( fh, disp_f )
    3200                    CALL posix_read( fh, buf, nr_bytes_f )
    3201 #endif
    3202                    disp_f     = disp
    3203                    nr_bytes_f = nr_bytes
    3204                    i_f = i
    3205                    j_f = j
    3206                 ENDIF
    3207 
    32083205             ENDDO
    32093206          ENDDO
    3210        ENDIF
    3211 
    3212 
    3213     ELSE
    3214 
    3215        message_string = 'surface array "' // TRIM( name ) // '" not found in restart file'
    3216        CALL message( 'rrd_mpi_io_surface', 'PA0722', 3, 2, 0, 6, 0 )
     3207
     3208
     3209#else
     3210          CALL posix_lseek( fh, array_position )
     3211          CALL posix_read( fh, array_out, SIZE(array_out) )
     3212
     3213          DO  i = nxl, nxr
     3214             DO  j = nys, nyn
     3215                data(m_start_index(j,i):m_end_index(j,i)) =                                        &
     3216                                                   array_out(e_start_index(j,i)+1:e_end_index(j,i))
     3217!
     3218!--             TODO: Test can be removed later.
     3219                IF ( e_end_index(j,i)-e_start_index(j,i)+1 /= NINT(array_out(e_start_index(j,i))) )&
     3220                THEN
     3221                   WRITE( 9, '(A,6I8)' ) 'Nr surface values does not match ', j, i,                &
     3222                                         e_start_index(j,i), e_end_index(j,i),                     &
     3223                                         e_end_index(j,i)-e_start_index(j,i)+1,                    &
     3224                                         NINT( array_out(e_start_index(j,i)) )
     3225                   FLUSH( 9 )
     3226                   CALL ABORT()
     3227                ENDIF
     3228             ENDDO
     3229          ENDDO
     3230#endif
     3231       ENDIF
    32173232
    32183233    ENDIF
     
    32263241       INTEGER(iwp) ::  i         !<
    32273242       INTEGER(iwp) ::  ie        !<
    3228 #if defined( __parallel )
    3229        INTEGER(iwp) ::  ierr      !<
    3230 #endif
    32313243       INTEGER(iwp) ::  is        !<
    32323244       INTEGER(iwp) ::  i_remote  !<
     
    32413253       INTEGER(KIND=MPI_ADDRESS_KIND) ::  rem_offs  !<
    32423254#else
    3243        INTEGER(idp) ::  rem_offs
    3244 #endif
    3245 
    3246        LOGICAL ::  write_done  !<
    3247 
    3248 
    3249 !
    3250 !--    In the current version, there is only 1 value per grid cell allowed.
    3251 !--    In this special case, the cyclical repetition can be done with the same method as for 2d-real
    3252 !--    array.
     3255       INTEGER(idp) ::  rem_offs                    !<
     3256#endif
     3257
     3258       REAL(wp), DIMENSION(:), ALLOCATABLE ::  c_data  !<
     3259
     3260
     3261!
     3262!--    ATTENTION: This version allows only 1 surface element per grid cell.
     3263!
     3264!--    Activate grid of the smaller prerun, i.e. grid variables like nxl, nxr, nys, nyn and other
     3265!--    values are set according to the prerun settings.
    32533266       CALL prerun_grid%activate_grid_from_this_class()
    32543267
    32553268       IF ( pe_active_for_read )  THEN
    3256           rmabuf_2d = -1.0
     3269
     3270          IF ( MAXVAL( m_end_index ) <= 0 )  THEN
     3271             CALL mainrun_grid%activate_grid_from_this_class()
     3272             IF ( debug_output )  THEN
     3273                CALL debug_message( 'PE inactive for reading restart or prerun data', 'start' )
     3274             ENDIF
     3275             RETURN
     3276          ENDIF
     3277
     3278          ALLOCATE( c_data(MAXVAL( m_end_index )) )
     3279
     3280!
     3281!--       Recursive CALL of rrd_mpi_io_surface.
     3282!--       rrd_mpi_io_surface is called with cyclic_fill_mode = .FALSE. on the smaller prerun grid.
     3283          cyclic_fill_mode = .FALSE.
     3284          CALL rrd_mpi_io_surface( name, c_data )
     3285          cyclic_fill_mode = .TRUE.
     3286
    32573287          DO  i = nxl, nxr
    32583288             DO  j = nys, nyn
    3259 
    3260                 IF ( m_global_start(j,i) > 0 )  THEN
    3261                    disp     = array_position+(m_global_start(j,i)-1) * wp
    3262                    nr_words = m_end_index(j,i)-m_start_index(j,i)+1
    3263                    nr_bytes = nr_words * wp
    3264                 ENDIF
    3265                 IF ( disp >= 0  .AND.  disp_f == -1 )  THEN   ! First entry
    3266                    disp_f     = disp
    3267                    nr_bytes_f = 0
    3268                    write_done = .TRUE.
    3269                 ENDIF
    3270                 IF( write_done )  THEN
    3271                    i_f = i
    3272                    j_f = j
    3273                    write_done = .FALSE.
    3274                 ENDIF
    3275 
    3276                 IF ( j == nyn  .AND.  i == nxr )  THEN        ! Last entry
    3277                    disp_n = -1
    3278                    IF (  nr_bytes > 0 )  THEN
    3279                       nr_bytes_f = nr_bytes_f+nr_bytes
    3280                    ENDIF
    3281                 ELSEIF ( j == nyn )  THEN                     ! Next x
    3282                    IF ( m_global_start(nys,i+1) > 0  .AND.  disp > 0 )  THEN
    3283                       disp_n = array_position + ( m_global_start(nys,i+1) - 1 ) * wp
    3284                    ELSE
    3285                       CYCLE
    3286                    ENDIF
    3287                 ELSE
    3288                    IF ( m_global_start(j+1,i) > 0  .AND.  disp > 0 )  THEN
    3289                       disp_n = array_position + ( m_global_start(j+1,i) - 1 ) * wp
    3290                    ELSE
    3291                       CYCLE
    3292                    ENDIF
    3293                 ENDIF
    3294 
    3295 
    3296                 IF ( disp + nr_bytes == disp_n )  THEN        ! Contiguous block
    3297                    nr_bytes_f = nr_bytes_f + nr_bytes
    3298                 ELSE                                          ! Read
    3299 #if defined( __parallel )
    3300                    CALL MPI_FILE_SEEK( fhs, disp_f, MPI_SEEK_SET, ierr )
    3301                    nr_words = nr_bytes_f / wp
    3302                    CALL MPI_FILE_READ( fhs, rmabuf_2d(j_f,i_f), nr_words, MPI_REAL, status, ierr )
    3303 #else
    3304                    CALL posix_lseek( fh, disp_f )
    3305                    CALL posix_read( fh, rmabuf_2d(j_f:,i_f:), nr_bytes_f )
    3306 #endif
    3307 
    3308                    disp_f     = disp
    3309                    nr_bytes_f = nr_bytes
    3310                    write_done = .TRUE.
    3311                 ENDIF
    3312 
     3289                rmabuf_2d(j,i) = c_data(c_start_index(j,i))
    33133290             ENDDO
    33143291          ENDDO
    33153292
    33163293       ENDIF
    3317 
     3294!
     3295!--    Activate grid of the mainrun, i.e. grid variables like nxl, nxr, nys, nyn and other values
     3296!--    are set according to the mainrun settings.
    33183297       CALL mainrun_grid%activate_grid_from_this_class()
    33193298
     
    33243303#endif
    33253304
    3326        IF ( .NOT. pe_active_for_read )  THEN
    3327 
    3328           is = nxl
    3329           ie = nxr
    3330           js = nys
    3331           je = nyn
    3332 
    3333        ELSE
    3334 
    3335           is = nxl
    3336           ie = nxr
    3337           js = prerun_grid%nys+1
    3338           je = nyn
    3339 
    3340           DO  i = is, ie
    3341              DO  j = js, je
    3342                 i_remote = MOD(i,nx_on_file+1)
    3343                 j_remote = MOD(j,ny_on_file+1)
    3344                 rem_pe   = remote_pe(i_remote,j_remote)
    3345                 rem_offs = rma_offset(i_remote,j_remote)
    3346                 nval     = 1
    3347 
    3348 #if defined( __parallel )
    3349                 IF ( rem_pe /= myid )  THEN
    3350                    CALL MPI_GET( data(m_start_index(j,i)), nval, MPI_REAL, rem_pe, rem_offs, nval, &
    3351                                  MPI_REAL, rmawin_2d, ierr)
    3352                 ELSE
    3353                    data(m_start_index(j,i)) = rmabuf_2d(j_remote,i_remote)
    3354                 ENDIF
    3355 #else
    3356                 data(m_start_index(j,i)) = array_2d(i_remote,j_remote)
    3357 #endif
    3358              ENDDO
    3359           ENDDO
    3360           is = prerun_grid%nxr+1
    3361           ie = nxr
    3362           js = nys
    3363           je = nyn
    3364 
    3365        ENDIF
     3305!
     3306!--   After reading surface data on the small grid, map these data in a cyclic way to all respective
     3307!--   grid points of the main run.
     3308      is = nxl
     3309      ie = nxr
     3310      js = nys
     3311      je = nyn
    33663312
    33673313       DO  i = is, ie
    33683314          DO  j = js, je
    3369              i_remote = MOD(i,nx_on_file+1)
    3370              j_remote = MOD(j,ny_on_file+1)
     3315             i_remote = MOD( i, nx_on_file+1 )
     3316             j_remote = MOD( j, ny_on_file+1 )
    33713317             rem_pe   = remote_pe(i_remote,j_remote)
    33723318             rem_offs = rma_offset(i_remote,j_remote)
     
    33753321#if defined( __parallel )
    33763322             IF ( rem_pe /= myid )  THEN
    3377                 CALL MPI_GET( data(m_start_index(j,i)), nval, MPI_REAL, rem_pe, rem_offs, nval,    &
     3323                CALL MPI_GET( data(o_start_index(j,i)), nval, MPI_REAL, rem_pe, rem_offs, nval,    &
    33783324                              MPI_REAL, rmawin_2d, ierr)
    33793325             ELSE
    3380                 data(m_start_index(j,i)) = rmabuf_2d(j_remote,i_remote)
     3326                data(o_start_index(j,i)) = rmabuf_2d(j_remote,i_remote)
    33813327             ENDIF
    33823328#else
    3383              data(m_start_index(j,i)) = array_2d(i_remote,j_remote)
     3329             data(o_start_index(j,i)) = array_2d(i_remote,j_remote)
    33843330#endif
    33853331          ENDDO
     
    33913337       CALL MPI_WIN_FENCE( 0, rmawin_2d, ierr )
    33923338#endif
     3339
     3340       IF ( ALLOCATED( c_data ) )  DEALLOCATE( c_data )
    33933341
    33943342    END SUBROUTINE rrd_mpi_io_surface_cyclic_fill
     
    35393487    array_position = prt_nr_bytes
    35403488
    3541     write(9,*) 'array_position after particle ',array_position,prt_nr_bytes,rs
    3542 
    35433489    DEALLOCATE( prt_data )
    35443490
     
    35563502    IMPLICIT NONE
    35573503
    3558     CHARACTER(LEN=*), INTENT(IN)       ::  name            !<
    3559 
    3560 #if defined( __parallel )
    3561     INTEGER(KIND=rd_offset_kind)       ::  disp            !<
    3562 #endif
    3563     INTEGER(iwp), OPTIONAL             ::  first_index     !<
    3564 #if defined( __parallel )
    3565     INTEGER(iwp)                       ::  i               !<
    3566 #endif
    3567     INTEGER(iwp)                       ::  lo_first_index  !<
    3568     INTEGER(KIND=rd_offset_kind)       ::  offset          !<
    3569 
    3570 #if defined( __parallel )
    3571     INTEGER, DIMENSION(rd_status_size) ::  status          !<
    3572 #endif
    3573 
    3574     REAL(wp), INTENT(IN), DIMENSION(:), TARGET ::  data    !<
    3575 
    3576 
    3577     offset = 0
     3504    CHARACTER(LEN=*), INTENT(IN) ::  name  !<
     3505
     3506    INTEGER(iwp), OPTIONAL ::  first_index     !<
     3507    INTEGER(iwp) ::  i               !<
     3508    INTEGER(iwp) ::  j               !<
     3509    INTEGER(iwp) ::  lo_first_index  !<
     3510#if defined( __parallel )
     3511    INTEGER(iwp) ::  buf_start       !<
     3512    INTEGER(iwp) ::  ie              !<
     3513    INTEGER(iwp) ::  is              !<
     3514    INTEGER(iwp) ::  ind_gb          !<
     3515    INTEGER(iwp) ::  ind_out         !<
     3516    INTEGER(iwp) ::  n               !<
     3517    INTEGER(iwp) ::  n_trans         !<
     3518#endif
     3519
     3520#if defined( __parallel )
     3521    INTEGER(KIND=MPI_ADDRESS_KIND) ::  disp    !< displacement in RMA window
     3522    INTEGER(KIND=rd_offset_kind)   ::  offset  !<
     3523
     3524    INTEGER(iwp), DIMENSION(0:numprocs-1)   ::  lo_index  !<
     3525    INTEGER(iwp), DIMENSION(rd_status_size) ::  status    !<
     3526#endif
     3527
     3528    REAL(wp), INTENT(IN), DIMENSION(:), TARGET ::  data  !<
     3529#if defined( __parallel )
     3530    REAL(wp), DIMENSION(:), ALLOCATABLE ::  get_buffer  !<
     3531#endif
     3532
     3533
    35783534    lo_first_index = 1
    35793535
     
    35813537       lo_first_index = first_index
    35823538    ENDIF
     3539
    35833540!
    35843541!-- In case of 2d-data, name is written only once
     
    35863543
    35873544       IF ( header_array_index == max_nr_arrays )  THEN
    3588           STOP '+++ maximum number of 2d/3d-array entries in restart file header exceeded'
     3545          message_string = 'maximum number of 2d/3d-array entries in restart file header exceeded'
     3546          CALL message( 'wrd_mpi_io_surface', 'PA0585', 1, 2, 0, 6, 0 )
    35893547       ENDIF
    35903548
     
    35963554
    35973555#if defined( __parallel )
    3598     IF ( sm_io%is_sm_active() )  THEN
    3599        DO  i = 1, nr_val
    3600           array_1d(i+local_start) = data(i)
     3556    offset = 0
     3557
     3558    ALLOCATE( get_buffer(SUM( transfer_index(4,:) )) )
     3559!
     3560!-- Copy from input array (data) to RMA window to allow the target PEs to get the appropiate data.
     3561!-- At this point, a dummy surface element is added. This makes sure that every x-y grid cell owns
     3562!-- at least one surface element. This way, bookkeeping becomes much easier.
     3563    lo_index = thread_values
     3564    DO  i = nxl, nxr
     3565       DO  j = nys, nyn
     3566          is = lo_index(target_thread(j,i)) + 1
     3567          ie = is + m_end_index(j,i) - m_start_index(j,i)
     3568!
     3569!--       Store number of surface elements in dummy additional surface element
     3570          array_1d(is-1)  = e_end_index(j,i) - e_start_index(j,i) + 1
     3571          array_1d(is:ie) = data(m_start_index(j,i):m_end_index(j,i))
     3572          lo_index(target_thread(j,i)) = lo_index(target_thread(j,i)) +                            &
     3573                                         e_end_index(j,i) - e_start_index(j,i) + 1
    36013574       ENDDO
    3602     ELSE
    3603 !       array_1d => data                           !kk Did not work in all cases    why???
    3604        ALLOCATE( array_1d( SIZE( data ) ) )
    3605        array_1d = data
    3606     ENDIF
    3607 
    3608     CALL sm_io%sm_node_barrier()  ! Has no effect if I/O on limited number of cores is inactive
     3575    ENDDO
     3576!
     3577!-- On target PE, get data from source PEs which are assigned for output on this PE.
     3578    CALL MPI_WIN_FENCE( 0, win_surf, ierr )
     3579
     3580    buf_start = 1
     3581    DO  n = 0, numprocs-1
     3582       n_trans = transfer_index(4,n)
     3583       IF ( n_trans > 0 )  THEN
     3584          disp = transfer_index(3,n) - 1
     3585          CALL MPI_GET( get_buffer(buf_start), n_trans, MPI_REAL, n, disp, n_trans, MPI_REAL,      &
     3586                        win_surf, ierr )
     3587          buf_start = buf_start + n_trans
     3588       ENDIF
     3589    ENDDO
     3590
     3591    CALL MPI_WIN_FENCE( 0, win_surf, ierr )
     3592!
     3593!-- Copy data to output buffer. Here, the outpuf buffer matches the indices global_start and
     3594!-- global_end.
     3595    ind_gb  = 1
     3596    DO  i = 1, SIZE( local_indices, 2 )
     3597       ind_out = local_indices(1,i)
     3598       DO  j = 1, local_indices(2,i)
     3599          array_out(ind_out) = get_buffer(ind_gb)
     3600          ind_out = ind_out+1
     3601          ind_gb  = ind_gb+1
     3602       ENDDO
     3603    ENDDO
     3604
     3605    DEALLOCATE( get_buffer )
     3606
     3607!
     3608!-- Write data to disk.
     3609    CALL sm_io%sm_node_barrier()  ! has no effect if I/O on limited number of cores is inactive
    36093610    IF ( sm_io%iam_io_pe )  THEN
    3610        IF ( all_pes_write )  THEN
    3611           CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_surf, 'native', MPI_INFO_NULL,  &
    3612                                   ierr )
    3613           CALL MPI_FILE_WRITE_ALL( fh, array_1d, nr_iope, MPI_REAL, status, ierr )
    3614        ELSE
    3615           CALL MPI_FILE_SET_VIEW( fh, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
    3616           IF ( nr_val > 0 )  THEN
    3617              disp = array_position + 8 * ( glo_start - 1 )
    3618              CALL MPI_FILE_SEEK( fh, disp, MPI_SEEK_SET, ierr )
    3619              CALL MPI_FILE_WRITE( fh, array_1d, nr_iope, MPI_REAL, status, ierr )
    3620           ENDIF
    3621        ENDIF
     3611       CALL MPI_FILE_SET_VIEW( fh, array_position, MPI_REAL, ft_surf, 'native', MPI_INFO_NULL,     &
     3612                               ierr )
     3613       CALL MPI_FILE_WRITE_ALL( fh, array_out, SIZE( array_out ), MPI_REAL, status, ierr )
    36223614    ENDIF
    36233615    CALL sm_io%sm_node_barrier()
    3624     IF( .NOT. sm_io%is_sm_active() )  DEALLOCATE( array_1d )
    3625 #else
     3616#else
     3617    DO  i = nxl, nxr
     3618       DO  j = nys, nyn
     3619          array_out(e_start_index(j,i)) = e_end_index(j,i) - e_start_index(j,i) + 1
     3620          array_out(e_start_index(j,i)+1:e_end_index(j,i)) =                                       &
     3621                                                          data(m_start_index(j,i):m_end_index(j,i))
     3622       ENDDO
     3623    ENDDO
     3624
    36263625    CALL posix_lseek( fh, array_position )
    3627     CALL posix_write( fh, data, nr_val )
    3628 #endif
     3626    CALL posix_write( fh, array_out, SIZE(array_out) )
     3627#endif
     3628
    36293629    array_position = array_position + total_number_of_surface_values * wp
    36303630
    3631 !    IF ( lo_first_index == 1 )  THEN
    3632 !       IF ( debug_level >= 2 .AND. nr_val  > 0 )  WRITE(9,*) 'w_surf_1 ', TRIM( name ), ' ', nr_val, SUM( data(1:nr_val) )
    3633 !    ELSE
    3634 !       IF ( debug_level >= 2 .AND. nr_val  > 0 ) WRITE(9,*) 'w_surf_n ', TRIM( name ), ' ', &
    3635 !                                                            lo_first_index, nr_val, SUM( data(1:nr_val) )
    3636 !    ENDIF
    3637 
    36383631 END SUBROUTINE wrd_mpi_io_surface
    3639 
    36403632
    36413633
     
    36903682    IF ( wr_flag  .AND.  sm_io%iam_io_pe )  THEN
    36913683
    3692        tgh%nr_int    = header_int_index - 1
    3693        tgh%nr_char   = header_char_index - 1
    3694        tgh%nr_real   = header_real_index - 1
    3695        tgh%nr_arrays = header_array_index - 1
    3696        tgh%total_nx  = iog%nx + 1
    3697        tgh%total_ny  = iog%ny + 1
     3684       tgh%nr_int      = header_int_index - 1
     3685       tgh%nr_char     = header_char_index - 1
     3686       tgh%nr_real     = header_real_index - 1
     3687       tgh%nr_arrays   = header_array_index - 1
     3688       tgh%total_nx    = iog%nx + 1
     3689       tgh%total_ny    = iog%ny + 1
     3690       tgh%pes_along_x = pdims(1)
     3691       tgh%pes_along_y = pdims(2)
    36983692       IF ( include_total_domain_boundaries )  THEN   ! Not sure, if LOGICAL interpretation is the same for all compilers,
    36993693          tgh%i_outer_bound = 1                       ! therefore store as INTEGER in general header
     
    38023796!-- Close MPI-IO files
    38033797#if defined( __parallel )
    3804 !
    3805 !-- Restart file has been opened with comm2d
    3806     IF ( fhs /= -1 )  THEN
    3807        CALL MPI_FILE_CLOSE( fhs, ierr )
    3808     ENDIF
    38093798!
    38103799!-- Free RMA windows
     
    38163805#endif
    38173806
    3818     IF (.NOT. pe_active_for_read )  RETURN
     3807    IF ( ALLOCATED( e_start_index ) )   DEALLOCATE( e_start_index  )
     3808    IF ( ALLOCATED( e_end_index )   )   DEALLOCATE( e_end_index    )
     3809    IF ( ALLOCATED( m_start_index ) )   DEALLOCATE( m_start_index  )
     3810    IF ( ALLOCATED( m_end_index )   )   DEALLOCATE( m_end_index    )
     3811    IF ( ALLOCATED( target_thread ) )   DEALLOCATE( target_thread  )
     3812    IF ( ALLOCATED( thread_index )  )   DEALLOCATE( thread_index   )
     3813    IF ( ALLOCATED( thread_values ) )   DEALLOCATE( thread_values  )
     3814    IF ( ALLOCATED( transfer_index ) )  DEALLOCATE( transfer_index  )
     3815
     3816    IF ( .NOT. pe_active_for_read )  RETURN
    38193817!
    38203818!-- TODO: better explain the following message
     
    38603858!> data is not time critical (data size is comparably small), it will be read by all cores.
    38613859!--------------------------------------------------------------------------------------------------!
    3862  SUBROUTINE rd_mpi_io_surface_filetypes( start_index, end_index, data_to_write, global_start )
     3860 RECURSIVE SUBROUTINE rd_mpi_io_surface_filetypes( start_index, end_index, data_to_write,          &
     3861                                                   global_start, global_end )
    38633862
    38643863    IMPLICIT NONE
    38653864
    3866     INTEGER(iwp)                          ::  i           !<  loop index
    3867     INTEGER(iwp)                          ::  j           !<  loop index
    3868     INTEGER(KIND=rd_offset_kind)          ::  offset      !<
    3869 
    3870     INTEGER(iwp), DIMENSION(1)            ::  dims1       !<
    3871     INTEGER(iwp), DIMENSION(1)            ::  lize1       !<
    3872     INTEGER(iwp), DIMENSION(1)            ::  start1      !<
    3873 
    3874     INTEGER(iwp), DIMENSION(0:numprocs-1) ::  all_nr_val  !< number of values for all PEs
    3875     INTEGER(iwp), DIMENSION(0:numprocs-1) ::  lo_nr_val   !< local number of values in x and y direction
    3876 
    3877 
    3878     INTEGER, INTENT(INOUT), DIMENSION(nys:nyn,nxl:nxr)  ::  end_index     !<
    3879     INTEGER, INTENT(OUT), DIMENSION(nys:nyn,nxl:nxr)    ::  global_start  !<
    3880     INTEGER, INTENT(INOUT), DIMENSION(nys:nyn,nxl:nxr)  ::  start_index   !<
    3881 
    3882     LOGICAL, INTENT(OUT) ::  data_to_write  !< returns, if surface data have to be written
    3883 
    3884 !
    3885 !-- Actions during reading
    3886     IF ( rd_flag )  THEN
    3887 !
    3888 !--    Set start index and end index for the mainrun grid.
    3889 !--    ATTENTION: This works only for horizontal surfaces with one vale per grid cell!!!
    3890        IF ( cyclic_fill_mode )  THEN
    3891           DO  i = nxl, nxr
    3892              DO  j = nys, nyn
    3893                 start_index (j,i) = (i-nxl) * nny + j - nys + 1
    3894                 end_index (j,i)   = start_index(j,i)
    3895              ENDDO
     3865    INTEGER(iwp) ::  e_lo_start        !<
     3866    INTEGER(iwp) ::  i                 !<  loop index
     3867    INTEGER(iwp) ::  j                 !<  loop index
     3868    INTEGER(iwp) ::  index_offset      !<
     3869    INTEGER(iwp) ::  last_end_index    !<
     3870    INTEGER(iwp) ::  lo_start          !<
     3871    INTEGER(iwp) ::  nr_surfcells_pe   !<
     3872    INTEGER(iwp) ::  rest_cells_pe     !<
     3873    INTEGER(iwp) ::  rest_bound        !<
     3874#if defined( __parallel )
     3875    INTEGER(iwp) ::  io_end_index      !<
     3876    INTEGER(iwp) ::  io_start_index    !<
     3877    INTEGER(iwp) ::  n                 !<  loop index
     3878    INTEGER(iwp) ::  nr_previous       !<
     3879#endif
     3880
     3881    INTEGER(iwp), DIMENSION(0:numprocs-1,2) ::  nr_surfcells_all_s  !<
     3882    INTEGER(iwp), DIMENSION(0:numprocs-1,2) ::  nr_surfcells_all_r  !<
     3883#if defined( __parallel )
     3884    INTEGER(iwp), DIMENSION(1)              ::  dims1               !< global dimension for MPI_TYPE_CREATE_SUBARRAY
     3885    INTEGER(iwp), DIMENSION(1)              ::  lsize1              !< local size for MPI_TYPE_CREATE_SUBARRAY
     3886    INTEGER(iwp), DIMENSION(0:numprocs-1)   ::  nr_cells_to_thread  !<
     3887    INTEGER(iwp), DIMENSION(0:pdims(1))     ::  nr_surf_cells_x     !<
     3888    INTEGER(iwp), DIMENSION(0:pdims(1))     ::  nr_surf_cells_x_s   !<
     3889    INTEGER(iwp), DIMENSION(0:numprocs-1)   ::  nr_values_to_thread !<
     3890    INTEGER(iwp), DIMENSION(1)              ::  start1              !< start index for MPI_TYPE_CREATE_SUBARRAY
     3891    INTEGER(iwp), DIMENSION(nxl:nxr)        ::  sum_y               !<
     3892#endif
     3893
     3894    INTEGER(iwp), INTENT(INOUT), DIMENSION(nys:nyn,nxl:nxr) ::  end_index     !< local end indx
     3895    INTEGER(iwp), INTENT(INOUT), DIMENSION(nys:nyn,nxl:nxr) ::  global_start  !< global start index
     3896    INTEGER(iwp), INTENT(INOUT), DIMENSION(nys:nyn,nxl:nxr) ::  global_end    !< global end index
     3897    INTEGER(iwp), INTENT(INOUT), DIMENSION(nys:nyn,nxl:nxr) ::  start_index   !< local start index
     3898#if defined( __parallel )
     3899    INTEGER(iwp), DIMENSION(0:myidy,nxl:nxr)    ::  nr_previous_y      !<
     3900    INTEGER(iwp), DIMENSION(0:pdims(2),nxl:nxr) ::  nr_surf_cells_y    !<
     3901    INTEGER(iwp), DIMENSION(0:pdims(2),nxl:nxr) ::  nr_surf_cells_y_s  !<
     3902    INTEGER(iwp), DIMENSION(4,0:numprocs-1)     ::  transfer_index_s   !<
     3903#endif
     3904
     3905    LOGICAL, INTENT(OUT) ::  data_to_write      !< returns .TRUE., if surface data have been written
     3906    LOGICAL              ::  only_dummy_values  !< only dummy values, i.e. no data to write
     3907
     3908
     3909!
     3910!-- Clean up previous calls.
     3911#if defined( __parallel )
     3912    IF ( win_surf /= -1 )  THEN
     3913       CALL MPI_WIN_FREE( win_surf, ierr )
     3914       DEALLOCATE( array_1d )
     3915       win_surf = -1
     3916    ENDIF
     3917    IF ( ft_surf /= -1  .AND.  sm_io%iam_io_pe )  THEN
     3918       CALL MPI_TYPE_FREE( ft_surf, ierr )
     3919    ENDIF
     3920    ft_surf = -1
     3921    IF ( sm_io%is_sm_active() )  THEN
     3922       IF ( win_out /= -1 )  THEN
     3923          CALL MPI_WIN_FREE( win_out, ierr )
     3924          win_out = -1
     3925       ENDIF
     3926    ELSE
     3927       IF ( ASSOCIATED( array_out ) )  DEALLOCATE( array_out )
     3928    ENDIF
     3929#else
     3930    IF ( ASSOCIATED( array_out ) )  DEALLOCATE( array_out )
     3931#endif
     3932
     3933    IF ( cyclic_fill_mode )  THEN
     3934       CALL cyclic_fill_surface_filetype
     3935       RETURN
     3936    ELSE
     3937       IF ( .NOT. ALLOCATED( e_end_index )    )  ALLOCATE( e_end_index(nys:nyn,nxl:nxr)   )
     3938       IF ( .NOT. ALLOCATED( e_start_index )  )  ALLOCATE( e_start_index(nys:nyn,nxl:nxr) )
     3939       IF ( .NOT. ALLOCATED( m_end_index )    )  ALLOCATE( m_end_index(nys:nyn,nxl:nxr)   )
     3940       IF ( .NOT. ALLOCATED( m_start_index )  )  ALLOCATE( m_start_index(nys:nyn,nxl:nxr) )
     3941       IF ( .NOT. ALLOCATED( target_thread )  )  ALLOCATE( target_thread(nys:nyn,nxl:nxr) )
     3942       IF ( .NOT. ALLOCATED( thread_index )   )  ALLOCATE( thread_index(0:numprocs-1)     )
     3943       IF ( .NOT. ALLOCATED( thread_values )  )  ALLOCATE( thread_values(0:numprocs-1)    )
     3944       IF ( .NOT. ALLOCATED( transfer_index ) )  ALLOCATE( transfer_index(4,0:numprocs-1) )
     3945    ENDIF
     3946
     3947    IF ( wr_flag)  THEN
     3948!
     3949!--    Add one dummy value at every grid box.
     3950!--    This allows to use MPI_FILE_WRITE_ALL and MPI_FILE_READ_ALL with subarray file type.
     3951       index_offset   = 0
     3952       last_end_index = 0
     3953       DO  i = nxl, nxr
     3954          DO  j = nys, nyn
     3955             e_start_index(j,i) = start_index (j,i) + index_offset
     3956             IF ( end_index (j,i) - start_index(j,i) < 0 )  THEN
     3957                e_end_index (j,i) = last_end_index+1
     3958                last_end_index    = last_end_index+1
     3959             ELSE
     3960                e_end_index (j,i) = end_index(j,i) + index_offset + 1
     3961                last_end_index    = e_end_index (j,i)
     3962             ENDIF
     3963             index_offset = index_offset + 1
     3964           ENDDO
     3965       ENDDO
     3966#if defined( __parallel )
     3967!
     3968!--    Compute indices for global, PE independent 1-d surface element array.
     3969       nr_surf_cells_y_s = 0
     3970!
     3971!--    Count number of surface elements in y-direction for every x.
     3972       DO  i = nxl, nxr
     3973          nr_surf_cells_y_s(myidy,i) = SUM( e_end_index (:,i) - e_start_index (:,i) + 1 )
     3974       ENDDO
     3975!
     3976!--    Distribute these values to all PEs along y.
     3977       CALL MPI_ALLREDUCE( nr_surf_cells_y_s, nr_surf_cells_y, SIZE( nr_surf_cells_y ),            &
     3978                           MPI_INTEGER, MPI_SUM, comm1dy, ierr )
     3979!
     3980!--    Sum all surface elements along y for individual x PEs
     3981       nr_surf_cells_x_s        = 0
     3982       nr_surf_cells_x_s(myidx) = SUM( nr_surf_cells_y )
     3983!
     3984!--    Distribute to all PEs along x.
     3985       CALL MPI_ALLREDUCE( nr_surf_cells_x_s, nr_surf_cells_x, SIZE( nr_surf_cells_x ),            &
     3986                           MPI_INTEGER, MPI_SUM, comm1dx, ierr )
     3987       DO  i = nxl, nxr
     3988          nr_previous_y(:,i) = 0
     3989          DO  n = 1, myidy
     3990             nr_previous_y(n,i) = nr_previous_y(n-1,i) + nr_surf_cells_y(n-1,i)
    38963991          ENDDO
    3897        ENDIF
    3898 
    3899        IF ( .NOT. ALLOCATED( m_start_index )  )  ALLOCATE( m_start_index(nys:nyn,nxl:nxr)  )
    3900        IF ( .NOT. ALLOCATED( m_end_index )    )  ALLOCATE( m_end_index(nys:nyn,nxl:nxr)    )
    3901        IF ( .NOT. ALLOCATED( m_global_start ) )  ALLOCATE( m_global_start(nys:nyn,nxl:nxr) )
    3902 !
    3903 !--    Save arrays for later reading
    3904        m_start_index  = start_index
    3905        m_end_index    = end_index
    3906        m_global_start = global_start
    3907        nr_val = MAXVAL( end_index )
    3908 
    3909     ENDIF
    3910 
    3911     IF ( .NOT. pe_active_for_read )  RETURN
    3912 
    3913     IF ( cyclic_fill_mode )  CALL prerun_grid%activate_grid_from_this_class()
    3914 
    3915     offset = 0
    3916     lo_nr_val= 0
    3917     lo_nr_val(myid) = MAXVAL( end_index )
    3918 #if defined( __parallel )
    3919     CALL MPI_ALLREDUCE( lo_nr_val, all_nr_val, numprocs, MPI_INTEGER, MPI_SUM, comm2d, ierr )
    3920     IF ( ft_surf /= -1  .AND.  sm_io%iam_io_pe )  THEN
    3921        CALL MPI_TYPE_FREE( ft_surf, ierr )    ! If set, free last surface filetype
    3922     ENDIF
    3923 
    3924     IF ( win_surf /= -1 )  THEN
    3925        IF ( sm_io%is_sm_active() )  THEN
    3926           CALL MPI_WIN_FREE( win_surf, ierr )
    3927        ENDIF
    3928        win_surf = -1
    3929     ENDIF
    3930 
    3931     IF ( sm_io%is_sm_active() .AND. rd_flag )  THEN
    3932        IF ( fhs == -1 )  THEN
    3933           CALL MPI_FILE_OPEN( comm2d, TRIM( io_file_name ), MPI_MODE_RDONLY, MPI_INFO_NULL, fhs,   &
    3934                               ierr )
    3935        ENDIF
     3992       ENDDO
     3993
     3994       sum_y(nxl) = SUM( nr_surf_cells_y(:,nxl) )
     3995       DO  i = nxl, nxr
     3996          IF ( i > nxl )  THEN
     3997             sum_y(i) = sum_y(i-1) + SUM( nr_surf_cells_y(:,i) )
     3998          ENDIF
     3999       ENDDO
     4000
     4001       nr_previous = 0
     4002       IF ( myidx >= 1 )  THEN
     4003          nr_previous = SUM(nr_surf_cells_x(0:myidx-1))
     4004       ENDIF
     4005
     4006       global_start(nys,nxl) = 1 + nr_previous + nr_previous_y(myidy,nxl)
     4007       DO  j = nys+1, nyn
     4008          global_start(j,nxl) = global_start(j-1,nxl) + e_end_index(j-1,nxl) -                     &
     4009                                e_start_index(j-1,nxl) + 1
     4010       ENDDO
     4011
     4012       DO  i = nxl+1, nxr
     4013          global_start(nys,i) = 1 + nr_previous + nr_previous_y(myidy,i) + sum_y(i-1)
     4014          DO  j = nys+1, nyn
     4015             global_start(j,i) = global_start(j-1,i) + e_end_index(j-1,i) - e_start_index(j-1,i) + 1
     4016          ENDDO
     4017       ENDDO
     4018#else
     4019       global_start = e_start_index
     4020#endif
     4021       DO  i = nxl, nxr
     4022          DO  j = nys, nyn
     4023             global_end(j,i) = global_start(j,i) + e_end_index (j,i) - e_start_index (j,i)
     4024          ENDDO
     4025       ENDDO
     4026
    39364027    ELSE
    3937        fhs = fh
    3938     ENDIF
    3939 #else
    3940     all_nr_val(myid) = lo_nr_val(myid)
    3941 #endif
    3942     nr_val = lo_nr_val(myid)
     4028!
     4029!--    In case of read, compute e_start_index and e_end_index for current processor grid.
     4030!--    This data contains one extra value for every i and j.
     4031       e_lo_start = 1
     4032       lo_start   = 1
     4033       DO  i = nxl, nxr
     4034          DO  j = nys, nyn
     4035             e_start_index(j,i) = e_lo_start
     4036             e_end_index(j,i)   = e_lo_start + global_end(j,i) - global_start(j,i)
     4037             e_lo_start         = e_lo_start + global_end(j,i) - global_start(j,i) + 1
     4038             start_index(j,i)   = lo_start
     4039             end_index(j,i)     = lo_start + global_end(j,i) - global_start(j,i) - 1
     4040             lo_start           = lo_start + global_end(j,i) - global_start(j,i)
     4041          ENDDO
     4042       ENDDO
     4043
     4044    ENDIF
     4045
     4046    nr_surfcells_all_s = 0
     4047    nr_surfcells_all_s(myid,1) = MAXVAL( e_end_index ) ! don't split surface elements of one gridbox
     4048    nr_surfcells_all_s(myid,2) = MAXVAL( e_end_index - e_start_index )
     4049
     4050#if defined( __parallel )
     4051    CALL MPI_ALLREDUCE( nr_surfcells_all_s, nr_surfcells_all_r, SIZE( nr_surfcells_all_s ),        &
     4052                        MPI_INTEGER, MPI_SUM, comm2d, ierr )
     4053#else
     4054    nr_surfcells_all_r = nr_surfcells_all_s
     4055#endif
    39434056
    39444057    total_number_of_surface_values = 0
     
    39474060          glo_start = total_number_of_surface_values + 1
    39484061       ENDIF
    3949        total_number_of_surface_values = total_number_of_surface_values + all_nr_val(i)
     4062       total_number_of_surface_values = total_number_of_surface_values + nr_surfcells_all_r(i,1)
    39504063    ENDDO
    3951 
    3952 !
    3953 !-- Actions during reading
    3954     IF ( rd_flag )  THEN
    3955 
    3956 #if defined( __parallel )
    3957        CALL MPI_FILE_SET_VIEW( fhs, offset, MPI_BYTE, MPI_BYTE, 'native', MPI_INFO_NULL, ierr )
    3958 #endif
    3959     ENDIF
    3960 
    3961     IF ( cyclic_fill_mode )  CALL mainrun_grid%activate_grid_from_this_class()
    3962 
    3963 !
    3964 !-- Actions during writing
    3965     IF ( wr_flag )  THEN
    3966 !
    3967 !--    Create surface filetype
    3968        ft_surf      = -1
    3969        global_start = start_index + glo_start - 1
    3970 
    3971        WHERE ( end_index < start_index )
    3972           global_start = -1
    3973        ENDWHERE
    3974 
    3975 #if defined( __parallel )
    3976        IF ( sm_io%is_sm_active() )  THEN
    3977           IF ( sm_io%iam_io_pe )  THEN
    3978 !
    3979 !--          Calculate number of values of all PEs of an I/O group
    3980              nr_iope = 0
    3981              DO  i = myid, myid+sm_io%sh_npes-1
    3982                 nr_iope = nr_iope + all_nr_val(i)
     4064    only_dummy_values = ( MAXVAL( nr_surfcells_all_r(:,2) ) <= 0 )
     4065
     4066!
     4067!-- Compute indices of equally distributed surface elements.
     4068!-- Number of surface values scheduled for ouput on this PE:
     4069    nr_surfcells_pe  = total_number_of_surface_values  / numprocs
     4070    rest_cells_pe    = MOD( total_number_of_surface_values, numprocs )
     4071    rest_bound       = rest_cells_pe * ( nr_surfcells_pe + 1 )
     4072    m_start_index    = start_index
     4073    m_end_index      = end_index
     4074
     4075!
     4076!-- Compute number of elements on source PE, which have to be send to the corresponding target PE.
     4077#if defined( __parallel )
     4078    nr_cells_to_thread  = 0
     4079    nr_values_to_thread = 0
     4080    DO  i = nxl, nxr
     4081       DO  j = nys, nyn
     4082          IF ( rest_cells_pe == 0 )  THEN
     4083             target_thread(j,i) = ( global_start(j,i) - 1 ) / nr_surfcells_pe
     4084          ELSE
     4085             IF ( global_start(j,i) <= rest_bound )  THEN
     4086                target_thread(j,i) = ( global_start(j,i) - 1 ) / ( nr_surfcells_pe + 1 )
     4087             ELSE
     4088                target_thread(j,i) = ( global_start(j,i) - rest_bound - 1 ) / nr_surfcells_pe
     4089                target_thread(j,i) = target_thread(j,i) + rest_cells_pe
     4090             ENDIF
     4091!
     4092!--          TODO: Test output, to be removed later.
     4093             IF ( target_thread(j,i) >= numprocs )  THEN
     4094                WRITE( 9,'(A,8I8)' )  'target_thread ', j, i, target_thread(j,i),                  &
     4095                                      global_start(j,i) , nr_surfcells_pe
     4096                FLUSH( 9 )
     4097                CALL MPI_ABORT( comm2d, 1, ierr )
     4098             ENDIF
     4099          ENDIF
     4100          nr_cells_to_thread(target_thread(j,i))  = nr_cells_to_thread(target_thread(j,i)) + 1
     4101          nr_values_to_thread(target_thread(j,i)) = nr_values_to_thread(target_thread(j,i)) +      &
     4102                                                    e_end_index(j,i) - e_start_index(j,i) + 1
     4103       ENDDO
     4104    ENDDO
     4105
     4106!
     4107!-- Compute start index in the transfer buffer on the source side for the corresponding target PE.
     4108    thread_index(0)  = 1
     4109    thread_values(0) = 1
     4110    DO  n = 1, numprocs-1
     4111       thread_index(n)  = thread_index(n-1) + nr_cells_to_thread(n-1)
     4112       thread_values(n) = thread_values(n-1) + nr_values_to_thread(n-1)
     4113    ENDDO
     4114!
     4115!-- Buffer distribution on the source side.
     4116    DO  n = 0, numprocs-1
     4117       transfer_index_s(1,n) = thread_index(n)
     4118       transfer_index_s(2,n) = nr_cells_to_thread(n)
     4119       transfer_index_s(3,n) = thread_values(n)
     4120       transfer_index_s(4,n) = nr_values_to_thread(n)
     4121    ENDDO
     4122
     4123    CALL MPI_ALLTOALL( transfer_index_s, 4, MPI_INTEGER, transfer_index, 4, MPI_INTEGER, comm2d,   &
     4124                       ierr)
     4125!
     4126!-- Buffer distribution on the target side side.
     4127    CALL get_remote_indices()
     4128!
     4129!-- Create surface element file type.
     4130    IF ( total_number_of_surface_values > 0 .AND. .NOT. only_dummy_values)  THEN
     4131        data_to_write = .TRUE.
     4132    ELSE
     4133        data_to_write = .FALSE.
     4134    ENDIF
     4135
     4136    CALL MPI_ALLREDUCE( global_end(nyn,nxr), dims1(1), 1, MPI_INTEGER, MPI_MAX, comm2d, ierr )
     4137    start1(1) = MINVAL( local_indices(1,:) ) - 1
     4138    IF ( sm_io%is_sm_active() )  THEN
     4139       CALL MPI_ALLREDUCE( SUM( local_indices(2,:) ), lsize1(1), 1, MPI_INTEGER, MPI_SUM,          &
     4140                           sm_io%comm_shared, ierr )
     4141    ELSE
     4142       lsize1(1) = SUM( local_indices(2,:) )
     4143    ENDIF
     4144
     4145    IF ( sm_io%iam_io_pe )  THEN
     4146       IF ( total_number_of_surface_values > 0 )  THEN
     4147           CALL MPI_TYPE_CREATE_SUBARRAY( 1, dims1, lsize1, start1, MPI_ORDER_FORTRAN, MPI_REAL,   &
     4148                                          ft_surf, ierr )
     4149           CALL MPI_TYPE_COMMIT( ft_surf, ierr )
     4150       ENDIF
     4151    ENDIF
     4152!
     4153!-- Allocate rma window to supply surface data to other PEs.
     4154    CALL rd_alloc_rma_mem( array_1d, SUM( nr_values_to_thread ), win_surf )
     4155!
     4156!-- Allocate shared array on IO-PE to supply data for MPI-IO (write or read).
     4157    IF ( sm_io%is_sm_active() )  THEN
     4158       IF ( sm_io%iam_io_pe )  THEN
     4159          io_start_index = start1(1) + 1
     4160          io_end_index   = start1(1) + lsize1(1)
     4161       ENDIF
     4162       CALL MPI_BCAST( io_start_index, 1, MPI_INTEGER, 0, sm_io%comm_shared, ierr )
     4163       CALL MPI_BCAST( io_end_index,   1, MPI_INTEGER, 0, sm_io%comm_shared, ierr )
     4164       CALL sm_io%sm_allocate_shared( array_out, io_start_index, io_end_index, win_out )
     4165    ELSE
     4166       ALLOCATE( array_out(start1(1)+1:start1(1)+lsize1(1)) )
     4167    ENDIF
     4168#else
     4169    IF ( total_number_of_surface_values > 0  .AND.  .NOT. only_dummy_values )  THEN
     4170        data_to_write = .TRUE.
     4171    ELSE
     4172        data_to_write = .FALSE.
     4173    ENDIF
     4174    ALLOCATE( array_out(1:total_number_of_surface_values) )
     4175#endif
     4176
     4177 CONTAINS
     4178
     4179    SUBROUTINE cyclic_fill_surface_filetype
     4180
     4181       INTEGER(iwp) ::  i  !<  loop index
     4182       INTEGER(iwp) ::  j  !<  loop index
     4183
     4184
     4185       IF ( .NOT. ALLOCATED( o_start_index ) )  ALLOCATE( o_start_index(nys:nyn,nxl:nxr) )
     4186       IF ( .NOT. ALLOCATED( o_end_index )   )  ALLOCATE( o_end_index(nys:nyn,nxl:nxr)   )
     4187
     4188       lo_start   = 1
     4189       DO  i = nxl, nxr
     4190           DO  j = nys, nyn
     4191               o_start_index(j,i) = lo_start
     4192               o_end_index(j,i)   = lo_start
     4193               lo_start           = lo_start + 1
     4194           ENDDO
     4195       ENDDO
     4196       start_index = o_start_index
     4197       end_index   = o_end_index
     4198
     4199       IF ( MAXVAL( global_end-global_start ) > 1 )  THEN
     4200          message_string = 'cylic-fill method does not allow more than one surface element ' //    &
     4201                           'per grid box'
     4202          CALL message( 'cyclic_fill_surface_filetype', 'PA0742', 1, 2, 0, 6, 0 )
     4203       ENDIF
     4204!
     4205!--    Activate grid of the smaller prerun, i.e. grid variables like nxl, nxr, nys, nyn and others
     4206!--    are set according to the prerun layout.
     4207       CALL prerun_grid%activate_grid_from_this_class()
     4208
     4209       IF ( pe_active_for_read )  THEN
     4210
     4211          IF ( .NOT. ALLOCATED( c_global_start ) )  ALLOCATE( c_global_start(nys:nyn,nxl:nxr) )
     4212          IF ( .NOT. ALLOCATED( c_global_end )   )  ALLOCATE( c_global_end(nys:nyn,nxl:nxr)   )
     4213          IF ( .NOT. ALLOCATED( c_start_index )  )  ALLOCATE( c_start_index(nys:nyn,nxl:nxr)  )
     4214          IF ( .NOT. ALLOCATED( c_end_index )    )  ALLOCATE( c_end_index(nys:nyn,nxl:nxr)    )
     4215
     4216          DO  i = nxl, nxr
     4217             DO  j = nys, nyn
     4218                c_global_start(j,i) = global_start(j,i)
     4219                c_global_end(j,i)   = global_end(j,i)
    39834220             ENDDO
    3984           ELSE
    3985              local_start = 0
    3986              DO  i = myid-sm_io%sh_rank, myid-1
    3987                 local_start = local_start + all_nr_val(i)
    3988              ENDDO
     4221          ENDDO
     4222!
     4223!--       Recursive call of rd_mpi_io_surface_filetypes.
     4224!--       Prerun data are read, but they are treated as if they are mainrun data, just on a smaller
     4225!--       grid.
     4226          cyclic_fill_mode = .FALSE.
     4227          CALL rd_mpi_io_surface_filetypes( c_start_index, c_end_index, data_to_write,             &
     4228                                            c_global_start, c_global_end )
     4229          cyclic_fill_mode = .TRUE.
     4230
     4231       ENDIF
     4232!
     4233!--    Activate grid of the mainrun, i.e. grid variables like nxl, nxr, nys, nyn and others
     4234!--    are set according to the mainrun layout.
     4235       CALL mainrun_grid%activate_grid_from_this_class()
     4236
     4237#if defined( __parallel )
     4238       CALL MPI_BCAST( data_to_write, 1, MPI_LOGICAL, 0, comm2d, ierr )
     4239#endif
     4240
     4241    END SUBROUTINE cyclic_fill_surface_filetype
     4242
     4243#if defined( __parallel )
     4244!
     4245!-- Get the indices of the surface elements inside the RMA window on the remote PE.
     4246!-- This information is required to fetch the surface element data on remote PEs
     4247!-- in rrd_mpi_io_surface and wrd_mpi_io_surface.
     4248    SUBROUTINE get_remote_indices
     4249
     4250       INTEGER(iwp) ::  buf_start  !<
     4251       INTEGER(iwp) ::  i          !<
     4252       INTEGER(iwp) ::  j          !<
     4253       INTEGER(iwp) ::  n          !<
     4254       INTEGER(iwp) ::  n_trans    !<
     4255       INTEGER(iwp) ::  win_ind    !<
     4256
     4257       INTEGER(KIND=MPI_ADDRESS_KIND) ::  disp     !< displacement in RMA window
     4258       INTEGER(KIND=MPI_ADDRESS_KIND) ::  winsize  !< size of RMA window
     4259
     4260       INTEGER(iwp), DIMENSION(0:numprocs-1) ::  lo_index         !<
     4261
     4262       INTEGER(iwp), POINTER, DIMENSION(:,:) ::  surf_val_index  !<
     4263
     4264
     4265       IF ( ALLOCATED( local_indices ) )  DEALLOCATE( local_indices )
     4266       ALLOCATE( local_indices(2,MAX( SUM( transfer_index(2,:) ), 2 )))
     4267
     4268       local_indices(1,:) = 0
     4269       local_indices(2,:) = 0
     4270
     4271       winsize = MAX( 2 * SUM( nr_cells_to_thread ), 2 )
     4272
     4273       ALLOCATE( surf_val_index(2,winsize) )
     4274       winsize = winsize * iwp
     4275       CALL MPI_WIN_CREATE( surf_val_index, winsize, iwp, MPI_INFO_NULL, comm2d, win_ind, ierr )
     4276       CALL MPI_WIN_FENCE( 0, win_ind, ierr )
     4277
     4278       lo_index = thread_index
     4279       DO  i = nxl, nxr
     4280          DO  j = nys, nyn
     4281             surf_val_index(1,lo_index(target_thread(j,i))) = global_start(j,i)
     4282             surf_val_index(2,lo_index(target_thread(j,i))) = global_end(j,i) - global_start(j,i)  &
     4283                                                              + 1
     4284             lo_index(target_thread(j,i)) = lo_index(target_thread(j,i)) + 1
     4285          ENDDO
     4286       ENDDO
     4287
     4288       CALL MPI_WIN_FENCE( 0, win_ind, ierr )
     4289
     4290       buf_start = 1
     4291       DO  n = 0, numprocs-1
     4292          n_trans = transfer_index(2,n)
     4293          IF ( n_trans > 0 )  THEN
     4294             disp = 2 * ( transfer_index(1,n) - 1 )
     4295             CALL MPI_GET( local_indices(1,buf_start), 2*n_trans, MPI_INTEGER, n, disp, 2*n_trans, &
     4296                           MPI_INTEGER, win_ind, ierr )
     4297             buf_start = buf_start + n_trans
    39894298          ENDIF
    3990 !
    3991 !--       Get the size of shared memory window on all PEs
    3992           CALL MPI_BCAST( nr_iope, 1, MPI_INTEGER, 0, sm_io%comm_shared, ierr )
    3993           CALL sm_io%sm_allocate_shared( array_1d, 1, MAX( 1, nr_iope ), win_surf )
    3994        ELSE
    3995           nr_iope = nr_val
    3996        ENDIF
    3997 #else
    3998        nr_iope = nr_val
    3999 #endif
    4000 
    4001 !
    4002 !--    Check, if surface data exist on this PE
    4003        data_to_write = .TRUE.
    4004        IF ( total_number_of_surface_values == 0 )  THEN
    4005           data_to_write = .FALSE.
    4006           RETURN
    4007        ENDIF
    4008 
    4009        IF ( sm_io%iam_io_pe )  THEN
    4010 
    4011           all_pes_write = ( MINVAL( all_nr_val ) > 0 )
    4012 
    4013           IF ( all_pes_write )  THEN
    4014              dims1(1)  = total_number_of_surface_values
    4015              lize1(1)  = nr_iope
    4016              start1(1) = glo_start-1
    4017 
    4018 #if defined( __parallel )
    4019              IF ( total_number_of_surface_values > 0 )  THEN
    4020                  CALL MPI_TYPE_CREATE_SUBARRAY( 1, dims1, lize1, start1, MPI_ORDER_FORTRAN,        &
    4021                                                 MPI_REAL, ft_surf, ierr )
    4022                  CALL MPI_TYPE_COMMIT( ft_surf, ierr )
    4023              ENDIF
    4024 #endif
     4299       ENDDO
     4300
     4301       CALL MPI_WIN_FENCE( 0, win_ind, ierr )
     4302
     4303       buf_start = 1
     4304       DO  n = 0, numprocs-1
     4305          n_trans = transfer_index(2,n)
     4306          IF ( n_trans > 0 )  THEN
     4307             disp = transfer_index(1,n) - 1
     4308             buf_start = buf_start + n_trans
    40254309          ENDIF
    4026        ENDIF
    4027 
    4028     ENDIF
     4310       ENDDO
     4311
     4312       CALL MPI_WIN_FREE( win_ind, ierr )
     4313
     4314       DEALLOCATE( surf_val_index )
     4315
     4316    END SUBROUTINE get_remote_indices
     4317
     4318!--------------------------------------------------------------------------------------------------!
     4319! Description:
     4320! ------------
     4321!> Allocate memory and create window for one-sided communication (1-d INTEGER array)
     4322!--------------------------------------------------------------------------------------------------!
     4323    SUBROUTINE rd_alloc_rma_mem( array, idim, win )
     4324
     4325       IMPLICIT NONE
     4326
     4327       INTEGER(iwp), INTENT(IN)       ::  idim     !< Dimension of this 1-D array
     4328       INTEGER                        ::  ierr     !< MPI error code
     4329       INTEGER(iwp), INTENT(OUT)      ::  win      !< MPI window
     4330       INTEGER(KIND=MPI_ADDRESS_KIND) ::  winsize  !< size of RMA window
     4331
     4332       REAL(wp), DIMENSION(:), POINTER, INTENT(INOUT) ::  array  !< array to access RMA window locally
     4333
     4334
     4335       winsize = MAX( idim, 2 )
     4336       ALLOCATE( array(winsize) )
     4337       winsize = winsize * wp
     4338       CALL MPI_WIN_CREATE( array, winsize, wp, MPI_INFO_NULL, comm2d, win, ierr )
     4339       array = -1
     4340       CALL MPI_WIN_FENCE( 0, win, ierr )
     4341
     4342    END SUBROUTINE rd_alloc_rma_mem
     4343#endif
    40294344
    40304345 END SUBROUTINE rd_mpi_io_surface_filetypes
    4031 
    40324346
    40334347
     
    40794393          iog%nnx = iog%nnx + nbgp
    40804394       ENDIF
    4081        IF ( myidx == npex-1  .OR.  npex == -1 )  THEN   ! npex == 1 if -D__parallel not set
     4395       IF ( myidx == pdims(1)-1 )  THEN
    40824396          iog%nxr = iog%nxr + nbgp
    40834397          iog%nnx = iog%nnx + nbgp
     
    40874401          iog%nny = iog%nny + nbgp
    40884402       ENDIF
    4089        IF ( myidy == npey-1  .OR.  npey == -1 )  THEN   ! npey == 1 if -D__parallel not set
     4403       IF ( myidy == pdims(2)-1 )  THEN
    40904404          iog%nyn = iog%nyn + nbgp
    40914405          iog%nny = iog%nny + nbgp
     
    42514565          iog%nnx = iog%nnx + nbgp
    42524566       ENDIF
    4253        IF ( myidx == npex-1  .OR.  npex == -1 )  THEN   ! npex == 1 if -D__parallel not set
     4567       IF ( myidx == pdims(1)-1 )  THEN
    42544568          iog%nxr = iog%nxr + nbgp
    42554569          iog%nnx = iog%nnx + nbgp
     
    42594573          iog%nny = iog%nny + nbgp
    42604574       ENDIF
    4261        IF ( myidy == npey-1  .OR.  npey == -1 )  THEN   ! npey == 1 if -D__parallel not set
     4575       IF ( myidy == pdims(2)-1 )  THEN
    42624576          iog%nyn = iog%nyn + nbgp
    42634577          iog%nny = iog%nny + nbgp
     
    43264640!> to a single file that contains the global arrays. It is not required for the serial mode.
    43274641!--------------------------------------------------------------------------------------------------!
    4328 #if defined( __parallel )
    43294642 SUBROUTINE rd_mpi_io_create_filetypes_3dsoil( nzb_soil, nzt_soil )
    43304643
     
    43344647    INTEGER, INTENT(IN)   ::  nzt_soil  !<
    43354648
     4649#if defined( __parallel )
    43364650    INTEGER, DIMENSION(3) ::  dims3     !<
    43374651    INTEGER, DIMENSION(3) ::  lize3     !<
     
    43674681       CALL MPI_TYPE_COMMIT( ft_3dsoil, ierr )
    43684682    ENDIF
     4683#else
     4684    ALLOCATE( array_3d_soil(nzb_soil:nzt_soil,iog%nxl:iog%nxr,iog%nys:iog%nyn) )
     4685    sm_io%io_grid = iog
     4686#endif
    43694687
    43704688 END SUBROUTINE rd_mpi_io_create_filetypes_3dsoil
    4371 #endif
    4372 
    4373 
    43744689
    43754690!--------------------------------------------------------------------------------------------------!
     
    43814696
    43824697    IMPLICIT NONE
    4383 
    43844698
    43854699#if defined( __parallel )
     
    44014715
    44024716    ENDIF
     4717
    44034718!
    44044719!-- Free last surface filetype
     
    44154730    IF ( sm_io%iam_io_pe .AND. ft_3di4 /= -1 )  THEN
    44164731       CALL MPI_TYPE_FREE( ft_3di4, ierr )
     4732       ft_3di4 = -1
     4733    ENDIF
     4734    IF ( sm_io%iam_io_pe .AND. ft_3di8 /= -1 )  THEN
    44174735       CALL MPI_TYPE_FREE( ft_3di8, ierr )
     4736       ft_3di8 = -1
    44184737    ENDIF
    44194738
    44204739    IF ( sm_io%is_sm_active() .AND.  win_3di4 /= -1 )  THEN
    44214740       CALL sm_io%sm_free_shared( win_3di4 )
     4741       win_3di4 = -1
     4742    ENDIF
     4743    IF ( sm_io%is_sm_active() .AND.  win_3di8 /= -1 )  THEN
    44224744       CALL sm_io%sm_free_shared( win_3di8 )
     4745       win_3di8 = -1
     4746    ENDIF
     4747
     4748    IF ( win_start /= -1 )  THEN
     4749       CALL sm_io%sm_free_shared( win_start)
     4750       CALL sm_io%sm_free_shared( win_end)
     4751       CALL sm_io%sm_free_shared( win_glost)
     4752       win_start = -1
     4753       win_end   = -1
     4754       win_glost = -1
    44234755    ENDIF
    44244756
     
    44264758    win_surf = -1
    44274759#else
    4428     IF ( ASSOCIATED(array_2d)   )  DEALLOCATE( array_2d )
    4429     IF ( ASSOCIATED(array_2di)  )  DEALLOCATE( array_2di )
    4430     IF ( ASSOCIATED(array_3d)   )  DEALLOCATE( array_3d )
    4431     IF ( ASSOCIATED(array_3di4) )  DEALLOCATE( array_3di4 )
    4432     IF ( ASSOCIATED(array_3di8) )  DEALLOCATE( array_3di8 )
     4760    IF ( ASSOCIATED( array_2d )   )  DEALLOCATE( array_2d )
     4761    IF ( ASSOCIATED( array_2di )  )  DEALLOCATE( array_2di )
     4762    IF ( ASSOCIATED( array_3d )   )  DEALLOCATE( array_3d )
     4763    IF ( ASSOCIATED( array_3di4 ) )  DEALLOCATE( array_3di4 )
     4764    IF ( ASSOCIATED( array_3di8 ) )  DEALLOCATE( array_3di8 )
    44334765#endif
    44344766
  • palm/trunk/SOURCE/shared_memory_io_mod.f90

    r4828 r4893  
    2525! -----------------
    2626! $Id$
     27! revised output of surface data via MPI-IO for better performance
     28!
     29! 4828 2021-01-05 11:21:41Z Giersch
    2730! additions for output of particle time series
    2831!
     
    3942! unused variable removed
    4043!
    41 !
    4244! Additions for cyclic fill mode
    43 !
    4445!
    4546! File re-formatted to follow the PALM coding standard
     
    4748!
    4849! Initial version (Klaus Ketelsen)
    49 !
    50 !
    5150!
    5251! Description:
     
    10099              comm1dy,                                                                             &
    101100              comm2d,                                                                              &
     101              comm_palm,                                                                           &
    102102              ierr,                                                                                &
    103103              myid,                                                                                &
     
    120120#endif
    121121
     122    USE transpose_indices,                                                                         &
     123        ONLY:  nxl_y, nxl_z, nxr_y, nxr_z, nys_x, nys_z, nyn_x, nyn_z, nzb_x, nzb_y, nzt_x, nzt_y
     124
    122125    IMPLICIT NONE
    123126
     
    149152    END TYPE domain_decomposition_grid_features
    150153
     154    TYPE, PUBLIC ::  sm_remote_array
     155
     156       TYPE(C_PTR)  ::  rem_ptr  !<
     157       INTEGER(iwp) ::  d1e      !<
     158       INTEGER(iwp) ::  d1s      !<
     159       INTEGER(iwp) ::  d2e      !<
     160       INTEGER(iwp) ::  d2s      !<
     161       INTEGER(iwp) ::  d3e      !<
     162       INTEGER(iwp) ::  d3s      !<
     163       INTEGER(iwp) ::  d4e      !<
     164       INTEGER(iwp) ::  d4s      !<
     165
     166    END TYPE sm_remote_array
     167
    151168!
    152169!-- Class definition for shared memory instances.
     
    164181       INTEGER(iwp), PUBLIC ::  sh_rank       !<
    165182
    166        LOGICAL, PUBLIC ::  iam_io_pe = .TRUE.  !< This PE is an IO-PE
    167183!
    168184!--    Variables for the I/O virtual grid
    169        INTEGER(iwp), PUBLIC ::  comm_io  !< Communicator for all IO processes
     185       INTEGER(iwp), PUBLIC ::  comm_io  !< communicator for all IO processes
    170186       INTEGER(iwp), PUBLIC ::  io_npes  !<
    171187       INTEGER(iwp), PUBLIC ::  io_rank  !<
    172188!
    173189!--    Variables for the node local communicator
    174        INTEGER(iwp) ::  comm_node          !< Communicator for all processes of current node
     190       INTEGER(iwp) ::  comm_node          !< communicator for all processes of current node
    175191       INTEGER(iwp) ::  io_pe_global_rank  !<
    176192       INTEGER(iwp) ::  n_npes             !<
    177193       INTEGER(iwp) ::  n_rank             !<
    178194
    179        TYPE(domain_decomposition_grid_features), PUBLIC ::  io_grid  !< io grid features, depending on reading from prerun or restart run
    180 
     195       LOGICAL, PUBLIC ::  is_root_pe          !<
     196       LOGICAL, PUBLIC ::  iam_io_pe = .TRUE.  !< this PE is an IO-PE
     197
     198       TYPE(domain_decomposition_grid_features), PUBLIC ::  io_grid  !< io grid features, depending on reading from prerun or main run
    181199
    182200       CONTAINS
     
    191209          PROCEDURE, PASS(this), PUBLIC ::  sm_node_barrier
    192210#if defined( __parallel )
     211          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_1d_32
    193212          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_1d_64
    194           PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_1d_32
    195213          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_1di
     214          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_2d_32
    196215          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_2d_64
    197           PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_2d_32
    198216          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_2di
     217          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_3d_32
    199218          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_3d_64
    200           PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_3d_32
     219          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_4d_32
     220          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_4d_64
    201221          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_3di_32
    202222          PROCEDURE, PASS(this), PUBLIC ::  sm_allocate_shared_3di_64
     223          PROCEDURE, PASS(this), PUBLIC ::  sm_all_allocate_shared_3d_64
    203224
    204225          GENERIC, PUBLIC ::  sm_allocate_shared =>                                                &
    205                                               sm_allocate_shared_1d_64,  sm_allocate_shared_1d_32, &
    206                                               sm_allocate_shared_2d_64,  sm_allocate_shared_2d_32, &
    207                                               sm_allocate_shared_2di,    sm_allocate_shared_3d_64, &
    208                                               sm_allocate_shared_3d_32,  sm_allocate_shared_1di,   &
    209                                               sm_allocate_shared_3di_32, sm_allocate_shared_3di_64
     226                                             sm_allocate_shared_1d_64,  sm_allocate_shared_1d_32,  &
     227                                             sm_allocate_shared_2d_64,  sm_allocate_shared_2d_32,  &
     228                                             sm_allocate_shared_2di,    sm_allocate_shared_3d_64,  &
     229                                             sm_allocate_shared_4d_64,  sm_allocate_shared_4d_32,  &
     230                                             sm_allocate_shared_3d_32,  sm_allocate_shared_1di,    &
     231                                             sm_allocate_shared_3di_32, sm_allocate_shared_3di_64
     232
     233          GENERIC, PUBLIC ::  sm_all_allocate_shared => sm_all_allocate_shared_3d_64
    210234#endif
    211235    END TYPE sm_class
     
    226250
    227251    CLASS(sm_class), INTENT(INOUT) ::  this        !< pointer to access internal variables of this call
    228     INTEGER, INTENT(IN), OPTIONAL ::  comm_input  !< main model communicator (comm2d) can optional be set
     252    INTEGER(iwp), INTENT(IN), OPTIONAL ::  comm_input  !< main model communicator (comm2d) can optional be set
    229253
    230254#if defined( __parallel )
    231     INTEGER ::  color
    232     INTEGER ::  max_n_npes  !< maximum number of PEs/node
     255    INTEGER ::  color              !<
     256    INTEGER ::  max_npes_per_node  !< maximum number of PEs/node
    233257#endif
    234258
     
    237261    this%nr_io_pe_per_node = 2
    238262
     263#if defined( __parallel )
    239264    IF ( PRESENT( comm_input ) )  THEN
    240265       this%comm_model = comm_input
     
    248273    IF ( this%no_shared_memory_in_this_run )  THEN
    249274       this%iam_io_pe = .TRUE.
     275       this%sh_rank   = 0
     276       this%sh_npes   = 1
    250277       RETURN
    251278    ENDIF
    252279
    253 #if defined( __parallel )
    254 !
    255 !-- Determine, how many MPI threads are running on a node
     280!
     281!-- Determine, how many PEs are running on a node.
    256282    this%iam_io_pe = .FALSE.
    257283    CALL MPI_COMM_SPLIT_TYPE( this%comm_model, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL,             &
     
    260286    CALL MPI_COMM_RANK( this%comm_node, this%n_rank, ierr )
    261287
    262     CALL MPI_ALLREDUCE( this%n_npes, max_n_npes, 1, MPI_INTEGER, MPI_MAX, this%comm_model, ierr )
     288    CALL MPI_ALLREDUCE( this%n_npes, max_npes_per_node, 1, MPI_INTEGER, MPI_MAX, this%comm_model,  &
     289                        ierr )
    263290!
    264291!-- Decide, if the configuration can run with shared-memory IO
    265     IF ( max_n_npes > 64 )  THEN
     292    IF ( max_npes_per_node > 64 )  THEN
    266293!
    267294!--    Special configuration on the HLRN-IV system with 4 shared memory blocks/node
    268295       this%nr_io_pe_per_node = 4
    269296
    270     ELSEIF ( max_n_npes <= 32 )  THEN
    271 !
    272 !--    No shared memory IO with less than 32 threads/node
     297    ELSEIF ( max_npes_per_node <= 3 )  THEN
     298!
     299!--    No shared memory IO with less than 3 MPI tasks/node
    273300       this%no_shared_memory_in_this_run = .TRUE.
    274301       this%iam_io_pe = .TRUE.
     
    277304
    278305!
    279 !-- No shared memory IO with small setups
    280     IF ( nx < 24  .OR.  ny < 24 )  THEN
     306!-- No shared memory IO with small setups.
     307    IF ( nx < 16  .OR.  ny < 16 )  THEN
    281308       this%no_shared_memory_in_this_run = .TRUE.
    282309       this%iam_io_pe = .TRUE.
     
    299326!
    300327!-- Setup the communicator across the nodes depending on the shared memory rank.
    301 !-- All threads with shared memory rank 0 will be I/O threads.
     328!-- All PEs with shared memory rank 0 will be I/O PEs.
    302329    color = this%sh_rank
    303330    CALL MPI_COMM_SPLIT( this%comm_model, color, 0, this%comm_io, ierr )
     
    316343    ENDIF
    317344    CALL MPI_BCAST( this%io_pe_global_rank, 1, MPI_INTEGER, 0, this%comm_shared, ierr )
    318 
    319345#else
    320     this%iam_io_pe = .TRUE.
     346    this%iam_io_pe  = .TRUE.
     347    this%comm_model = comm2d
     348    this%sh_rank    = 0
     349    this%sh_npes    = 1
     350    this%no_shared_memory_in_this_run = .TRUE.
    321351#endif
    322 
    323 !      write(9,'(a,8i7)') ' end of sm_init_comm ',this%sh_rank,this%sh_npes,this%io_rank,this%io_npes,this%io_pe_global_rank
    324 !      write(9,*) 'This process is IO Process ',this%iam_io_pe
    325352
    326353#if defined( __parallel )
     
    362389    CALL MPI_ALLREDUCE( local_dim_s, local_dim_r, SIZE( local_dim_s ), MPI_INTEGER, MPI_SUM,       &
    363390                        this%comm_node, ierr )
    364     sh_group_size = ( max_n_npes + this%nr_io_pe_per_node - 1 ) / this%nr_io_pe_per_node
     391    sh_group_size = ( max_npes_per_node + this%nr_io_pe_per_node - 1 ) / this%nr_io_pe_per_node
    365392
    366393    pe       = 0
     
    417444 END SUBROUTINE sm_init_comm
    418445
    419 
    420446!
    421447!-- Initializing setup for output of particle time series.
     
    428454
    429455#if defined( __parallel )
    430     INTEGER(iwp) ::  color             !<
    431     INTEGER(iwp) ::  ierr              !<
    432     INTEGER(iwp) ::  max_n_npes        !< maximum number of PEs/node
     456    INTEGER(iwp) ::  color              !<
     457    INTEGER(iwp) ::  ierr               !<
     458    INTEGER(iwp) ::  max_npes_per_node  !< maximum number of PEs/node
    433459#endif
    434460
     
    451477#if defined( __parallel )
    452478!
    453 !-- Determine, how many MPI threads are running on a node
     479!-- Determine, how many PEs are running on a node.
    454480    this%iam_io_pe = .FALSE.
    455481    CALL MPI_COMM_SPLIT_TYPE( this%comm_model, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL,             &
     
    458484    CALL MPI_COMM_RANK( this%comm_node, this%n_rank, ierr )
    459485
    460     CALL MPI_ALLREDUCE( this%n_npes, max_n_npes, 1, MPI_INTEGER, MPI_MAX, this%comm_model, ierr )
    461 
     486    CALL MPI_ALLREDUCE( this%n_npes, max_npes_per_node, 1, MPI_INTEGER, MPI_MAX, this%comm_model,  &
     487                        ierr )
    462488!
    463489!-- TODO: better explanation
     
    465491!-- even better to use the complete node for MPI shared memory (this%nr_io_pe_per_node = 1).
    466492!-  In the latter case, the access to the MPI shared memory buffer is slower, the number of
    467 !-- particles to move between threads will be much smaller.
    468     IF ( max_n_npes > 64 )  THEN
     493!-- particles to move between PEs will be much smaller.
     494    IF ( max_npes_per_node > 64 )  THEN
    469495!
    470496!--    Special configuration on the HLRN-IV system with 4 shared memory blocks/node
     
    526552!
    527553!-- Setup the communicator across the nodes depending on the shared memory rank.
    528 !-- All threads with shared memory rank 0 will be I/O threads.
     554!-- All PEs with shared memory rank 0 will be I/O PEs.
    529555    color = this%sh_rank
    530556    CALL MPI_COMM_SPLIT( this%comm_model, color, 0, this%comm_io, ierr )
     
    573599! Description:
    574600! ------------
    575 !> Allocate shared 1d-REAL (64 Bit) array on ALL threads
     601!> Allocate shared 1d-REAL (64 bit) array on PE 0 and pass address to all PEs.
    576602!--------------------------------------------------------------------------------------------------!
    577603 SUBROUTINE sm_allocate_shared_1d_64( this, p1, d1, d2, win )
     
    590616    INTEGER(KIND=MPI_ADDRESS_KIND)  ::  wsize
    591617
    592     INTEGER, DIMENSION(1)           ::  buf_shape
     618    INTEGER(iwp), DIMENSION(1)      ::  buf_shape
    593619
    594620    REAL(dp), DIMENSION(:), POINTER ::  buf
     
    601627    IF ( this%no_shared_memory_in_this_run )  RETURN
    602628!
    603 !-- Allocate shared memory on node rank 0 threads.
     629!-- Allocate shared memory on node rank 0 PEs.
    604630    IF ( this%sh_rank == pe_from )  THEN
    605631       wsize = d2 - d1 + 1
     
    629655! Description:
    630656! ------------
    631 !> Allocate shared 1d-REAL (32 Bit) array on ALL threads
     657!> Allocate shared 1d-REAL (32 bit) array on PE 0 and pass address to all PEs
    632658!--------------------------------------------------------------------------------------------------!
    633659 SUBROUTINE sm_allocate_shared_1d_32( this, p1, d1, d2, win )
     
    646672    INTEGER(KIND=MPI_ADDRESS_KIND)  ::  wsize
    647673
    648     INTEGER, DIMENSION(1)           ::  buf_shape
     674    INTEGER(iwp), DIMENSION(1)      ::  buf_shape
    649675
    650676    REAL(sp), DIMENSION(:), POINTER ::  buf
     
    657683    IF ( this%no_shared_memory_in_this_run )  RETURN
    658684!
    659 !-- Allocate shared memory on node rank 0 threads.
     685!-- Allocate shared memory on node rank 0 PEs.
    660686    IF ( this%sh_rank == pe_from )  THEN
    661687       wsize = d2 - d1 + 1
     
    685711! Description:
    686712! ------------
    687 !> Allocate shared 1d-INTEGER array on ALL threads
     713!> Allocate shared 1d-INTEGER array on PE 0 and pass address to all PEs.
    688714!--------------------------------------------------------------------------------------------------!
    689715 SUBROUTINE sm_allocate_shared_1di( this, p1, d1, d2, win )
     
    702728    INTEGER(KIND=MPI_ADDRESS_KIND)  ::  wsize
    703729
    704     INTEGER, DIMENSION(1)           ::  buf_shape
     730    INTEGER(iwp), DIMENSION(1)          ::  buf_shape
    705731
    706732    INTEGER(iwp), DIMENSION(:), POINTER ::  buf
     
    713739    IF ( this%no_shared_memory_in_this_run )  RETURN
    714740!
    715 !-- Allocate shared memory on node rank 0 threads.
     741!-- Allocate shared memory on node rank 0 PEs.
    716742    IF ( this%sh_rank == pe_from )  THEN
    717743       wsize = d2 - d1 + 1
     
    741767! Description:
    742768! ------------
    743 !> Allocate shared 2d-REAL array on ALL threads (64 Bit)
     769!> Allocate shared 2d-REAL array (64 bit) on PE 0 and pass address to all PEs.
    744770!--------------------------------------------------------------------------------------------------!
    745771 SUBROUTINE sm_allocate_shared_2d_64( this, p2, n_nxlg, n_nxrg, n_nysg, n_nyng, win )
     
    771797    IF ( this%no_shared_memory_in_this_run )  RETURN
    772798!
    773 !-- Allocate shared memory on node rank 0 threads.
     799!-- Allocate shared memory on node rank 0 PEs.
    774800    IF ( this%sh_rank == pe_from )  THEN
    775801       wsize = ( n_nyng - n_nysg + 1 ) * ( n_nxrg - n_nxlg + 1 )
     
    801827! Description:
    802828! ------------
    803 !> Allocate shared 2d-REAL (32 Bit) array on ALL threads
     829!> Allocate shared 2d-REAL (32 Bit) array on PE 0 and pass address to all PEs.
    804830!--------------------------------------------------------------------------------------------------!
    805831 SUBROUTINE sm_allocate_shared_2d_32( this, p2, n_nxlg, n_nxrg, n_nysg, n_nyng, win )
     
    831857    IF ( this%no_shared_memory_in_this_run )  RETURN
    832858!
    833 !-- Allocate shared memory on node rank 0 threads.
     859!-- Allocate shared memory on node rank 0 PEs.
    834860    IF ( this%sh_rank == pe_from )  THEN
    835861       wsize = ( n_nyng - n_nysg + 1 ) * ( n_nxrg - n_nxlg + 1 )
     
    861887! Description:
    862888! ------------
    863 !> Allocate shared 2d-INTEGER array on ALL threads
     889!> Allocate shared 2d-INTEGER array on PE 0 and pass address to all PEs.
    864890!--------------------------------------------------------------------------------------------------!
    865891 SUBROUTINE sm_allocate_shared_2di( this, p2i, n_nxlg, n_nxrg, n_nysg, n_nyng, win )
     
    891917    IF ( this%no_shared_memory_in_this_run )  RETURN
    892918!
    893 !-- Allocate shared memory on node rank 0 threads.
     919!-- Allocate shared memory on node rank 0 PEs.
    894920    IF ( this%sh_rank == pe_from )  THEN
    895921       wsize = ( n_nyng - n_nysg + 1 ) * ( n_nxrg - n_nxlg + 1 )
     
    921947! Description:
    922948! ------------
    923 !> Allocate shared 3d-REAL (64 Bit) array on ALL threads
     949!> Allocate shared 3d-REAL (64 bit) array on PE 0 and pass address to all PEs.
    924950!--------------------------------------------------------------------------------------------------!
    925951 SUBROUTINE sm_allocate_shared_3d_64( this, p3, d1s, d1e, d2s, d2e, d3s, d3e, win )
     
    929955    CLASS(sm_class), INTENT(inout)      ::  this         !<
    930956
    931     INTEGER                             ::  disp_unit    !<
    932     INTEGER, INTENT(IN)                 ::  d1e          !<
    933     INTEGER, INTENT(IN)                 ::  d1s          !<
    934     INTEGER, INTENT(IN)                 ::  d2e          !<
    935     INTEGER, INTENT(IN)                 ::  d2s          !<
    936     INTEGER, INTENT(IN)                 ::  d3e          !<
    937     INTEGER, INTENT(IN)                 ::  d3s          !<
    938     INTEGER, SAVE                       ::  pe_from = 0  !<
    939     INTEGER, INTENT(OUT)                ::  win          !<
     957    INTEGER(iwp)                        ::  disp_unit    !<
     958    INTEGER(iwp), INTENT(IN)            ::  d1e          !<
     959    INTEGER(iwp), INTENT(IN)            ::  d1s          !<
     960    INTEGER(iwp), INTENT(IN)            ::  d2e          !<
     961    INTEGER(iwp), INTENT(IN)            ::  d2s          !<
     962    INTEGER(iwp), INTENT(IN)            ::  d3e          !<
     963    INTEGER(iwp), INTENT(IN)            ::  d3s          !<
     964    INTEGER(iwp), SAVE                  ::  pe_from = 0  !<
     965    INTEGER(iwp), INTENT(OUT)           ::  win          !<
    940966
    941967    INTEGER(KIND=MPI_ADDRESS_KIND)      ::  rem_size     !<
    942968    INTEGER(KIND=MPI_ADDRESS_KIND)      ::  wsize        !<
    943969
    944     INTEGER, DIMENSION(3)               ::  buf_shape    !<
     970    INTEGER(iwp), DIMENSION(3)          ::  buf_shape    !<
    945971
    946972    REAL(dp), DIMENSION(:,:,:), POINTER ::  buf          !<
     
    953979    IF ( this%no_shared_memory_in_this_run )  RETURN
    954980!
    955 !-- Allocate shared memory on node rank 0 threads.
     981!-- Allocate shared memory on node rank 0 PEs.
    956982    IF ( this%sh_rank == pe_from )  THEN
    957983       wsize = ( d3e - d3s + 1 ) * ( d2e - d2s + 1 ) * ( d1e - d1s + 1 )
     
    9841010! Description:
    9851011! ------------
    986 !> Allocate shared 3d-REAL (32 Bit) array on ALL threads
     1012!> Allocate shared 3d-REAL (32 bit) array on PE 0 and pass address to all PEs.
    9871013!--------------------------------------------------------------------------------------------------!
    9881014 SUBROUTINE sm_allocate_shared_3d_32( this, p3, d1s, d1e, d2s, d2e, d3s, d3e, win )
     
    9921018    CLASS(sm_class), INTENT(inout)      ::  this
    9931019
    994     INTEGER                             ::  disp_unit
    995     INTEGER, INTENT(IN)                 ::  d1e
    996     INTEGER, INTENT(IN)                 ::  d1s
    997     INTEGER, INTENT(IN)                 ::  d2e
    998     INTEGER, INTENT(IN)                 ::  d2s
    999     INTEGER, INTENT(IN)                 ::  d3e
    1000     INTEGER, INTENT(IN)                 ::  d3s
    1001     INTEGER, SAVE                       ::  pe_from = 0
    1002     INTEGER, INTENT(OUT)                ::  win
     1020    INTEGER(iwp)                        ::  disp_unit
     1021    INTEGER(iwp), INTENT(IN)            ::  d1e
     1022    INTEGER(iwp), INTENT(IN)            ::  d1s
     1023    INTEGER(iwp), INTENT(IN)            ::  d2e
     1024    INTEGER(iwp), INTENT(IN)            ::  d2s
     1025    INTEGER(iwp), INTENT(IN)            ::  d3e
     1026    INTEGER(iwp), INTENT(IN)            ::  d3s
     1027    INTEGER(iwp), SAVE                  ::  pe_from = 0
     1028    INTEGER(iwp), INTENT(OUT)           ::  win
    10031029
    10041030    INTEGER(KIND=MPI_ADDRESS_KIND)      ::  rem_size
    10051031    INTEGER(KIND=MPI_ADDRESS_KIND)      ::  wsize
    10061032
    1007     INTEGER, DIMENSION(3)               ::  buf_shape
     1033    INTEGER(iwp), DIMENSION(3)          ::  buf_shape
    10081034
    10091035    REAL(sp), DIMENSION(:,:,:), POINTER ::  buf
     
    10161042    IF ( this%no_shared_memory_in_this_run )  RETURN
    10171043!
    1018 !-- Allocate shared memory on node rank 0 threads.
     1044!-- Allocate shared memory on node rank 0 PEs.
    10191045    IF ( this%sh_rank == pe_from )  THEN
    10201046       wsize = ( d3e - d3s + 1 ) * ( d2e - d2s + 1 ) * ( d1e - d1s + 1 )
     
    10471073! Description:
    10481074! ------------
    1049 !> Allocate shared 3d-REAL (32 bit) array on ALL threads
     1075!> Allocate shared 4d-REAL (64 bit) array on PE 0 and pass address to all PEs.
     1076!--------------------------------------------------------------------------------------------------!
     1077 SUBROUTINE sm_allocate_shared_4d_64( this, p3, d1s, d1e, d2s, d2e, d3s, d3e, d4s, d4e, win )
     1078
     1079    IMPLICIT NONE
     1080
     1081    CLASS(sm_class), INTENT(inout)        ::  this         !<
     1082
     1083    INTEGER                               ::  disp_unit    !<
     1084    INTEGER(iwp), INTENT(IN)              ::  d1e          !<
     1085    INTEGER(iwp), INTENT(IN)              ::  d1s          !<
     1086    INTEGER(iwp), INTENT(IN)              ::  d2e          !<
     1087    INTEGER(iwp), INTENT(IN)              ::  d2s          !<
     1088    INTEGER(iwp), INTENT(IN)              ::  d3e          !<
     1089    INTEGER(iwp), INTENT(IN)              ::  d3s          !<
     1090    INTEGER(iwp), INTENT(IN)              ::  d4e          !<
     1091    INTEGER(iwp), INTENT(IN)              ::  d4s          !<
     1092    INTEGER(iwp), SAVE                    ::  pe_from = 0  !<
     1093    INTEGER(iwp), INTENT(OUT)             ::  win          !<
     1094
     1095    INTEGER(KIND=MPI_ADDRESS_KIND)        ::  rem_size     !<
     1096    INTEGER(KIND=MPI_ADDRESS_KIND)        ::  wsize        !<
     1097
     1098    INTEGER(iwp), DIMENSION(4)            ::  buf_shape    !<
     1099
     1100    REAL(dp), DIMENSION(:,:,:,:), POINTER ::  buf          !<
     1101    REAL(dp), DIMENSION(:,:,:,:), POINTER ::  p3           !<
     1102
     1103    TYPE(C_PTR), SAVE                     ::  base_ptr     !<
     1104    TYPE(C_PTR), SAVE                     ::  rem_ptr      !<
     1105
     1106
     1107    IF ( this%no_shared_memory_in_this_run )  RETURN
     1108!
     1109!-- Allocate shared memory on node rank 0 PEs.
     1110    IF ( this%sh_rank == pe_from )  THEN
     1111       wsize = (d4e - d4s +1) * ( d3e - d3s + 1 ) * ( d2e - d2s + 1 ) * ( d1e - d1s + 1 )
     1112    ELSE
     1113       wsize = 1
     1114    ENDIF
     1115
     1116    wsize = wsize * dp ! Please note, size is always in bytes, independently of the displacement
     1117                       ! unit
     1118
     1119    CALL MPI_WIN_ALLOCATE_SHARED( wsize, dp, MPI_INFO_NULL, this%comm_shared, base_ptr, win, ierr )
     1120!
     1121!-- Get C-pointer of the memory located on node-rank pe_from (sh_rank == pe_from)
     1122    CALL MPI_WIN_SHARED_QUERY( win, pe_from, rem_size, disp_unit, rem_ptr, ierr )
     1123!
     1124!-- Convert C- to Fortran-pointer
     1125    buf_shape(4) = d4e - d4s + 1
     1126    buf_shape(3) = d3e - d3s + 1
     1127    buf_shape(2) = d2e - d2s + 1
     1128    buf_shape(1) = d1e - d1s + 1
     1129    CALL C_F_POINTER( rem_ptr, buf, buf_shape )
     1130    p3(d1s:,d2s:,d3s:,d4s:) => buf
     1131!
     1132!-- Allocate shared memory in round robin on all PEs of a node.
     1133    pe_from = MOD( pe_from, this%sh_npes )
     1134
     1135 END SUBROUTINE sm_allocate_shared_4d_64
     1136
     1137
     1138!--------------------------------------------------------------------------------------------------!
     1139! Description:
     1140! ------------
     1141!> Allocate shared 4d-REAL (32 bit) array on PE 0 and pass address to all PEs.
     1142!--------------------------------------------------------------------------------------------------!
     1143 SUBROUTINE sm_allocate_shared_4d_32( this, p3, d1s, d1e, d2s, d2e, d3s, d3e, d4s, d4e, win )
     1144
     1145    IMPLICIT NONE
     1146
     1147    CLASS(sm_class), INTENT(inout)        ::  this         !<
     1148
     1149    INTEGER                               ::  disp_unit    !<
     1150    INTEGER(iwp), INTENT(IN)              ::  d1e          !<
     1151    INTEGER(iwp), INTENT(IN)              ::  d1s          !<
     1152    INTEGER(iwp), INTENT(IN)              ::  d2e          !<
     1153    INTEGER(iwp), INTENT(IN)              ::  d2s          !<
     1154    INTEGER(iwp), INTENT(IN)              ::  d3e          !<
     1155    INTEGER(iwp), INTENT(IN)              ::  d3s          !<
     1156    INTEGER(iwp), INTENT(IN)              ::  d4e          !<
     1157    INTEGER(iwp), INTENT(IN)              ::  d4s          !<
     1158    INTEGER(iwp), SAVE                    ::  pe_from = 0  !<
     1159    INTEGER(iwp), INTENT(OUT)             ::  win          !<
     1160
     1161    INTEGER(KIND=MPI_ADDRESS_KIND)        ::  rem_size     !<
     1162    INTEGER(KIND=MPI_ADDRESS_KIND)        ::  wsize        !<
     1163
     1164    INTEGER(iwp), DIMENSION(4)            ::  buf_shape    !<
     1165
     1166    REAL(sp), DIMENSION(:,:,:,:), POINTER ::  buf          !<
     1167    REAL(sp), DIMENSION(:,:,:,:), POINTER ::  p3           !<
     1168
     1169    TYPE(C_PTR), SAVE                     ::  base_ptr     !<
     1170    TYPE(C_PTR), SAVE                     ::  rem_ptr      !<
     1171
     1172
     1173    IF ( this%no_shared_memory_in_this_run )  RETURN
     1174!
     1175!-- Allocate shared memory on node rank 0 PEs.
     1176    IF ( this%sh_rank == pe_from )  THEN
     1177       wsize = (d4e - d4s +1) * ( d3e - d3s + 1 ) * ( d2e - d2s + 1 ) * ( d1e - d1s + 1 )
     1178    ELSE
     1179       wsize = 1
     1180    ENDIF
     1181
     1182    wsize = wsize * sp ! Please note, size is always in bytes, independently of the displacement
     1183                       ! unit
     1184
     1185    CALL MPI_WIN_ALLOCATE_SHARED( wsize, dp, MPI_INFO_NULL, this%comm_shared, base_ptr, win, ierr )
     1186!
     1187!-- Get C-pointer of the memory located on node-rank pe_from (sh_rank == pe_from)
     1188    CALL MPI_WIN_SHARED_QUERY( win, pe_from, rem_size, disp_unit, rem_ptr, ierr )
     1189!
     1190!-- Convert C- to Fortran-pointer
     1191    buf_shape(4) = d4e - d4s + 1
     1192    buf_shape(3) = d3e - d3s + 1
     1193    buf_shape(2) = d2e - d2s + 1
     1194    buf_shape(1) = d1e - d1s + 1
     1195    CALL C_F_POINTER( rem_ptr, buf, buf_shape )
     1196    p3(d1s:,d2s:,d3s:,d4s:) => buf
     1197!
     1198!-- Allocate shared memory in round robin on all PEs of a node.
     1199    pe_from = MOD( pe_from, this%sh_npes )
     1200
     1201 END SUBROUTINE sm_allocate_shared_4d_32
     1202
     1203
     1204!--------------------------------------------------------------------------------------------------!
     1205! Description:
     1206! ------------
     1207!> Allocate shared 3d-INTEGER (32 bit) array on PE 0 and pass address to all PEs.
    10501208!--------------------------------------------------------------------------------------------------!
    10511209 SUBROUTINE sm_allocate_shared_3di_32( this, p3, d1s, d1e, d2s, d2e, d3s, d3e, win )
     
    10531211    IMPLICIT NONE
    10541212
    1055     CLASS(sm_class), INTENT(inout)      ::  this
    1056 
    1057     INTEGER                             ::  disp_unit
    1058     INTEGER, INTENT(IN)                 ::  d1e
    1059     INTEGER, INTENT(IN)                 ::  d1s
    1060     INTEGER, INTENT(IN)                 ::  d2e
    1061     INTEGER, INTENT(IN)                 ::  d2s
    1062     INTEGER, INTENT(IN)                 ::  d3e
    1063     INTEGER, INTENT(IN)                 ::  d3s
    1064     INTEGER, SAVE                       ::  pe_from = 0
    1065     INTEGER, INTENT(OUT)                ::  win
    1066 
    1067     INTEGER(KIND=MPI_ADDRESS_KIND)      ::  rem_size
    1068     INTEGER(KIND=MPI_ADDRESS_KIND)      ::  wsize
    1069 
    1070     INTEGER, DIMENSION(3)               ::  buf_shape
     1213    CLASS(sm_class), INTENT(inout)          ::  this
     1214
     1215    INTEGER                                 ::  disp_unit
     1216    INTEGER(iwp), INTENT(IN)                ::  d1e
     1217    INTEGER(iwp), INTENT(IN)                ::  d1s
     1218    INTEGER(iwp), INTENT(IN)                ::  d2e
     1219    INTEGER(iwp), INTENT(IN)                ::  d2s
     1220    INTEGER(iwp), INTENT(IN)                ::  d3e
     1221    INTEGER(iwp), INTENT(IN)                ::  d3s
     1222    INTEGER(iwp), SAVE                      ::  pe_from = 0
     1223    INTEGER(iwp), INTENT(OUT)               ::  win
     1224
     1225    INTEGER(KIND=MPI_ADDRESS_KIND)          ::  rem_size
     1226    INTEGER(KIND=MPI_ADDRESS_KIND)          ::  wsize
     1227
     1228    INTEGER(iwp), DIMENSION(3)              ::  buf_shape
    10711229
    10721230    INTEGER(isp), DIMENSION(:,:,:), POINTER ::  buf
    10731231    INTEGER(isp), DIMENSION(:,:,:), POINTER ::  p3
    10741232
    1075     TYPE(C_PTR), SAVE                   ::  base_ptr
    1076     TYPE(C_PTR), SAVE                   ::  rem_ptr
    1077 
    1078 
    1079     IF ( this%no_shared_memory_in_this_run )  RETURN
    1080 !
    1081 !-- Allocate shared memory on node rank 0 threads.
     1233    TYPE(C_PTR), SAVE                       ::  base_ptr
     1234    TYPE(C_PTR), SAVE                       ::  rem_ptr
     1235
     1236
     1237    IF ( this%no_shared_memory_in_this_run )  RETURN
     1238!
     1239!-- Allocate shared memory on node rank 0 PEs.
    10821240    IF ( this%sh_rank == pe_from )  THEN
    10831241       wsize = ( d3e - d3s + 1 ) * ( d2e - d2s + 1 ) * ( d1e - d1s + 1 )
     
    11101268! Description:
    11111269! ------------
    1112 !> Allocate shared 3d-REAL (64 bit) array on ALL threads
     1270!> Allocate shared 3d-INTEGER (64 bit) array on PE 0 and pass address to all PEs.
    11131271!--------------------------------------------------------------------------------------------------!
    11141272 SUBROUTINE sm_allocate_shared_3di_64( this, p3, d1s, d1e, d2s, d2e, d3s, d3e, win )
     
    11161274    IMPLICIT NONE
    11171275
    1118     CLASS(sm_class), INTENT(inout)      ::  this         !<
    1119 
    1120     INTEGER                             ::  disp_unit    !<
    1121     INTEGER, INTENT(IN)                 ::  d1e          !<
    1122     INTEGER, INTENT(IN)                 ::  d1s          !<
    1123     INTEGER, INTENT(IN)                 ::  d2e          !<
    1124     INTEGER, INTENT(IN)                 ::  d2s          !<
    1125     INTEGER, INTENT(IN)                 ::  d3e          !<
    1126     INTEGER, INTENT(IN)                 ::  d3s          !<
    1127     INTEGER, SAVE                       ::  pe_from = 0  !<
    1128     INTEGER, INTENT(OUT)                ::  win          !<
    1129 
    1130     INTEGER(KIND=MPI_ADDRESS_KIND)      ::  rem_size     !<
    1131     INTEGER(KIND=MPI_ADDRESS_KIND)      ::  wsize        !<
    1132 
    1133     INTEGER, DIMENSION(3)               ::  buf_shape    !<
     1276    CLASS(sm_class), INTENT(inout)          ::  this         !<
     1277
     1278    INTEGER                                 ::  disp_unit    !<
     1279    INTEGER(iwp), INTENT(IN)                ::  d1e          !<
     1280    INTEGER(iwp), INTENT(IN)                ::  d1s          !<
     1281    INTEGER(iwp), INTENT(IN)                ::  d2e          !<
     1282    INTEGER(iwp), INTENT(IN)                ::  d2s          !<
     1283    INTEGER(iwp), INTENT(IN)                ::  d3e          !<
     1284    INTEGER(iwp), INTENT(IN)                ::  d3s          !<
     1285    INTEGER(iwp), SAVE                      ::  pe_from = 0  !<
     1286    INTEGER(iwp), INTENT(OUT)               ::  win          !<
     1287
     1288    INTEGER(KIND=MPI_ADDRESS_KIND)          ::  rem_size     !<
     1289    INTEGER(KIND=MPI_ADDRESS_KIND)          ::  wsize        !<
     1290
     1291    INTEGER(iwp), DIMENSION(3)              ::  buf_shape    !<
    11341292
    11351293    INTEGER(idp), DIMENSION(:,:,:), POINTER ::  buf          !<
    11361294    INTEGER(idp), DIMENSION(:,:,:), POINTER ::  p3           !<
    11371295
    1138     TYPE(C_PTR), SAVE                   ::  base_ptr     !<
    1139     TYPE(C_PTR), SAVE                   ::  rem_ptr      !<
    1140 
    1141 
    1142     IF ( this%no_shared_memory_in_this_run )  RETURN
    1143 !
    1144 !-- Allocate shared memory on node rank 0 threads.
     1296    TYPE(C_PTR), SAVE                       ::  base_ptr     !<
     1297    TYPE(C_PTR), SAVE                       ::  rem_ptr      !<
     1298
     1299
     1300    IF ( this%no_shared_memory_in_this_run )  RETURN
     1301!
     1302!-- Allocate shared memory on node rank 0 PEs.
    11451303    IF ( this%sh_rank == pe_from )  THEN
    11461304       wsize = ( d3e - d3s + 1 ) * ( d2e - d2s + 1 ) * ( d1e - d1s + 1 )
     
    11691327 END SUBROUTINE sm_allocate_shared_3di_64
    11701328
     1329
     1330!--------------------------------------------------------------------------------------------------!
     1331! Description:
     1332! ------------
     1333!> Allocate shared 3d-REAL (64 Bit) array on ALL PEs.
     1334!>
     1335!> Every PE allocates the local part of a node-shared array.
     1336!> The C-Pointer of this array and the local limits are broadcasted to all PEs of the node
     1337!> The information is store in an array of type sm_remote_array and can be retrieved
     1338!> by sm_remote_array to access remote data.
     1339!--------------------------------------------------------------------------------------------------!
     1340 SUBROUTINE sm_all_allocate_shared_3d_64( this, p3, d1s, d1e, d2s, d2e, d3s, d3e, remote_arrays, win )
     1341
     1342    IMPLICIT NONE
     1343
     1344    CLASS(sm_class), INTENT(inout)      ::  this         !< class pointer
     1345    REAL(dp), DIMENSION(:,:,:), POINTER ::  p3           !< return local array pointer
     1346
     1347    INTEGER(iwp), INTENT(IN)            ::  d1e          !< end index dimension 1
     1348    INTEGER(iwp), INTENT(IN)            ::  d1s          !< start index dimension 1
     1349    INTEGER(iwp), INTENT(IN)            ::  d2e          !< end index dimension 2
     1350    INTEGER(iwp), INTENT(IN)            ::  d2s          !< start index dimension 2
     1351    INTEGER(iwp), INTENT(IN)            ::  d3e          !< end index dimension 3
     1352    INTEGER(iwp), INTENT(IN)            ::  d3s          !< start index dimension 3
     1353    INTEGER(iwp), INTENT(OUT)           ::  win          !< MPI Window
     1354
     1355    INTEGER(iwp), DIMENSION(3)          ::  buf_shape    !<
     1356    INTEGER(iwp)                        ::  disp_unit    !<
     1357    INTEGER(iwp)                        ::  i            !<
     1358    INTEGER(iwp), SAVE                  ::  pe_from = 0  !<
     1359
     1360    INTEGER(KIND=MPI_ADDRESS_KIND)      ::  rem_size     !<
     1361    INTEGER(KIND=MPI_ADDRESS_KIND)      ::  wsize        !<
     1362
     1363    REAL(dp), DIMENSION(:,:,:), POINTER ::  buf          !<
     1364
     1365    TYPE(sm_remote_array),INTENT(INOUT), DIMENSION(0:this%sh_npes-1) :: remote_arrays !< info about all remote arrays
     1366
     1367    TYPE(C_PTR), SAVE                   ::  base_ptr     !<
     1368
     1369    INTEGER(iwp),DIMENSION(6,0:this%sh_npes-1)              ::  all_indices_s
     1370    INTEGER(iwp),DIMENSION(6,0:this%sh_npes-1)              ::  all_indices
     1371
     1372
     1373    IF ( this%no_shared_memory_in_this_run )  RETURN
     1374
     1375    all_indices_s = 0
     1376
     1377
     1378    wsize = ( d3e - d3s + 1 ) * ( d2e - d2s + 1 ) * ( d1e - d1s + 1 )
     1379
     1380    wsize = wsize * dp   ! Please note, size is always in bytes, independently of the displacement unit
     1381
     1382    CALL MPI_WIN_ALLOCATE_SHARED( wsize, dp, MPI_INFO_NULL, this%comm_shared, base_ptr, win, ierr )
     1383!
     1384!-- Get C-pointer of the memory located on node-rank pe_from (sh_rank == pe_from)
     1385
     1386    all_indices_s(1,this%sh_rank) = d1s
     1387    all_indices_s(2,this%sh_rank) = d1e
     1388    all_indices_s(3,this%sh_rank) = d2s
     1389    all_indices_s(4,this%sh_rank) = d2e
     1390    all_indices_s(5,this%sh_rank) = d3s
     1391    all_indices_s(6,this%sh_rank) = d3e
     1392
     1393    CALL MPI_ALLREDUCE (all_indices_s ,all_indices, SIZE(all_indices_s), MPI_INTEGER, MPI_SUM, this%comm_shared, ierr)
     1394
     1395    DO i=0,this%sh_npes-1
     1396       CALL MPI_WIN_SHARED_QUERY( win, i, rem_size, disp_unit, remote_arrays(i)%rem_ptr, ierr )
     1397       remote_arrays(i)%d1s = all_indices(1,i)
     1398       remote_arrays(i)%d1e = all_indices(2,i)
     1399       remote_arrays(i)%d2s = all_indices(3,i)
     1400       remote_arrays(i)%d2e = all_indices(4,i)
     1401       remote_arrays(i)%d3s = all_indices(5,i)
     1402       remote_arrays(i)%d3e = all_indices(6,i)
     1403    END DO
     1404
     1405!
     1406!-- Convert C- to Fortran-pointer
     1407    buf_shape(3) = d3e - d3s + 1
     1408    buf_shape(2) = d2e - d2s + 1
     1409    buf_shape(1) = d1e - d1s + 1
     1410    CALL C_F_POINTER( remote_arrays(this%sh_rank)%rem_ptr, buf, buf_shape )
     1411    p3(d1s:,d2s:,d3s:) => buf
     1412!
     1413!-- Allocate shared memory in round robin on all PEs of a node.
     1414    pe_from = MOD( pe_from, this%sh_npes )
     1415
     1416 END SUBROUTINE sm_all_allocate_shared_3d_64
    11711417#endif
    11721418
     
    12431489!> ...
    12441490!--------------------------------------------------------------------------------------------------!
    1245  SUBROUTINE sm_node_barrier( this )
    1246 
    1247     IMPLICIT NONE
     1491 SUBROUTINE sm_node_barrier( this, win )
     1492
     1493    IMPLICIT NONE
     1494
     1495    INTEGER(iwp), OPTIONAL         ::  win   !<
    12481496
    12491497    CLASS(sm_class), INTENT(inout) ::  this  !<
     
    12541502#if defined( __parallel )
    12551503    CALL MPI_BARRIER( this%comm_shared, ierr )
     1504    IF ( PRESENT(win) )  THEN
     1505       CALL MPI_WIN_FENCE(0, win, ierr )
     1506    ENDIF
    12561507#endif
    12571508
  • palm/trunk/SOURCE/surface_data_output_mod.f90

    r4892 r4893  
    2525! -----------------
    2626! $Id$
     27! revised output of surface data via MPI-IO for better performance
     28!
     29! 4892 2021-03-02 11:53:58Z suehring
    2730! Remove outdated error message.
    2831!
     
    45784581
    45794582    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  end_index           !< end index of surface data at (j,i)
    4580     INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  global_start_index  !< index array for surface data (MPI-IO)
     4583    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  global_end_index    !< end index array for surface data (MPI-IO)
     4584    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  global_start_index  !< start index array for surface data (MPI-IO)
    45814585    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  num_surf            !< number of surface data at (j,i)
    45824586    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  start_index         !< start index of surface data at (j,i)
     
    46004604    ALLOCATE( surf_in(1:surfaces%ns) )
    46014605
    4602     CALL rd_mpi_io_check_array( 'surfaces%start_index', found = array_found )
    4603     IF ( array_found )  CALL rrd_mpi_io( 'surfaces%start_index', start_index )
    4604 
    4605     CALL rd_mpi_io_check_array( 'surfaces%end_index', found = array_found )
    4606     IF ( array_found )  CALL rrd_mpi_io( 'surfaces%end_index', end_index )
    4607 
    46084606    CALL rd_mpi_io_check_array( 'surfaces%global_start_index', found = array_found )
    46094607    IF ( array_found )  CALL rrd_mpi_io( 'surfaces%global_start_index', global_start_index )
    46104608
    4611     CALL rd_mpi_io_surface_filetypes( start_index, end_index, ldum, global_start_index )
     4609    CALL rd_mpi_io_check_array( 'surfaces%global_end_index', found = array_found )
     4610    IF ( array_found )  CALL rrd_mpi_io( 'surfaces%global_end_index', global_end_index )
     4611
     4612    CALL rd_mpi_io_surface_filetypes( start_index, end_index, ldum, global_start_index,            &
     4613                                      global_end_index )
    46124614
    46134615    DO  nv = 1, dosurf_no(1)
     
    47284730
    47294731    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  end_index           !< end index of surface data at (j,i)
    4730     INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  global_start_index  !< index array for surface data (MPI-IO)
     4732    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  global_end_index    !< end index array for surface data (MPI-IO)
     4733    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  global_start_index  !< start index array for surface data (MPI-IO)
    47314734    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  num_surf            !< number of surface data at (j,i)
    47324735    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE ::  start_index         !< start index of surface data at (j,i)
     
    48114814
    48124815          CALL rd_mpi_io_surface_filetypes( start_index, end_index, surface_data_to_write,         &
    4813                                             global_start_index )
    4814           CALL wrd_mpi_io( 'surfaces%start_index',        start_index        )
    4815           CALL wrd_mpi_io( 'surfaces%end_index',          end_index          )
     4816                                            global_start_index, global_end_index )
     4817
    48164818          CALL wrd_mpi_io( 'surfaces%global_start_index', global_start_index )
     4819          CALL wrd_mpi_io( 'surfaces%global_end_index', global_end_index )
    48174820
    48184821          DO  nv = 1, dosurf_no(1)
  • palm/trunk/SOURCE/surface_mod.f90

    r4882 r4893  
    2525! -----------------
    2626! $Id$
     27! revised output of surface data via MPI-IO for better performance
     28!
     29! 4882 2021-02-19 22:49:44Z forkel
    2730! removed lsp in subroutine nitialize_top
    28 !
    2931!
    3032! 4881 2021-02-19 22:05:08Z forkel
    3133! removed constant_top_csflux option
    32 !
    3334!
    3435! 4877 2021-02-17 16:17:35Z suehring
     
    14841485!> Allocating memory for upward and downward-facing horizontal surface types, except for top fluxes.
    14851486!--------------------------------------------------------------------------------------------------!
    1486  SUBROUTINE allocate_surface_attributes_h( surfaces, nys_l, nyn_l, nxl_l, nxr_l )
     1487 SUBROUTINE allocate_surface_attributes_h( surfaces, nys_l, nyn_l, nxl_l, nxr_l,                   &
     1488                                           no_allocate_index_arrays )
    14871489
    14881490    IMPLICIT NONE
     
    14931495    INTEGER(iwp) ::  nxr_l  !< east bound of local 2d array start/end_index, is equal to nyn, except for restart-array
    14941496
     1497    LOGICAL ::  allocate_index_arrays
     1498    LOGICAL, INTENT(IN), OPTIONAL  :: no_allocate_index_arrays
     1499
    14951500    TYPE(surf_type) ::  surfaces  !< respective surface type
    14961501
     1502
     1503    IF ( PRESENT( no_allocate_index_arrays ) )  THEN
     1504       allocate_index_arrays = .NOT. no_allocate_index_arrays
     1505    ELSE
     1506       allocate_index_arrays = .TRUE.
     1507    ENDIF
    14971508!
    14981509!-- Allocate arrays for start and end index of horizontal surface type for each (j,i)-grid point.
    1499 !-- This is required e.g. in diffion_x, which is called for each (j,i). In order to find the
     1510!-- This is required e.g. in diffusion_x, which is called for each (j,i). In order to find the
    15001511!-- location where the respective flux is store within the surface-type, start- and end-index are
    15011512!-- stored for each (j,i). For example, each (j,i) can have several entries where fluxes for
     
    15031514!-- surfaces might exist for given (j,i). If no surface of respective type exist at current (j,i),
    15041515!-- set indicies such that loop in diffusion routines will not be entered.
    1505     ALLOCATE ( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
    1506     ALLOCATE ( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l) )
    1507     surfaces%start_index = 0
    1508     surfaces%end_index   = -1
     1516    IF ( allocate_index_arrays )  THEN
     1517       ALLOCATE( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
     1518       ALLOCATE( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l)   )
     1519       surfaces%start_index = 0
     1520       surfaces%end_index   = -1
     1521    ENDIF
    15091522!
    15101523!-- Indices to locate surface element
     
    20512064!> Allocating memory for vertical surface types.
    20522065!--------------------------------------------------------------------------------------------------!
    2053  SUBROUTINE allocate_surface_attributes_v( surfaces, nys_l, nyn_l, nxl_l, nxr_l )
     2066 SUBROUTINE allocate_surface_attributes_v( surfaces, nys_l, nyn_l, nxl_l, nxr_l,                   &
     2067                                           no_allocate_index_arrays )
    20542068
    20552069    IMPLICIT NONE
     
    20602074    INTEGER(iwp) ::  nxr_l  !< east bound of local 2d array start/end_index, is equal to nyn, except for restart-array
    20612075
     2076    LOGICAL ::  allocate_index_arrays
     2077    LOGICAL, INTENT(IN), OPTIONAL ::  no_allocate_index_arrays
     2078
    20622079    TYPE(surf_type) ::  surfaces !< respective surface type
    20632080
     2081
     2082    IF ( PRESENT( no_allocate_index_arrays ) )  THEN
     2083       allocate_index_arrays = .NOT. no_allocate_index_arrays
     2084    ELSE
     2085       allocate_index_arrays = .TRUE.
     2086    ENDIF
     2087
    20642088!
    20652089!-- Allocate arrays for start and end index of vertical surface type for each (j,i)-grid point. This
    2066 !-- is required in diffion_x, which is called for each (j,i). In order to find the location where
     2090!-- is required in diffusion_x, which is called for each (j,i). In order to find the location where
    20672091!-- the respective flux is store within the surface-type, start- and end-index are stored for each
    20682092!-- (j,i). For example, each (j,i) can have several entries where fluxes for vertical surfaces might
    20692093!-- be stored. In the flat case, where no vertical walls exit, set indicies such that loop in
    20702094!-- diffusion routines will not be entered.
    2071     ALLOCATE ( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
    2072     ALLOCATE ( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l) )
    2073     surfaces%start_index = 0
    2074     surfaces%end_index   = -1
     2095    IF ( allocate_index_arrays )  THEN
     2096       ALLOCATE( surfaces%start_index(nys_l:nyn_l,nxl_l:nxr_l) )
     2097       ALLOCATE( surfaces%end_index(nys_l:nyn_l,nxl_l:nxr_l)   )
     2098       surfaces%start_index = 0
     2099       surfaces%end_index   = -1
     2100    ENDIF
    20752101!
    20762102!-- Indices to locate surface element.
     
    31413167    INTEGER(iwp), DIMENSION(0:3) ::  start_index_v  !< start index for vertical surface elements on gathered surface array
    31423168
    3143     INTEGER(iwp),DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index  !< index for surface data (MPI-IO)
     3169    INTEGER(iwp),DIMENSION(nys:nyn,nxl:nxr) ::  global_end_index    !< end index for surface data (MPI-IO)
     3170    INTEGER(iwp),DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index  !< start index for surface data (MPI-IO)
    31443171
    31453172    LOGICAL ::  surface_data_to_write  !< switch for MPI-I/O if PE has surface data to write
     
    39884015
    39894016          CALL rd_mpi_io_surface_filetypes( surf_h(l)%start_index, surf_h(l)%end_index,            &
    3990                                             surface_data_to_write, global_start_index )
     4017                                            surface_data_to_write, global_start_index,             &
     4018                                            global_end_index )
    39914019          IF ( .NOT. surface_data_to_write )  CYCLE
    39924020
    39934021          ns_h_on_file(l) = total_number_of_surface_values
    39944022
    3995           CALL wrd_mpi_io( 'surf_h(' // dum // ')%start_index', surf_h(l)%start_index )
    3996           CALL wrd_mpi_io( 'surf_h(' // dum // ')%end_index', surf_h(l)%end_index )
    39974023          CALL wrd_mpi_io( 'global_start_index_h_' // dum, global_start_index )
     4024          CALL wrd_mpi_io( 'global_end_index_h_' // dum, global_end_index )
    39984025
    39994026          IF ( ALLOCATED ( surf_h(l)%us ) )  THEN
     
    41214148
    41224149          CALL rd_mpi_io_surface_filetypes( surf_v(l)%start_index, surf_v(l)%end_index,            &
    4123                                             surface_data_to_write, global_start_index )
     4150                                            surface_data_to_write, global_start_index,             &
     4151                                            global_end_index )
     4152
     4153          IF ( .NOT. surface_data_to_write )  CYCLE
    41244154
    41254155          ns_v_on_file(l) = total_number_of_surface_values
    41264156
    4127           CALL wrd_mpi_io( 'surf_v(' // dum // ')%start_index', surf_v(l)%start_index )
    4128           CALL wrd_mpi_io( 'surf_v(' // dum // ')%end_index', surf_v(l)%end_index )
    41294157          CALL wrd_mpi_io( 'global_start_index_v_' // dum, global_start_index )
     4158          CALL wrd_mpi_io( 'global_end_index_v_' // dum, global_end_index )
    41304159
    41314160          IF ( .NOT. surface_data_to_write )  CYCLE
     
    53815410    INTEGER(iwp) ::  mm !< loop index for surface types - file array
    53825411
    5383     INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index  !< index for surface data (MPI-IO)
    5384 
    5385     LOGICAL ::  ldum            !< dummy variable
     5412    INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_end_index    !< end index for surface data (MPI-IO)
     5413    INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index  !< start index for surface data (MPI-IO)
     5414
     5415    LOGICAL ::  data_to_read    !< cycle in l loop, if no values to read
    53865416    LOGICAL ::  surf_match_def  !< flag indicating that surface element is of default type
    53875417    LOGICAL ::  surf_match_lsm  !< flag indicating that surface element is of natural type
     
    54015431       IF ( ns_h_on_file(l) == 0 )  CYCLE  !< No data of this surface type on file
    54025432
     5433       WRITE( dum, '(I1)')  l
     5434
    54035435       IF ( ALLOCATED( surf_h(l)%start_index ) )  CALL deallocate_surface_attributes_h( surf_h(l) )
    5404        surf_h(l)%ns = ns_h_on_file(l)
    5405        CALL allocate_surface_attributes_h( surf_h(l), nys, nyn, nxl, nxr )
    5406 
    5407        WRITE( dum, '(I1)') l
    5408 
    5409        CALL rrd_mpi_io( 'surf_h(' // dum // ')%start_index',  surf_h(l)%start_index )
    5410        CALL rrd_mpi_io( 'surf_h(' // dum // ')%end_index',  surf_h(l)%end_index )
     5436
     5437       ALLOCATE( surf_h(l)%start_index(nys:nyn,nxl:nxr) )
     5438       ALLOCATE( surf_h(l)%end_index(nys:nyn,nxl:nxr) )
     5439       surf_h(l)%start_index = 0
     5440       surf_h(l)%end_index   = -1
     5441
    54115442       CALL rrd_mpi_io( 'global_start_index_h_' // dum , global_start_index )
    5412 
    5413        CALL rd_mpi_io_surface_filetypes( surf_h(l)%start_index, surf_h(l)%end_index, ldum,         &
    5414                                          global_start_index )
     5443       CALL rrd_mpi_io( 'global_end_index_h_' // dum , global_end_index )
     5444
     5445       CALL rd_mpi_io_surface_filetypes( surf_h(l)%start_index, surf_h(l)%end_index, data_to_read, &
     5446                                         global_start_index, global_end_index )
     5447
     5448       surf_h(l)%ns = MAX( 2, MAXVAL( surf_h(l)%end_index ) )
     5449
     5450       CALL allocate_surface_attributes_h( surf_h(l), nys, nyn, nxl, nxr,                          &
     5451                                           no_allocate_index_arrays = .TRUE. )
     5452       IF ( .NOT. data_to_read )  CYCLE
    54155453
    54165454       IF ( ALLOCATED ( surf_h(l)%us ) )  THEN
     
    55395577
    55405578       IF ( ALLOCATED( surf_v(l)%start_index ) )  CALL deallocate_surface_attributes_v( surf_v(l) )
    5541        surf_v(l)%ns = ns_v_on_file(l)
    5542        CALL allocate_surface_attributes_v( surf_v(l), nys, nyn, nxl, nxr )
    5543 
    5544        WRITE( dum, '(I1)' ) l
    5545 
    5546        CALL rrd_mpi_io( 'surf_v(' // dum // ')%start_index', surf_v(l)%start_index )
    5547        CALL rrd_mpi_io( 'surf_v(' // dum // ')%end_index', surf_v(l)%end_index )
     5579
     5580       ALLOCATE( surf_v(l)%start_index(nys:nyn,nxl:nxr) )
     5581       ALLOCATE( surf_v(l)%end_index(nys:nyn,nxl:nxr) )
     5582       surf_v(l)%start_index = 0
     5583       surf_v(l)%end_index   = -1
     5584
     5585       WRITE( dum, '(I1)' )  l
     5586
    55485587       CALL rrd_mpi_io( 'global_start_index_v_' // dum , global_start_index )
    5549 
    5550        CALL rd_mpi_io_surface_filetypes( surf_v(l)%start_index, surf_v(l)%end_index, ldum,         &
    5551                                          global_start_index )
     5588       CALL rrd_mpi_io( 'global_end_index_v_' // dum , global_end_index )
     5589
     5590       CALL rd_mpi_io_surface_filetypes( surf_v(l)%start_index, surf_v(l)%end_index, data_to_read, &
     5591                                         global_start_index, global_end_index )
     5592       IF ( .NOT. data_to_read )  CYCLE
     5593
     5594       surf_v(l)%ns = MAX( 2, MAXVAL( surf_v(l)%end_index ) )
     5595       CALL allocate_surface_attributes_v( surf_v(l), nys, nyn, nxl, nxr,                          &
     5596                                           no_allocate_index_arrays = .TRUE.  )
    55525597
    55535598       IF ( ALLOCATED ( surf_v(l)%us ) )  THEN
  • palm/trunk/SOURCE/urban_surface_mod.f90

    r4872 r4893  
    2727! -----------------
    2828! $Id$
     29! revised output of surface data via MPI-IO for better performance
     30!
     31! 4872 2021-02-12 15:49:02Z raasch
    2932! internal switch removed from namelist
    3033!
     
    66086611 SUBROUTINE usm_rrd_local_mpi
    66096612
    6610 
    66116613    CHARACTER(LEN=1) ::  dum  !< dummy string to create input-variable name
    66126614
    66136615    INTEGER(iwp) ::  l  !< loop index for surface types
    66146616
    6615     INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_start
    6616 
    6617     LOGICAL ::  ldum  !< dummy variable
     6617    INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_end_index
     6618    INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index
     6619
     6620    LOGICAL ::  data_to_read  !< dummy variable
     6621
    66186622
    66196623    DO  l = 0, 1
     
    66216625       WRITE( dum, '(I1)' )  l
    66226626
    6623        CALL rrd_mpi_io( 'usm_start_index_h_' //dum,  surf_usm_h(l)%start_index )
    6624        CALL rrd_mpi_io( 'usm_end_index_h_' //dum, surf_usm_h(l)%end_index )
    6625        CALL rrd_mpi_io( 'usm_global_start_h_' //dum, global_start )
    6626 
    6627        CALL rd_mpi_io_surface_filetypes( surf_usm_h(l)%start_index, surf_usm_h(l)%end_index, ldum, &
    6628                                          global_start )
    6629 
    6630        IF ( MAXVAL( surf_usm_h(l)%end_index ) <= 0 )  CYCLE
    6631 
    6632        IF ( .NOT.  ALLOCATED( t_surf_wall_h_1(l)%val ) )                                             &
     6627       CALL rrd_mpi_io( 'usm_global_start_h_' //dum, global_start_index )
     6628       CALL rrd_mpi_io( 'usm_global_end_h_' //dum, global_end_index )
     6629
     6630       CALL rd_mpi_io_surface_filetypes( surf_usm_h(l)%start_index, surf_usm_h(l)%end_index,       &
     6631                                         data_to_read, global_start_index, global_end_index )
     6632       IF ( .NOT. data_to_read )  CYCLE
     6633
     6634       IF ( .NOT. ALLOCATED( t_surf_wall_h_1(l)%val ) )                                            &
    66336635          ALLOCATE( t_surf_wall_h_1(l)%val(1:surf_usm_h(l)%ns) )
    66346636       CALL rrd_mpi_io_surface( 't_surf_wall_h(' // dum // ')', t_surf_wall_h_1(l)%val )
    66356637
    6636        IF ( .NOT.  ALLOCATED( t_surf_window_h_1(l)%val ) )                                           &
     6638       IF ( .NOT. ALLOCATED( t_surf_window_h_1(l)%val ) )                                          &
    66376639          ALLOCATE( t_surf_window_h_1(l)%val(1:surf_usm_h(l)%ns) )
    66386640       CALL rrd_mpi_io_surface( 't_surf_window_h(' // dum // ')', t_surf_window_h_1(l)%val )
    66396641
    6640        IF ( .NOT.  ALLOCATED( t_surf_green_h_1(l)%val ) )                                            &
     6642       IF ( .NOT. ALLOCATED( t_surf_green_h_1(l)%val ) )                                           &
    66416643          ALLOCATE( t_surf_green_h_1(l)%val(1:surf_usm_h(l)%ns) )
    66426644       CALL rrd_mpi_io_surface( 't_surf_green_h(' // dum // ')', t_surf_green_h_1(l)%val )
    66436645
    6644        IF ( .NOT.  ALLOCATED( m_liq_usm_h_1(l)%val ) )                                             &
     6646       IF ( .NOT. ALLOCATED( m_liq_usm_h_1(l)%val ) )                                              &
    66456647          ALLOCATE( m_liq_usm_h_1(l)%val(1:surf_usm_h(l)%ns) )
    66466648       CALL rrd_mpi_io_surface( 'm_liq_usm_h(' // dum // ')', m_liq_usm_h_1(l)%val )
    66476649
    66486650       IF ( indoor_model )  THEN
    6649           IF ( .NOT.  ALLOCATED( surf_usm_h(l)%waste_heat ) )                                      &
     6651          IF ( .NOT. ALLOCATED( surf_usm_h(l)%waste_heat ) )                                       &
    66506652             ALLOCATE( surf_usm_h(l)%waste_heat(1:surf_usm_h(l)%ns) )
    66516653          CALL rrd_mpi_io_surface( 'waste_heat_h(' // dum // ')', surf_usm_h(l)%waste_heat )
    6652           IF ( .NOT.  ALLOCATED( surf_usm_h(l)%t_prev ) )                                          &
     6654          IF ( .NOT. ALLOCATED( surf_usm_h(l)%t_prev ) )                                           &
    66536655             ALLOCATE( surf_usm_h(l)%t_prev(1:surf_usm_h(l)%ns) )
    66546656          CALL rrd_mpi_io_surface( 't_prev_h(' // dum // ')', surf_usm_h(l)%t_prev )
     
    66606662       WRITE( dum, '(I1)' )  l
    66616663
    6662        CALL rrd_mpi_io( 'usm_start_index_v_' //dum, surf_usm_v(l)%start_index )
    6663        CALL rrd_mpi_io( 'usm_end_index_v_' // dum, surf_usm_v(l)%end_index )
    6664        CALL rrd_mpi_io( 'usm_global_start_v_' // dum, global_start )
    6665 
    6666        CALL rd_mpi_io_surface_filetypes( surf_usm_v(l)%start_index, surf_usm_v(l)%end_index, ldum, &
    6667                                          global_start )
    6668 
    6669        IF ( MAXVAL( surf_usm_v(l)%end_index ) <= 0 )  CYCLE
    6670 
    6671        IF ( .NOT.  ALLOCATED( t_surf_wall_v_1(l)%val ) )                                             &
     6664       CALL rrd_mpi_io( 'usm_global_start_v_' // dum, global_start_index )
     6665       CALL rrd_mpi_io( 'usm_global_end_v_' // dum, global_end_index )
     6666
     6667       CALL rd_mpi_io_surface_filetypes( surf_usm_v(l)%start_index, surf_usm_v(l)%end_index,       &
     6668                                         data_to_read, global_start_index, global_end_index )
     6669       IF ( .NOT. data_to_read )  CYCLE
     6670
     6671       IF ( .NOT. ALLOCATED( t_surf_wall_v_1(l)%val ) )                                            &
    66726672          ALLOCATE( t_surf_wall_v_1(l)%val(1:surf_usm_v(l)%ns) )
    66736673       CALL rrd_mpi_io_surface( 't_surf_wall_v(' // dum // ')', t_surf_wall_v_1(l)%val )
    66746674
    6675        IF ( .NOT.  ALLOCATED( t_surf_window_v_1(l)%val ) )                                           &
     6675       IF ( .NOT. ALLOCATED( t_surf_window_v_1(l)%val ) )                                          &
    66766676          ALLOCATE( t_surf_window_v_1(l)%val(1:surf_usm_v(l)%ns) )
    66776677       CALL rrd_mpi_io_surface( 't_surf_window_v(' // dum // ')', t_surf_window_v_1(l)%val )
    66786678
    6679        IF ( .NOT.  ALLOCATED( t_surf_green_v_1(l)%val ) )                                            &
     6679       IF ( .NOT. ALLOCATED( t_surf_green_v_1(l)%val ) )                                           &
    66806680          ALLOCATE( t_surf_green_v_1(l)%val(1:surf_usm_v(l)%ns) )
    66816681       CALL rrd_mpi_io_surface( 't_surf_green_v(' // dum // ')', t_surf_green_v_1(l)%val)
    66826682
    66836683       IF ( indoor_model )  THEN
    6684           IF ( .NOT.  ALLOCATED( surf_usm_v(l)%waste_heat ) )                                      &
     6684          IF ( .NOT. ALLOCATED( surf_usm_v(l)%waste_heat ) )                                       &
    66856685             ALLOCATE( surf_usm_v(l)%waste_heat(1:surf_usm_v(l)%ns) )
    66866686          CALL rrd_mpi_io_surface( 'waste_heat_v(' // dum // ')', surf_usm_v(l)%waste_heat )
    6687           IF ( .NOT.  ALLOCATED( surf_usm_v(l)%t_prev ) )                                          &
     6687          IF ( .NOT. ALLOCATED( surf_usm_v(l)%t_prev ) )                                           &
    66886688             ALLOCATE( surf_usm_v(l)%t_prev(1:surf_usm_v(l)%ns) )
    66896689          CALL rrd_mpi_io_surface( 't_prev_v(' // dum // ')', surf_usm_v(l)%t_prev )
     
    66966696       WRITE( dum, '(I1)' )  l
    66976697
    6698        CALL rrd_mpi_io( 'usm_start_index_h_2_' //dum,  surf_usm_h(l)%start_index )
    6699        CALL rrd_mpi_io( 'usm_end_index_h_2_' //dum, surf_usm_h(l)%end_index )
    6700        CALL rrd_mpi_io( 'usm_global_start_h_2_' //dum, global_start )
    6701 
    6702        CALL rd_mpi_io_surface_filetypes( surf_usm_h(l)%start_index, surf_usm_h(l)%end_index, ldum,          &
    6703                                          global_start )
    6704 
    6705        IF ( MAXVAL( surf_usm_h(l)%end_index ) <= 0 )  CYCLE
    6706 
    6707        IF ( .NOT.  ALLOCATED( t_wall_h_1(l)%val ) )                                                          &
     6698       CALL rrd_mpi_io( 'usm_global_start_h_2_' //dum, global_start_index )
     6699       CALL rrd_mpi_io( 'usm_global_end_h_2_' //dum, global_end_index )
     6700
     6701       CALL rd_mpi_io_surface_filetypes( surf_usm_h(l)%start_index, surf_usm_h(l)%end_index,       &
     6702                                         data_to_read, global_start_index, global_end_index )
     6703       IF ( .NOT. data_to_read )  CYCLE
     6704
     6705       IF ( .NOT. ALLOCATED( t_wall_h_1(l)%val ) )                                                 &
    67086706          ALLOCATE( t_wall_h_1(l)%val(nzb_wall:nzt_wall+1,1:surf_usm_h(l)%ns) )
    67096707       CALL rrd_mpi_io_surface( 't_wall_h(' // dum // ')', t_wall_h_1(l)%val )
    67106708
    6711        IF ( .NOT.  ALLOCATED( t_window_h_1(l)%val ) )                                                        &
     6709       IF ( .NOT. ALLOCATED( t_window_h_1(l)%val ) )                                               &
    67126710          ALLOCATE( t_window_h_1(l)%val(nzb_wall:nzt_wall+1,1:surf_usm_h(l)%ns) )
    67136711       CALL rrd_mpi_io_surface( 't_window_h(' // dum // ')', t_window_h_1(l)%val )
    67146712
    6715        IF ( .NOT.  ALLOCATED( t_green_h_1(l)%val ) )                                                         &
     6713       IF ( .NOT. ALLOCATED( t_green_h_1(l)%val ) )                                                &
    67166714          ALLOCATE( t_green_h_1(l)%val(nzb_wall:nzt_wall+1,1:surf_usm_h(l)%ns) )
    67176715       CALL rrd_mpi_io_surface( 't_green_h(' // dum // ')', t_green_h_1(l)%val )
     
    67236721       WRITE( dum, '(I1)' )  l
    67246722
    6725        CALL rrd_mpi_io( 'usm_start_index_v_2_' //dum, surf_usm_v(l)%start_index )
    6726        CALL rrd_mpi_io( 'usm_end_index_v_2_' // dum, surf_usm_v(l)%end_index )
    6727        CALL rrd_mpi_io( 'usm_global_start_v_2_' // dum, global_start )
    6728 
    6729        CALL rd_mpi_io_surface_filetypes( surf_usm_v(l)%start_index, surf_usm_v(l)%end_index, ldum, &
    6730                                          global_start )
    6731 
    6732        IF ( MAXVAL( surf_usm_v(l)%end_index ) <= 0 )  CYCLE
    6733 
    6734        IF ( .NOT. ALLOCATED( t_wall_v_1(l)%val ) )                                                   &
     6723       CALL rrd_mpi_io( 'usm_global_start_v_2_' // dum, global_start_index )
     6724       CALL rrd_mpi_io( 'usm_global_end_v_2_' // dum, global_end_index )
     6725
     6726       CALL rd_mpi_io_surface_filetypes( surf_usm_v(l)%start_index, surf_usm_v(l)%end_index,       &
     6727                                         data_to_read, global_start_index, global_end_index )
     6728       IF ( .NOT. data_to_read )  CYCLE
     6729
     6730       IF ( .NOT. ALLOCATED( t_wall_v_1(l)%val ) )                                                 &
    67356731          ALLOCATE ( t_wall_v_1(l)%val(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
    67366732       CALL rrd_mpi_io_surface( 't_wall_v(' // dum // ')', t_wall_v_1(l)%val )
    67376733
    6738        IF ( .NOT. ALLOCATED( t_window_v_1(l)%val ) )                                                 &
     6734       IF ( .NOT. ALLOCATED( t_window_v_1(l)%val ) )                                               &
    67396735          ALLOCATE ( t_window_v_1(l)%val(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
    67406736       CALL rrd_mpi_io_surface( 't_window_v(' // dum // ')', t_window_v_1(l)%val )
    67416737
    6742        IF ( .NOT. ALLOCATED( t_green_v_1(l)%val ) )                                                  &
     6738       IF ( .NOT. ALLOCATED( t_green_v_1(l)%val ) )                                                &
    67436739          ALLOCATE ( t_green_v_1(l)%val(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
    67446740       CALL rrd_mpi_io_surface( 't_green_v(' // dum // ')', t_green_v_1(l)%val )
     
    67476743
    67486744 END SUBROUTINE usm_rrd_local_mpi
     6745
    67496746
    67506747!--------------------------------------------------------------------------------------------------!
     
    74227419    INTEGER(iwp)  ::  l  !< index surface type orientation
    74237420
    7424     INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr)  ::  global_start_index  !< index for surface data (MPI-IO)
     7421    INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_end_index    !< end index for surface data (MPI-IO)
     7422    INTEGER(iwp), DIMENSION(nys:nyn,nxl:nxr) ::  global_start_index  !< start index for surface data (MPI-IO)
    74257423
    74267424    LOGICAL  ::  surface_data_to_write  !< switch for MPI-I/O if PE has surface data to write
     
    75487546          WRITE( dum, '(I1)')  l
    75497547
    7550           CALL rd_mpi_io_surface_filetypes( surf_usm_h(l)%start_index, surf_usm_h(l)%end_index,             &
    7551                                             surface_data_to_write, global_start_index )
    7552 
    7553           CALL wrd_mpi_io( 'usm_start_index_h_' // dum,  surf_usm_h(l)%start_index )
    7554           CALL wrd_mpi_io( 'usm_end_index_h_' // dum, surf_usm_h(l)%end_index )
     7548          CALL rd_mpi_io_surface_filetypes( surf_usm_h(l)%start_index, surf_usm_h(l)%end_index,    &
     7549                                            surface_data_to_write, global_start_index,             &
     7550                                            global_end_index )
     7551
    75557552          CALL wrd_mpi_io( 'usm_global_start_h_' // dum, global_start_index )
     7553          CALL wrd_mpi_io( 'usm_global_end_h_' // dum, global_end_index )
    75567554
    75577555          IF ( .NOT. surface_data_to_write )  CYCLE
     
    75747572
    75757573          CALL rd_mpi_io_surface_filetypes( surf_usm_v(l)%start_index, surf_usm_v(l)%end_index,    &
    7576                                             surface_data_to_write, global_start_index )
    7577 
    7578           CALL wrd_mpi_io( 'usm_start_index_v_' // dum, surf_usm_v(l)%start_index )
    7579           CALL wrd_mpi_io( 'usm_end_index_v_' // dum, surf_usm_v(l)%end_index )
     7574                                            surface_data_to_write, global_start_index,             &
     7575                                            global_end_index )
     7576
    75807577          CALL wrd_mpi_io( 'usm_global_start_v_' // dum, global_start_index )
     7578          CALL wrd_mpi_io( 'usm_global_end_v_' // dum,   global_end_index )
    75817579
    75827580          IF ( .NOT. surface_data_to_write )  CYCLE
     
    75967594          WRITE( dum, '(I1)')  l
    75977595
    7598           CALL rd_mpi_io_surface_filetypes( surf_usm_h(l)%start_index, surf_usm_h(l)%end_index,             &
    7599                                             surface_data_to_write, global_start_index )
    7600 
    7601           CALL wrd_mpi_io( 'usm_start_index_h_2_' // dum,  surf_usm_h(l)%start_index )
    7602           CALL wrd_mpi_io( 'usm_end_index_h_2_' // dum, surf_usm_h(l)%end_index )
     7596          CALL rd_mpi_io_surface_filetypes( surf_usm_h(l)%start_index, surf_usm_h(l)%end_index,    &
     7597                                            surface_data_to_write, global_start_index,             &
     7598                                            global_end_index )
     7599
    76037600          CALL wrd_mpi_io( 'usm_global_start_h_2_' // dum, global_start_index )
     7601          CALL wrd_mpi_io( 'usm_global_end_h_2_' // dum, global_end_index )
    76047602
    76057603          IF ( .NOT. surface_data_to_write )  CYCLE
     
    76167614
    76177615          CALL rd_mpi_io_surface_filetypes( surf_usm_v(l)%start_index, surf_usm_v(l)%end_index,    &
    7618                                             surface_data_to_write, global_start_index )
    7619 
    7620           CALL wrd_mpi_io( 'usm_start_index_v_2_' //dum, surf_usm_v(l)%start_index )
    7621           CALL wrd_mpi_io( 'usm_end_index_v_2_' // dum, surf_usm_v(l)%end_index )
     7616                                            surface_data_to_write, global_start_index,             &
     7617                                            global_end_index )
     7618
    76227619          CALL wrd_mpi_io( 'usm_global_start_v_2_' // dum, global_start_index )
     7620          CALL wrd_mpi_io( 'usm_global_end_v_2_' // dum, global_end_index )
    76237621
    76247622          IF ( .NOT. surface_data_to_write )  CYCLE
  • palm/trunk/SOURCE/write_restart_data_mod.f90

    r4848 r4893  
    2424! -----------------
    2525! $Id$
     26! version number update because of revised output of surface data via MPI-IO for better performance
     27!
     28! 4848 2021-01-21 15:51:51Z gronemeier
    2629! bugfix: removed syn_turb_gen from restart files
    2730!
     
    201204    INTEGER ::  i                                !< loop index
    202205
    203     binary_version_global = '5.2'
     206    binary_version_global = '5.3'
    204207
    205208    IF ( restart_data_format_output == 'fortran_binary' )  THEN
Note: See TracChangeset for help on using the changeset viewer.