regular_grid_spawn.F90

Go to the documentation of this file.

!SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
!SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
!SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
!SFX_LIC for details. version 1.
!     ################################################################
      SUBROUTINE regular_grid_spawn(U,KLUOUT,                               &
                                      KL1, KIMAX1,KJMAX1,PX1,PY1,PDX1,PDY1, &
                                      KXOR, KYOR, KDXRATIO, KDYRATIO,       &
                                      KXSIZE, KYSIZE,                       &
                                      KL2, KIMAX_C_ll,KJMAX_C_ll,PX2,PY2,PDX2,PDY2  )
!     ################################################################
!
!!****  *REGULAR_GRID_SPAWN* - routine to read in namelist the horizontal grid
!!
!!    PURPOSE
!!    -------
!!
!!**  METHOD
!!    ------
!!
!!    EXTERNAL
!!    --------
!!
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!
!!    REFERENCE
!!    ---------
!!
!!
!!    AUTHOR
!!    ------
!!      V. Masson   *Meteo France*      
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    01/2004 
!!        M.Moge    04/2015  Parallelization using routines from MNH/SURCOUCHE
!!        M.Moge    06/2015  bug fix for reproductibility using UPDATE_NHALO1D
!!        M.Moge    01/2016  bug fix for parallel execution with SPLIT2
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE modd_surf_par, ONLY : nundef
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
USE modd_surf_atm_n, ONLY : surf_atm_t
!
USE modi_abor1_sfx
#ifdef MNH_PARALLEL
USE mode_ll
USE mode_modeln_handler

USE mode_splitting_ll, ONLY : split2, def_splitting2
USE modd_var_ll, ONLY : nproc, ip, ysplitting
USE modd_structure_ll, ONLY : zone_ll, crspd_ll
USE modd_parameters, ONLY : jphext
USE mode_tools_ll, ONLY : intersection
USE mode_exchange_ll, ONLY : send_recv_field
USE modi_update_nhalo1d
#endif
!
IMPLICIT NONE
!
!*       0.1   Declarations of arguments
!              -------------------------
!
TYPE(surf_atm_t), INTENT(INOUT) :: U
INTEGER,                      INTENT(IN)    :: KLUOUT     ! output listing logical unit
INTEGER,                      INTENT(IN)    :: KL1        ! total number of points KIMAX1 * KJMAX1
INTEGER,                      INTENT(IN)    :: KIMAX1     ! number of points in x direction
INTEGER,                      INTENT(IN)    :: KJMAX1     ! number of points in y direction
REAL, DIMENSION(KL1),         INTENT(IN)    :: PX1        ! X coordinate of all points
REAL, DIMENSION(KL1),         INTENT(IN)    :: PY1        ! Y coordinate of all points
REAL, DIMENSION(KL1),         INTENT(IN)    :: PDX1       ! X mesh size of all points
REAL, DIMENSION(KL1),         INTENT(IN)    :: PDY1       ! Y mesh size of all points
INTEGER,                      INTENT(IN)    :: KXOR       ! position of modified bottom left point
INTEGER,                      INTENT(IN)    :: KYOR       ! according to initial grid
INTEGER,                      INTENT(IN)    :: KXSIZE     ! number of grid meshes in initial grid to be
INTEGER,                      INTENT(IN)    :: KYSIZE     ! covered by the modified grid
INTEGER,                      INTENT(IN)    :: KDXRATIO   ! resolution ratio between modified grid
INTEGER,                      INTENT(IN)    :: KDYRATIO   ! and initial grid
INTEGER,                      INTENT(IN)    :: KL2        ! total number of points KIMAX_C_ll * KJMAX_C_ll
INTEGER,                      INTENT(INOUT)    :: KIMAX_C_ll     ! number of points in x direction (glb on entry, lcl on exit)
INTEGER,                      INTENT(INOUT)    :: KJMAX_C_ll     ! number of points in y direction (glb on entry, lcl on exit)
REAL, DIMENSION(:),ALLOCATABLE,         INTENT(OUT)   :: PX2        ! X coordinate of all points
REAL, DIMENSION(:),ALLOCATABLE,         INTENT(OUT)   :: PY2        ! Y coordinate of all points
REAL, DIMENSION(:),ALLOCATABLE,         INTENT(OUT)   :: PDX2       ! X mesh size of all points
REAL, DIMENSION(:),ALLOCATABLE,         INTENT(OUT)   :: PDY2       ! Y mesh size of all points
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
!* coarse/father grid
!
REAL, DIMENSION(:),     ALLOCATABLE :: ZXM1     ! X coordinate of center of mesh (IIMAX1   points)
REAL, DIMENSION(:),     ALLOCATABLE :: ZYM1     ! Y coordinate of center of mesh (IJMAX1   points)
REAL, DIMENSION(:),     ALLOCATABLE :: ZXHAT1   ! X coordinate of left side      (IIMAX1+1 points)
REAL, DIMENSION(:),     ALLOCATABLE :: ZYHAT1   ! Y coordinate of bottom side    (IJMAX1+1 points)
REAL, DIMENSION(:,:,:), ALLOCATABLE :: ZXHAT1_3D, ZYHAT1_3D   ! ZXHAT1 and ZXHAT1 copied in a 3D field for the communications
!
!* fine/son grid
!
REAL, DIMENSION(:),   ALLOCATABLE :: ZXHAT2   ! X coordinate of left side      (IIMAX2 points)
REAL, DIMENSION(:),   ALLOCATABLE :: ZYHAT2   ! Y coordinate of bottom side    (IJMAX2 points)
REAL, DIMENSION(:),   ALLOCATABLE :: ZXHAT2_F_TMP
REAL, DIMENSION(:),   ALLOCATABLE :: ZYHAT2_F_TMP
REAL, DIMENSION(:,:,:),   ALLOCATABLE :: ZXHAT2_F, ZYHAT2_F   ! temporary 3D fields to communicate the values on the father grid to the local son subgrid
!
!* other variables
!
INTEGER     :: JL            ! loop counter
INTEGER     :: JI,JJ         ! loop controls relatively to modified grid
#ifndef MNH_PARALLEL
INTEGER     :: JIBOX,JJBOX   ! grid mesh relatively to initial grid
REAL        :: ZCOEF         ! ponderation coefficient for linear interpolation
#endif
REAL(KIND=JPRB) :: ZHOOK_HANDLE
INTEGER :: IMI
INTEGER :: IINFO_ll
INTEGER :: IXDOMAINS, IYDOMAINS               ! number of subdomains in X and Y directions
LOGICAL :: GPREM                              ! needed for DEF_SPLITTING2, true if NPROC is a prime number
INTEGER :: IXOR_F_ll, IYOR_F_ll               ! origin of local father subdomain in global coordinates
INTEGER :: IXDIM_C, IYDIM_C                   ! size of local son subdomain (in coarse/father grid)
INTEGER :: IXOR_C_ll, IYOR_C_ll               ! origin of local son subdomain (in global fine/son grid)
INTEGER :: IXEND_C_ll, IYEND_C_ll             ! end of local son subdomain (in global fine/son grid)
INTEGER :: IXOR_C_COARSE_ll, IYOR_C_COARSE_ll ! origin of local son subdomain (in global coarse/father grid)
INTEGER :: IIMAX_C     ! number of points in x direction in local portion of son model (in fine grid)
INTEGER :: IJMAX_C     ! number of points in y direction in local portion of son model (in fine grid)
REAL, DIMENSION(KDXRATIO) :: ZCOEFX   ! ponderation coefficients for linear interpolation
REAL, DIMENSION(KDYRATIO) :: ZCOEFY   ! ponderation coefficients for linear interpolation
!
! structures for the partitionning
!
#ifdef MNH_PARALLEL
TYPE(zone_ll), DIMENSION(NPROC) :: TZSPLITTING_C  !splitting of child model
TYPE(zone_ll), ALLOCATABLE, DIMENSION(:) :: TZCOARSEFATHER ! Coarse father grid splitting
TYPE(zone_ll), ALLOCATABLE, DIMENSION(:) :: TZCOARSESONSPLIT   ! coarse son grid intersection with local father subdomain : coordinates in the father grid
!
! structures for the communications
!
TYPE(zone_ll), ALLOCATABLE, DIMENSION(:)  :: TZSEND, TZRECV
TYPE(crspd_ll), POINTER      :: TZCRSPDSEND, TZCRSPDRECV
TYPE(crspd_ll), ALLOCATABLE, DIMENSION(:), TARGET :: TZCRSPDSENDTAB, TZCRSPDRECVTAB
#endif
!
INTEGER :: J
INTEGER :: INBMSG
INTEGER :: ICARD
INTEGER :: ICARDDIF
!
!------------------------------------------------------------------------------
!
!*       1.    Coherence tests
!              ---------------
!
!* tests
!
IF (lhook) CALL dr_hook('REGULAR_GRID_SPAWN',0,zhook_handle)
#ifdef MNH_PARALLEL
IF ( kxor+kxsize-1 > u%NIMAX_SURF_ll ) THEN
  WRITE(kluout,*) 'spawned domain is not contained in the input domain'
  WRITE(kluout,*) 'IXOR = ', kxor, ' IXSIZE = ', kxsize,&
                    ' with NIMAX(file) = ', u%NIMAX_SURF_ll
  CALL abor1_sfx('REGULAR_GRID_SPAWN: (1) SPAWNED DOMAIN NOT CONTAINED IN INPUT DOMAIN')
END IF
IF ( kyor+kysize-1 > u%NJMAX_SURF_ll ) THEN
  WRITE(kluout,*) 'spawned domain is not contained in the input domain'
  WRITE(kluout,*) 'IYOR = ', kyor, ' IYSIZE = ', kysize,&
                    ' with NJMAX(file) = ',  u%NJMAX_SURF_ll
  CALL abor1_sfx('REGULAR_GRID_SPAWN: (2) SPAWNED DOMAIN NOT CONTAINED IN INPUT DOMAIN')
END IF
#else
IF ( kxor+kxsize-1 > kimax1 ) THEN
  WRITE(kluout,*) 'spawned domain is not contained in the input domain'
  WRITE(kluout,*) 'IXOR = ', kxor, ' IXSIZE = ', kxsize,&
                    ' with NIMAX(file) = ', kimax1  
  CALL abor1_sfx('REGULAR_GRID_SPAWN: (1) SPAWNED DOMAIN NOT CONTAINED IN INPUT DOMAIN')
END IF
IF ( kyor+kysize-1 > kjmax1 ) THEN
  WRITE(kluout,*) 'spawned domain is not contained in the input domain'
  WRITE(kluout,*) 'IYOR = ', kyor, ' IYSIZE = ', kysize,&
                    ' with NJMAX(file) = ', kjmax1  
  CALL abor1_sfx('REGULAR_GRID_SPAWN: (2) SPAWNED DOMAIN NOT CONTAINED IN INPUT DOMAIN')
END IF
#endif

!
!------------------------------------------------------------------------------
!
!*       2.    Partitionning of the son subdomain
!              --------------------------------------------------------------
!
#ifdef MNH_PARALLEL
! get origin of local father subdomain in global coordinates
!
 CALL get_or_ll( "B", ixor_f_ll, iyor_f_ll )
!
! origin of local son subdomain in global father coordinates
!
!IXOR_C_COARSE_ll = MAX( IXOR_F_ll, KXOR+1 )
!IYOR_C_COARSE_ll = MAX( IYOR_F_ll, KYOR+1 )
ixor_c_coarse_ll = max( ixor_f_ll-1, kxor ) ! we have to add one point on the west and south sides -> hence the "- 1"
iyor_c_coarse_ll = max( iyor_f_ll-1, kyor ) ! we have to add one point on the west and south sides -> hence the "- 1"
!
ALLOCATE(tzcoarsefather(nproc))
ALLOCATE(tzcoarsesonsplit(nproc))
!
! compute father partitioning
!
 CALL split2( u%NIMAX_SURF_ll,  u%NJMAX_SURF_ll, 1, nproc,tzcoarsefather, ysplitting)
! we don't want the halo
DO j = 1, nproc
  tzcoarsefather(j)%NXOR = tzcoarsefather(j)%NXOR - jphext
  tzcoarsefather(j)%NYOR = tzcoarsefather(j)%NYOR - jphext
  tzcoarsefather(j)%NXEND = tzcoarsefather(j)%NXEND - jphext
  tzcoarsefather(j)%NYEND = tzcoarsefather(j)%NYEND - jphext
ENDDO
!
! partition son domain on father grid (with global coordinates on father grid)
!
! we have to add one point on the west and south sides -> hence the "- 1"
! Warning : we cannot just call SPLIT2(KXSIZE, KYSIZE, 1, NPROC, TZCOARSESONSPLIT, YSPLITTING) as it would not 
!           necessarily split the son domain the same way the father domain was splitted
!           example : if father domain is 30x40 and son domain is 6x5 (in father grid dimensions) then
!                     with NPROC = 2, SPLIT2 will split father domain along Y dimension -> 30x20 local domains
!                     but SPLIT2 will split son domain along X dimension -> 3x5 local domains.
!           therefore we have to use DEF_SPLITTING2 and force the decomposition in the call to SPLIT2
! we want the same domain partitioning for the child domain and for the father domain
 CALL def_splitting2(ixdomains,iydomains, u%NIMAX_SURF_ll, u%NJMAX_SURF_ll,nproc,gprem)
 CALL split2(kxsize, kysize, 1, nproc, tzcoarsesonsplit, ysplitting, ixdomains, iydomains)

! compute the local size of son grid
! KIMAX_C_ll, KJMAX_C_ll are the global sizes of son domain
iimax_c = ( tzcoarsesonsplit(ip)%NXEND - tzcoarsesonsplit(ip)%NXOR + 1 ) * kdxratio
ijmax_c = ( tzcoarsesonsplit(ip)%NYEND - tzcoarsesonsplit(ip)%NYOR + 1 ) * kdyratio
! get the coordinates of the son domain partition on father grid
DO j = 1, nproc
  tzcoarsesonsplit(j)%NXOR = tzcoarsesonsplit(j)%NXOR + kxor - jphext - 1
  tzcoarsesonsplit(j)%NYOR = tzcoarsesonsplit(j)%NYOR + kyor - jphext - 1
  tzcoarsesonsplit(j)%NXEND = tzcoarsesonsplit(j)%NXEND + kxor - jphext
  tzcoarsesonsplit(j)%NYEND = tzcoarsesonsplit(j)%NYEND + kyor - jphext
ENDDO
!
! compute the local size of son grid
! KIMAX_C_ll, KJMAX_C_ll are the global sizes of son domain
!
!CALL SPLIT2 ( KIMAX_C_ll, KJMAX_C_ll, 1, NPROC, TZSPLITTING_C, YSPLITTING )
!IXOR_C_ll = TZSPLITTING_C(IP)%NXOR - JPHEXT
!IXEND_C_ll = TZSPLITTING_C(IP)%NXEND - JPHEXT
!IYOR_C_ll = TZSPLITTING_C(IP)%NYOR - JPHEXT
!IYEND_C_ll = TZSPLITTING_C(IP)%NYEND - JPHEXT
!!
!IIMAX_C = IXEND_C_ll - IXOR_C_ll + 1
!IJMAX_C = IYEND_C_ll - IYOR_C_ll + 1
!IIMAX_C = ( TZCOARSESONSPLIT(IP)%NXEND - TZCOARSESONSPLIT(IP)%NXOR + 1 ) * KDXRATIO
!IJMAX_C = ( TZCOARSESONSPLIT(IP)%NYEND - TZCOARSESONSPLIT(IP)%NYOR + 1 ) * KDYRATIO
!
!------------------------------------------------------------------------------
!
!*       3.    Preparing the structures for the communications for the initialization of son fields using father fields
!              --------------------------------------------------------------
!
  !
  ! ######## initializing the structures for the SEND ########
  !
  ALLOCATE(tzsend(nproc))
  CALL intersection( tzcoarsesonsplit, nproc, tzcoarsefather(ip), tzsend)
  ! il faut initialiser le TAG de manière a avoir un meme tag unique pour le send et le recv :
  !   on concatene le num du proc qui envoie et le num du proc qui recoit
  DO j = 1, nproc
    IF ( tzsend(j)%NUMBER > 0 ) THEN
      IF (tzsend(j)%NUMBER == 1) THEN
        tzsend(j)%MSSGTAG = ip * 10 + 1
      ELSE
        tzsend(j)%MSSGTAG = ip * 10**(ceiling(log10(real(tzsend(j)%number)))) + tzsend(j)%number
      ENDIF
    ENDIF
  ENDDO
  ! switching to local coordinates
  DO j = 1, nproc
    IF ( tzsend(j)%NUMBER > 0 ) THEN
       tzsend(j)%NXOR = tzsend(j)%NXOR - ixor_f_ll + 1
       tzsend(j)%NXEND = tzsend(j)%NXEND - ixor_f_ll + 1
       tzsend(j)%NYOR = tzsend(j)%NYOR - iyor_f_ll + 1
       tzsend(j)%NYEND = tzsend(j)%NYEND - iyor_f_ll + 1
    ENDIF
  ENDDO
  ! we do not need the Z dimension
  DO j = 1, nproc
    IF ( tzsend(j)%NUMBER > 0 ) THEN
       tzsend(j)%NZOR = 1
       tzsend(j)%NZEND = 1
    ENDIF
  ENDDO
  ! switching from an array of CRSPD_ll to a CRSPD_ll pointer
  inbmsg = 0
  DO j = 1, nproc
    IF ( tzsend(j)%NUMBER > 0 ) THEN
      inbmsg = inbmsg+1
    ENDIF
  ENDDO
  IF ( inbmsg > 0 ) THEN
    ALLOCATE( tzcrspdsendtab(inbmsg) )
    icard = 0
    icarddif = 0
    DO j = 1, nproc
      IF ( tzsend(j)%NUMBER > 0 ) THEN
        icard = icard+1
        IF ( tzsend(icard)%NUMBER /= ip ) THEN
          icarddif = icarddif+1
        ENDIF
        tzcrspdsendtab(icard)%TELT = tzsend(j)
        IF ( icard == inbmsg ) THEN
          tzcrspdsendtab(icard)%TNEXT => null()
        ELSE
          tzcrspdsendtab(icard)%TNEXT => tzcrspdsendtab(icard+1)
        ENDIF
      ENDIF
    ENDDO
    DO j = 1, icard
      tzcrspdsendtab(j)%NCARD = icard
      tzcrspdsendtab(j)%NCARDDIF = icarddif
    ENDDO
  ELSE 
    !il faut tout de meme mettre un element de taille 0 dans TZCRSPDSENDTAB
    !sinon SEND_RECV_FIELD plante en 02
    ALLOCATE( tzcrspdsendtab(1) )
    icard = 0
    icarddif = 0
    tzcrspdsendtab(1)%TELT = tzsend(1)
    tzcrspdsendtab(1)%TNEXT => null()
    tzcrspdsendtab(1)%NCARD = 0
    tzcrspdsendtab(1)%NCARDDIF = 0
  ENDIF
  IF (icard > 0) THEN
    tzcrspdsend => tzcrspdsendtab(1)
  ELSE
    tzcrspdsend => null()
  ENDIF
  !
  ! ######## initializing the structures for the RECV ########
  !
  ALLOCATE(tzrecv(nproc))
  CALL intersection( tzcoarsefather, nproc, tzcoarsesonsplit(ip), tzrecv )
  ! il faut initialiser le TAG de manière a avoir un meme tag unique pour le send et le recv :
  !   on concatene le num du proc qui envoie et le num du proc qui recoit
  DO j = 1, nproc
    IF ( tzrecv(j)%NUMBER > 0 ) THEN
      IF (ip == 1) THEN
        tzrecv(j)%MSSGTAG = tzrecv(j)%NUMBER * 10 + 1
      ELSE
        tzrecv(j)%MSSGTAG = tzrecv(j)%NUMBER * 10**(ceiling(log10(real(ip)))) + ip
      ENDIF
    ENDIF
  ENDDO
  ! switching to local coordinates
  DO j = 1, nproc
    IF ( tzrecv(j)%NUMBER > 0 ) THEN
       tzrecv(j)%NXOR = tzrecv(j)%NXOR - tzcoarsesonsplit(ip)%NXOR + 1
       tzrecv(j)%NXEND = tzrecv(j)%NXEND - tzcoarsesonsplit(ip)%NXOR + 1
       tzrecv(j)%NYOR = tzrecv(j)%NYOR - tzcoarsesonsplit(ip)%NYOR + 1
       tzrecv(j)%NYEND = tzrecv(j)%NYEND - tzcoarsesonsplit(ip)%NYOR + 1
    ENDIF
  ENDDO
  ! we do not need the Z dimension
  DO j = 1, nproc
    IF ( tzrecv(j)%NUMBER > 0 ) THEN
       tzrecv(j)%NZOR = 1
       tzrecv(j)%NZEND = 1
    ENDIF
  ENDDO
  ! switching from an array of CRSPD_ll to a CRSPD_ll pointer
  inbmsg = 0
  DO j = 1, nproc
    IF ( tzrecv(j)%NUMBER > 0 ) THEN
      inbmsg = inbmsg+1
    ENDIF
  ENDDO
  IF ( inbmsg > 0 ) THEN
    ALLOCATE( tzcrspdrecvtab(inbmsg) )
    icard = 0
    icarddif = 0
    DO j = 1, nproc
      IF ( tzrecv(j)%NUMBER > 0 ) THEN
        icard = icard+1
        IF ( tzrecv(icard)%NUMBER /= ip ) THEN
          icarddif = icarddif+1
        ENDIF
        tzcrspdrecvtab(icard)%TELT = tzrecv(j)
        IF ( icard == inbmsg ) THEN
          tzcrspdrecvtab(icard)%TNEXT => null()
        ELSE
          tzcrspdrecvtab(icard)%TNEXT => tzcrspdrecvtab(icard+1)
        ENDIF
      ENDIF
    ENDDO
    DO j = 1, icard
      tzcrspdrecvtab(j)%NCARD = icard
      tzcrspdrecvtab(j)%NCARDDIF = icarddif
    ENDDO
  ELSE
    !il faut tout de meme mettre un element de taille 0 dans TZCRSPDRECVTAB
    !sinon SEND_RECV_FIELD plante en 02
    ALLOCATE( tzcrspdrecvtab(1) )
    icard = 0
    icarddif = 0
    tzcrspdrecvtab(1)%TELT = tzsend(1)
    tzcrspdrecvtab(1)%TNEXT => null()
    tzcrspdrecvtab(1)%NCARD = 0
    tzcrspdrecvtab(1)%NCARDDIF = 0
  ENDIF
  IF (icard > 0) THEN
    tzcrspdrecv => tzcrspdrecvtab(1)
  ELSE
    tzcrspdrecv => null()
  ENDIF
#else
iimax_c = kimax_c_ll
ijmax_c = kjmax_c_ll
#endif
!
!------------------------------------------------------------------------------
!
!*       4.    Center of mesh coordinate arrays for each direction separately
!              --------------------------------------------------------------
!
! allocate the fields on the local son grid
!
#ifdef MNH_PARALLEL
ALLOCATE(px2(iimax_c*ijmax_c))
ALLOCATE(py2(iimax_c*ijmax_c))
ALLOCATE(pdx2(iimax_c*ijmax_c))
ALLOCATE(pdy2(iimax_c*ijmax_c))
#else
ALLOCATE(px2(kl2))
ALLOCATE(py2(kl2))
ALLOCATE(pdx2(kl2))
ALLOCATE(pdy2(kl2))
#endif
!
ALLOCATE(zxhat2(iimax_c+1))
ALLOCATE(zyhat2(ijmax_c+1))
!
! allocate the fields on the local father grid
!
ALLOCATE(zxm1(kimax1))
ALLOCATE(zym1(kjmax1))
ALLOCATE(zxhat1(kimax1+1))
ALLOCATE(zyhat1(kjmax1+1))
!
zxm1(:) = px1(1:kimax1)
DO jl=1,kl1
  IF (mod(jl,kimax1)==0) zym1(jl/kimax1) = py1(jl)
END DO
!
!------------------------------------------------------------------------------
!
!*       5.    side of mesh coordinate arrays for each direction separately
!              ------------------------------------------------------------
!
!
IF (kimax1==1) THEN
  zxhat1(1) = zxm1(1) - 0.5 * pdx1(1)
  zxhat1(2) = zxm1(1) + 0.5 * pdx1(1)
ELSE
  zxhat1(1) = 1.5 * zxm1(1) - 0.5 * zxm1(2)
  DO ji=2,kimax1
    zxhat1(ji) = 0.5 * zxm1(ji-1) + 0.5 * zxm1(ji)
  END DO
  zxhat1(kimax1+1) = 1.5 * zxm1(kimax1) - 0.5 * zxm1(kimax1-1)
END IF
!
IF (kjmax1==1) THEN
  zyhat1(1) = zym1(1) - 0.5 * pdy1(1)
  zyhat1(2) = zym1(1) + 0.5 * pdy1(1)
ELSE
  zyhat1(1) = 1.5 * zym1(1) - 0.5 * zym1(2)
  DO jj=2,kjmax1
    zyhat1(jj) = 0.5 * zym1(jj-1) + 0.5 * zym1(jj)
  END DO
  zyhat1(kjmax1+1) = 1.5 * zym1(kjmax1) - 0.5 * zym1(kjmax1-1)
END IF
#ifdef MNH_PARALLEL
  !
  ! do the communication
  !
  ixdim_c = tzcoarsesonsplit(ip)%NXEND-tzcoarsesonsplit(ip)%NXOR+1
  iydim_c = tzcoarsesonsplit(ip)%NYEND-tzcoarsesonsplit(ip)%NYOR+1
  ALLOCATE(zxhat2_f(ixdim_c,iydim_c,1))
  ALLOCATE(zyhat2_f(ixdim_c,iydim_c,1))
  ALLOCATE(zxhat1_3d(kimax1,kjmax1,1))
  ALLOCATE(zyhat1_3d(kimax1,kjmax1,1))
  zxhat1_3d(:,:,:) = 0
  zyhat1_3d(:,:,:) = 0
  zxhat2_f(:,:,:) = 0
  zyhat2_f(:,:,:) = 0
  DO j=1, kjmax1
    zxhat1_3d(:,j,1) = zxhat1(1:kimax1)
  ENDDO
  DO j=1, kimax1
    zyhat1_3d(j,:,1) = zyhat1(1:kjmax1)
  ENDDO
  CALL send_recv_field( tzcrspdsend, tzcrspdrecv, zxhat1_3d, zxhat2_f, iinfo_ll)
  CALL send_recv_field( tzcrspdsend, tzcrspdrecv, zyhat1_3d, zyhat2_f, iinfo_ll)
!
! We have to copy the entries of ZXHAT1_3D and ZYHAT1_3D that are local to the current process,
! and that are therefore not communicated in SEND_RECV_FIELD, in ZXHAT2_F and ZYHAT2_F
!
  IF ( tzsend(ip)%NUMBER /= 0 ) THEN  !if there are entries of ZXHAT1_3D and ZYHAT1_3D that are local to the current process
!  DO J=TZSEND(IP)%NXOR-KXOR,TZSEND(IP)%NXEND-KXOR
    zxhat2_f( tzrecv(ip)%NXOR:tzrecv(ip)%NXEND, 1, 1) = zxhat1_3d( tzsend(ip)%NXOR:tzsend(ip)%NXEND, 1, 1)
!  ENDDO
!  DO J=TZSEND(IP)%NYOR-KYOR,TZSEND(IP)%NYEND-KYOR
    zyhat2_f( 1,tzrecv(ip)%NYOR:tzrecv(ip)%NYEND, 1) = zyhat1_3d( 1,tzsend(ip)%NYOR:tzsend(ip)%NYEND, 1)
!  ENDDO
  ENDIF
  ! 
  ! We need one halo point on the east and north sides of each local subdomain to do a proper interpolation
  !
  ALLOCATE( zxhat2_f_tmp(ixdim_c+1) )
  ALLOCATE( zyhat2_f_tmp(iydim_c+1) )
  zxhat2_f_tmp(:) = 0.
  zyhat2_f_tmp(:) = 0.
  zxhat2_f_tmp(1:ixdim_c) = zxhat2_f(:,1,1)
  zyhat2_f_tmp(1:iydim_c) = zyhat2_f(1,:,1)
! we want the same domain partitioning for the child domain and for the father domain
  CALL def_splitting2(ixdomains,iydomains, u%NIMAX_SURF_ll, u%NJMAX_SURF_ll,nproc,gprem)
  CALL split2(kxsize, kysize, 1, nproc, tzcoarsesonsplit, ysplitting, ixdomains, iydomains)
  CALL update_nhalo1d( 1, zxhat2_f_tmp, kxsize, kysize,tzcoarsesonsplit(ip)%NXOR, &
    tzcoarsesonsplit(ip)%NXEND,tzcoarsesonsplit(ip)%NYOR,tzcoarsesonsplit(ip)%NYEND, 'XX    ')
  CALL update_nhalo1d( 1, zyhat2_f_tmp, kxsize, kysize,tzcoarsesonsplit(ip)%NXOR, &
    tzcoarsesonsplit(ip)%NXEND,tzcoarsesonsplit(ip)%NYOR,tzcoarsesonsplit(ip)%NYEND, 'YY    ')
#endif
!
!------------------------------------------------------------------------------
!
!*       6.    Interpolation of coordinate arrays for each direction separately
!              ----------------------------------------------------------------
!
!* X coordinate array
!
#ifndef MNH_PARALLEL
DO ji=1,iimax_c
  jibox=(ji-1)/kdxratio + kxor
  zcoef= float(mod(ji-1,kdxratio))/float(kdxratio)
  zxhat2(ji)=(1.-zcoef)*zxhat1(jibox)+zcoef*zxhat1(jibox+1)
END DO
IF (iimax_c==1) THEN
  zxhat2(iimax_c+1) = zxhat2(iimax_c) + zxhat1(jibox+1) - zxhat1(jibox)
ELSE
  zxhat2(iimax_c+1) = 2. * zxhat2(iimax_c) - zxhat2(iimax_c-1)
END IF
#else
DO j=0,kdxratio-1
  zcoefx(j+1) = float(j)/float(kdxratio)
ENDDO
DO ji=1,ixdim_c-1
  DO jj=1,kdxratio
    zxhat2((ji-1)*kdxratio+jj)=(1.-zcoefx(jj))*zxhat2_f(ji,1,1)+zcoefx(jj)*zxhat2_f(ji+1,1,1)
  ENDDO
ENDDO
IF (iimax_c==1) THEN
  zxhat2(iimax_c+1) = zxhat2(iimax_c) + zxhat2_f(ji,1,1) - zxhat2_f(ji,1,1)
ELSE
   IF ( least_ll() ) THEN ! the east halo point does not have a correct value so have to do an extrapolation
      zxhat2(iimax_c+1) = 2. * zxhat2(iimax_c) - zxhat2(iimax_c-1)
   ELSE
      zxhat2(iimax_c+1)=(1.-zcoefx(1))*zxhat2_f_tmp(ixdim_c)+zcoefx(1)*zxhat2_f_tmp(ixdim_c+1)
   ENDIF
END IF
#endif
!
!* Y coordinate array
!
#ifndef MNH_PARALLEL
DO jj=1,ijmax_c
  jjbox=(jj-1)/kdyratio + kyor
  zcoef= float(mod(jj-1,kdyratio))/float(kdyratio)
  zyhat2(jj)=(1.-zcoef)*zyhat1(jjbox)+zcoef*zyhat1(jjbox+1)
END DO
IF (ijmax_c==1) THEN
  zyhat2(ijmax_c+1) = zyhat2(ijmax_c) + zyhat1(jjbox+1) - zyhat1(jjbox)
ELSE
  zyhat2(ijmax_c+1) = 2. * zyhat2(ijmax_c) - zyhat2(ijmax_c-1)
END IF
#else
DO j=0,kdyratio-1
  zcoefy(j+1) = float(j)/float(kdyratio)
ENDDO
DO ji=1,iydim_c-1
  DO jj=1,kdyratio
    zyhat2((ji-1)*kdyratio+jj)=(1.-zcoefy(jj))*zyhat2_f(1,ji,1)+zcoefy(jj)*zyhat2_f(1,ji+1,1)
  ENDDO
ENDDO
IF (ijmax_c==1) THEN
  zyhat2(ijmax_c+1) = zyhat2(ijmax_c) + zyhat2_f(1,ji,1) - zyhat2_f(1,ji,1)
ELSE
  IF ( lnorth_ll() ) THEN ! the east halo point does not have a correct value so have to do an extrapolation
    zyhat2(ijmax_c+1) = 2. * zyhat2(ijmax_c) - zyhat2(ijmax_c-1)
  ELSE
    zyhat2(ijmax_c+1)=(1.-zcoefy(1))*zyhat2_f_tmp(iydim_c)+zcoefy(1)*zyhat2_f_tmp(iydim_c+1)
  ENDIF
END IF
#endif
!---------------------------------------------------------------------------
DEALLOCATE(zxm1)
DEALLOCATE(zym1)
DEALLOCATE(zxhat1)
DEALLOCATE(zyhat1)
#ifdef MNH_PARALLEL
DEALLOCATE(zxhat1_3d)
DEALLOCATE(zyhat1_3d)
#endif
!------------------------------------------------------------------------------
!
!*       7.    Coordinate arrays of all points
!              -------------------------------
!
DO jj=1,ijmax_c
  DO ji=1,iimax_c
    jl = (jj-1) * iimax_c + ji
      px2(jl) = 0.5 * zxhat2(ji) + 0.5 * zxhat2(ji+1)
      pdx2(jl) = zxhat2(ji+1) - zxhat2(ji)
      py2(jl) = 0.5 * zyhat2(jj) + 0.5 * zyhat2(jj+1)
      pdy2(jl) = zyhat2(jj+1) - zyhat2(jj)
  END DO
END DO
!
#ifdef MNH_PARALLEL
kimax_c_ll = iimax_c
kjmax_c_ll = ijmax_c
#endif
!---------------------------------------------------------------------------
#ifdef MNH_PARALLEL
DEALLOCATE(zxhat2_f)
DEALLOCATE(zyhat2_f)
DEALLOCATE(tzcrspdsendtab)
DEALLOCATE(tzcrspdrecvtab)
DEALLOCATE(tzsend)
DEALLOCATE(tzrecv)
DEALLOCATE(tzcoarsefather)
DEALLOCATE(tzcoarsesonsplit)
#endif
DEALLOCATE(zxhat2)
DEALLOCATE(zyhat2)
IF (lhook) CALL dr_hook('REGULAR_GRID_SPAWN',1,zhook_handle)
!---------------------------------------------------------------------------
!
END SUBROUTINE regular_grid_spawn