diff --git a/MNH/get_halo.f90 b/MNH/get_halo.f90
index 1acdd951c781d41c7b7bfa6563ef40f1d530972a..8cf520eaa5c34eb42164bd0424283e6560ed5566 100644
--- a/MNH/get_halo.f90
+++ b/MNH/get_halo.f90
@@ -9,29 +9,41 @@
 !     ####################
 !
 INTERFACE
+   SUBROUTINE GET_HALO2(PSRC,TP_PSRC_HALO2_ll)
+     !
+     USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
+     REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
+     TYPE(HALO2LIST_ll), POINTER         :: TP_PSRC_HALO2_ll          ! halo2 for SRC
+     !
+   END SUBROUTINE GET_HALO2
+END INTERFACE
 !
-SUBROUTINE GET_HALO2(PSRC,TP_PSRC_HALO2_ll)
-!
-USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-TYPE(HALO2LIST_ll), POINTER         :: TP_PSRC_HALO2_ll          ! halo2 for SRC
-!
-END SUBROUTINE GET_HALO2
-!
-SUBROUTINE GET_HALO(PSRC,HDIR)
-!
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-CHARACTER(len=4), OPTIONAL :: HDIR ! to send only halo on X or Y direction
-!
-END SUBROUTINE GET_HALO
+INTERFACE 
+   SUBROUTINE GET_HALO(PSRC,HDIR)
+     !
+     REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
+     CHARACTER(len=4), OPTIONAL :: HDIR ! to send only halo on X or Y direction
+     !
+   END SUBROUTINE GET_HALO
+END INTERFACE
+! 
+INTERFACE            
+   SUBROUTINE GET_HALO_D(PSRC,HDIR)
+     !
+     REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
+     !$acc reflected (PSRC)
+     CHARACTER(len=4), OPTIONAL :: HDIR ! to send only halo on X or Y direction
+     !
+   END SUBROUTINE GET_HALO_D
+END INTERFACE
 !
+INTERFACE
 SUBROUTINE DEL_HALO2_ll(TPHALO2LIST)
 !
 USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
 TYPE(HALO2LIST_ll), POINTER :: TPHALO2LIST ! list of HALO2_lls
 !
 END SUBROUTINE DEL_HALO2_ll
-!
 END INTERFACE
 !
 END MODULE MODI_GET_HALO         
@@ -105,6 +117,50 @@ CALL CLEANLIST_ll(TZ_PSRC_ll)
 !
 END SUBROUTINE GET_HALO
 !-----------------------------------------------------------------------
+!-------------------------------------------------------------------------------
+!     #########################
+      SUBROUTINE GET_HALO_D(PSRC,HDIR)
+!     #########################
+!
+USE MODE_ll
+USE MODD_ARGSLIST_ll, ONLY : LIST_ll
+USE MODD_PARAMETERS, ONLY : JPHEXT
+!
+IMPLICIT NONE
+!
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
+!$acc reflected (PSRC)
+CHARACTER(len=4), OPTIONAL :: HDIR ! to send only halo on X or Y direction
+!
+TYPE(LIST_ll)     , POINTER      :: TZ_PSRC_ll               ! halo
+INTEGER                          :: IERROR                 ! error return code 
+INTEGER:: IIB,IJB    ! Begining useful area in x,y,z directions
+INTEGER:: IIE,IJE    ! End useful area in x,y,z directions
+!
+!
+NULLIFY( TZ_PSRC_ll)
+CALL GET_INDICE_ll(IIB,IJB,IIE,IJE)
+!
+!$acc update host (PSRC)
+! acc update host (PSRC(    :     , IJB:IJB  , : ))
+! acc update host (PSRC(    :     , IJE:IJE  , : ))
+! acc update host (PSRC( IIB:IIB  , IJB:IJE  , : ))
+! acc update host (PSRC( IIE:IIE  , IJB:IJE  , : ))
+
+CALL ADD3DFIELD_ll(TZ_PSRC_ll,PSRC)
+CALL UPDATE_HALO_ll(TZ_PSRC_ll,IERROR, HDIR=HDIR )
+CALL CLEANLIST_ll(TZ_PSRC_ll)
+
+! acc update device (PSRC)
+!$acc update device (PSRC(      :       ,      : IJB-1 , : ))
+!$acc update device (PSRC(      :       , IJE+1:       , : ))
+!$acc update device (PSRC(      :IIB-1  ,      :       , : ))
+!$acc update device (PSRC( IIE+1:       ,      :       , : ))
+
+!
+END SUBROUTINE GET_HALO_D
+!-----------------------------------------------------------------------
+!
 !
 !     ####################################
       SUBROUTINE DEL_HALO2_ll(TPHALO2LIST)
diff --git a/MNH/ppm.f90 b/MNH/ppm.f90
index 688423f7ffa5100f02bbfcff986ef0b0fd319ae9..f293b5b1daf89daae878f23929b7fb147f0d731b 100644
--- a/MNH/ppm.f90
+++ b/MNH/ppm.f90
@@ -10,73 +10,73 @@
 !
 INTERFACE
 !
-FUNCTION PPM_01_X(HLBCX, KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+FUNCTION PPM_01_X(HLBCX, KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX  ! X direction LBC type
 !
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_01_X
 !
-FUNCTION PPM_01_Y(HLBCY, KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+FUNCTION PPM_01_Y(HLBCY, KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY  ! Y direction LBC type
 !
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_01_Y
 !
-FUNCTION PPM_01_Z(KGRID, PSRC, ZCR, PRHO, PTSTEP) RESULT(PR)
+FUNCTION PPM_01_Z(KGRID, PSRC, PCR, PRHO, PTSTEP) RESULT(PR)
 !
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_01_Z
 !
-FUNCTION PPM_S0_X(HLBCX, KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+FUNCTION PPM_S0_X(HLBCX, KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX  ! X direction LBC type
 !
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_S0_X
 !
-FUNCTION PPM_S0_Y(HLBCY, KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+FUNCTION PPM_S0_Y(HLBCY, KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 !
 CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY ! Y direction LBC type
@@ -84,33 +84,33 @@ CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY ! Y direction LBC type
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_S0_Y
 !
-FUNCTION PPM_S0_Z(KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+FUNCTION PPM_S0_Z(KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 !
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_S0_Z
 !
-FUNCTION PPM_S1_X(HLBCX, KGRID, PSRC, ZCR, PRHO, PRHOT, &
+FUNCTION PPM_S1_X(HLBCX, KGRID, PSRC, PCR, PRHO, PRHOT, &
                         PTSTEP) RESULT(PR)
 !
 CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX  ! X direction LBC type
@@ -118,18 +118,18 @@ CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX  ! X direction LBC type
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHOT ! density at t+dt
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_S1_X
 !
-FUNCTION PPM_S1_Y(HLBCY, KGRID, PSRC, ZCR, PRHO, PRHOT, &
+FUNCTION PPM_S1_Y(HLBCY, KGRID, PSRC, PCR, PRHO, PRHOT, &
                         PTSTEP) RESULT(PR)
 !
 CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY  ! X direction LBC type
@@ -137,31 +137,31 @@ CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY  ! X direction LBC type
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHOT ! density at t+dt
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_S1_Y
 !
-FUNCTION PPM_S1_Z(KGRID, PSRC, ZCR, PRHO, PRHOT, PTSTEP) &
+FUNCTION PPM_S1_Z(KGRID, PSRC, PCR, PRHO, PRHOT, PTSTEP) &
                         RESULT(PR)
 !
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHOT ! density at t+dt
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 END FUNCTION PPM_S1_Z
 !
@@ -173,7 +173,7 @@ END MODULE MODI_PPM
 !-------------------------------------------------------------------------------
 !-------------------------------------------------------------------------------
 !     ########################################################################
-      FUNCTION PPM_01_X(HLBCX, KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+      FUNCTION PPM_01_X(HLBCX, KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 !     ########################################################################
 !!
@@ -193,7 +193,6 @@ USE MODI_SHUMAN
 USE MODI_GET_HALO
 !
 USE MODD_CONF
-!USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
 !BEG JUAN PPM_LL
 USE MODD_LUNIT
 !END JUAN PPM_LL
@@ -209,13 +208,13 @@ CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX  ! X direction LBC type
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 !*       0.2   Declarations of local variables :
 !
@@ -223,15 +222,15 @@ INTEGER:: IIB,IJB    ! Begining useful area in x,y,z directions
 INTEGER:: IIE,IJE    ! End useful area in x,y,z directions
 !
 ! terms used in parabolic interpolation, dmq, qL, qR, dq, q6
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZQL,ZQR
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZDQ,ZQ6
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZDMQ
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZQL,ZQR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZDQ,ZQ6
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZDMQ
 !
 ! extra variables for the initial guess of parabolae parameters
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZQL0,ZQR0,ZQ60
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZQL0,ZQR0,ZQ60
 !
 ! advection fluxes
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFPOS, ZFNEG
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZFPOS, ZFNEG
 !
 !BEG JUAN PPM_LL
 INTEGER                          :: ILUOUT,IRESP             ! for prints
@@ -266,7 +265,7 @@ GEAST = LEAST_ll()
 !
 !BEG JUAN PPM_LL
 !
-!*              initialise & update halo & halo2 for PSRC
+!*              initialise & update halo for PSRC
 !
 IF(NHALO /= 1) THEN
     CALL FMLOOK_ll(CLUOUT0,CLUOUT0,ILUOUT,IRESP)
@@ -279,9 +278,12 @@ IF(NHALO /= 1) THEN
 ENDIF
 !
 CALL GET_HALO(PSRC,HDIR="01_X")
-
-!$acc data region local (ZQL,ZQR,ZDQ,ZQ6,ZDMQ,ZQL0,ZQR0,ZQ60,ZFPOS,ZFNEG) copyin (psrc,zcr,prho) copyout(pr) !
+!
+#define JUAN_ACC_01_X
+#ifdef JUAN_ACC_01_X
+!$acc data region local (ZQL,ZQR,ZDQ,ZQ6,ZDMQ,ZQL0,ZQR0,ZQ60,ZFPOS,ZFNEG) copyin (psrc,pcr,prho) copyout(pr) 
 !$acc region 
+#endif
  DO K=1,IKU ; DO J = 1,IJU ; DO I=1,IIU
 PR(I,J,K)=PSRC(I,J,K)
 ZQL(I,J,K)=PSRC(I,J,K)
@@ -384,10 +386,10 @@ ENDDO ; ENDDO ; ENDDO
 !!$!
 !!$! and finally calculate fluxes for the advection
 !!$!
-!!$!  ZFPOS(i) = Fct[ ZQR(i-1),ZCR(i),ZDQ(i-1),ZQ6(i-1) ]
+!!$!  ZFPOS(i) = Fct[ ZQR(i-1),PCR(i),ZDQ(i-1),ZQ6(i-1) ]
 !!$!
-!!$   ZFPOS(IIB:IIE+1,IJS:IJN,:) = ZQR(IIB-1:IIE,IJS:IJN,:) - 0.5*ZCR(IIB:IIE+1,IJS:IJN,:) * &            
-!!$        (ZDQ(IIB-1:IIE,IJS:IJN,:) - (1.0 - 2.0*ZCR(IIB:IIE+1,IJS:IJN,:)/3.0)        &
+!!$   ZFPOS(IIB:IIE+1,IJS:IJN,:) = ZQR(IIB-1:IIE,IJS:IJN,:) - 0.5*PCR(IIB:IIE+1,IJS:IJN,:) * &            
+!!$        (ZDQ(IIB-1:IIE,IJS:IJN,:) - (1.0 - 2.0*PCR(IIB:IIE+1,IJS:IJN,:)/3.0)        &
 !!$        * ZQ6(IIB-1:IIE,IJS:IJN,:))
 !!$!
 !!$   CALL GET_HALO(ZFPOS,HDIR="Z1_X")
@@ -397,15 +399,15 @@ ENDDO ; ENDDO ; ENDDO
 !!$! PPOSX(IIB-1,:,:) is not important for the calc of advection so 
 !!$! we set it to 0
 !!$!
-!!$   ZFNEG(:,IJS:IJN,:) = ZQL(:,IJS:IJN,:) - 0.5*ZCR(:,IJS:IJN,:) *      &            
-!!$        ( ZDQ(:,IJS:IJN,:) + (1.0 + 2.0*ZCR(:,IJS:IJN,:)/3.0) * ZQ6(:,IJS:IJN,:) )
+!!$   ZFNEG(:,IJS:IJN,:) = ZQL(:,IJS:IJN,:) - 0.5*PCR(:,IJS:IJN,:) *      &            
+!!$        ( ZDQ(:,IJS:IJN,:) + (1.0 + 2.0*PCR(:,IJS:IJN,:)/3.0) * ZQ6(:,IJS:IJN,:) )
 !!$!
 !!$   CALL GET_HALO(ZFNEG,HDIR="Z1_X")
 !!$!
 !!$! advect the actual field in X direction by U*dt
 !!$!
-!!$   PR = DXF( ZCR*MXM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                             ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+!!$   PR = DXF( PCR*MXM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                             ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 !!$   CALL GET_HALO(PR,HDIR="01_X")   
 !!$   CALL MPPDB_CHECK3DM("PPM::PPM_01_X CYCL ::PR",PRECISION,PR)   
 !!$!
@@ -460,8 +462,13 @@ ENDDO ; ENDDO ; ENDDO
 !
 !  update ZDMQ HALO before next/further  utilisation 
 !
-!TEMPO_JUAN   CALL  GET_HALO(ZDMQ,HDIR="01_X")
-!TEMPO_JUAN   CALL MPPDB_CHECK3DM("PPM::PPM_01_X OPEN ::ZDMQ",PRECISION,ZDMQ)
+#define TEMPO_JUAN
+#ifdef TEMPO_JUAN  
+!$acc end region   
+CALL  GET_HALO_D(ZDMQ,HDIR="01_X") 
+CALL MPPDB_CHECK3DM("PPM::PPM_01_X OPEN ::ZDMQ",PRECISION,ZDMQ)
+!$acc region
+#endif   
 !
 ! calculate qL and qR
 !
@@ -470,7 +477,13 @@ ENDDO ; ENDDO ; ENDDO
    ZQL0(IIB:IIE+1,IJS:IJN,:) = 0.5*(PSRC(IIB:IIE+1,IJS:IJN,:) + PSRC(IIB-1:IIE,IJS:IJN,:)) - &
         (ZDMQ(IIB:IIE+1,IJS:IJN,:) - ZDMQ(IIB-1:IIE,IJS:IJN,:))/6.0
 !
-!TEMPO_JUAN  CALL  GET_HALO(ZQL0,HDIR="01_X")
+
+!
+#ifdef TEMPO_JUAN  
+!$acc end region  
+CALL  GET_HALO_D(ZQL0,HDIR="01_X") 
+!$acc region
+#endif
 !  
 !  WEST BOUND
 !
@@ -483,7 +496,11 @@ ENDDO ; ENDDO ; ENDDO
 !
    ZQR0(IIB-1:IIE,IJS:IJN,:) = ZQL0(IIB:IIE+1,IJS:IJN,:)
 !
-!TEMPO_JUAN   CALL  GET_HALO(ZQR0,HDIR="01_X")
+#ifdef TEMPO_JUAN   
+!$acc end region  
+CALL  GET_HALO_D(ZQR0,HDIR="01_X")
+!$acc region
+#endif
 !
 !  EAST BOUND
 !
@@ -530,13 +547,17 @@ ENDDO ; ENDDO ; ENDDO
 ! and finally calculate fluxes for the advection
 !
 !
-!  ZFPOS(i) = Fct[ ZQR(i-1),ZCR(i),ZDQ(i-1),ZQ6(i-1) ]
+!  ZFPOS(i) = Fct[ ZQR(i-1),PCR(i),ZDQ(i-1),ZQ6(i-1) ]
 !
-   ZFPOS(IIB:IIE+1,IJS:IJN,:) = ZQR(IIB-1:IIE,IJS:IJN,:) - 0.5*ZCR(IIB:IIE+1,IJS:IJN,:) * &            
-        (ZDQ(IIB-1:IIE,IJS:IJN,:) - (1.0 - 2.0*ZCR(IIB:IIE+1,IJS:IJN,:)/3.0)        &
+   ZFPOS(IIB:IIE+1,IJS:IJN,:) = ZQR(IIB-1:IIE,IJS:IJN,:) - 0.5*PCR(IIB:IIE+1,IJS:IJN,:) * &            
+        (ZDQ(IIB-1:IIE,IJS:IJN,:) - (1.0 - 2.0*PCR(IIB:IIE+1,IJS:IJN,:)/3.0)        &
         * ZQ6(IIB-1:IIE,IJS:IJN,:))
 !
-!TEMPO_JUAN   CALL GET_HALO(ZFPOSS,HDIR="01_X")
+#ifdef TEMPO_JUAN 
+!$acc end region     
+CALL GET_HALO_D(ZFPOS,HDIR="01_X")
+!$acc region 
+#endif
 !
 !
 !  WEST BOUND
@@ -545,24 +566,28 @@ ENDDO ; ENDDO ; ENDDO
 ! 
   IF (GWEST) THEN 
   ! acc region
-   ZFPOS(IIB,IJS:IJN,:) = (PSRC(IIB-1,IJS:IJN,:) - ZQR(IIB-1,IJS:IJN,:))*ZCR(IIB,IJS:IJN,:) + &
+   ZFPOS(IIB,IJS:IJN,:) = (PSRC(IIB-1,IJS:IJN,:) - ZQR(IIB-1,IJS:IJN,:))*PCR(IIB,IJS:IJN,:) + &
                     ZQR(IIB-1,IJS:IJN,:)
   ! acc end region
 ! PPOSX(IIB-1,:,:) is not important for the calc of advection so 
 ! we set it to 0
    ENDIF
 !
-   ZFNEG(:,IJS:IJN,:) = ZQL(:,IJS:IJN,:) - 0.5*ZCR(:,IJS:IJN,:) *      &            
-        ( ZDQ(:,IJS:IJN,:) + (1.0 + 2.0*ZCR(:,IJS:IJN,:)/3.0) * ZQ6(:,IJS:IJN,:) )
+   ZFNEG(:,IJS:IJN,:) = ZQL(:,IJS:IJN,:) - 0.5*PCR(:,IJS:IJN,:) *      &            
+        ( ZDQ(:,IJS:IJN,:) + (1.0 + 2.0*PCR(:,IJS:IJN,:)/3.0) * ZQ6(:,IJS:IJN,:) )
 !
-!TEMPO_JUAN   CALL GET_HALO(ZFNEG,HDIR="01_X")
+#ifdef TEMPO_JUAN  
+!$acc end region   
+CALL GET_HALO_D(ZFNEG,HDIR="01_X")
+!$acc region
+#endif
 !
 !  EAST BOUND
 !
 ! advection flux at open boundary when u(IIE+1) < 0
   IF (GEAST) THEN
   ! acc region
-   ZFNEG(IIE+1,IJS:IJN,:) = (ZQR(IIE,IJS:IJN,:)-PSRC(IIE+1,IJS:IJN,:))*ZCR(IIE+1,IJS:IJN,:) + &
+   ZFNEG(IIE+1,IJS:IJN,:) = (ZQR(IIE,IJS:IJN,:)-PSRC(IIE+1,IJS:IJN,:))*PCR(IIE+1,IJS:IJN,:) + &
                       ZQR(IIE,IJS:IJN,:)
   ! acc end region
   ENDIF
@@ -571,15 +596,16 @@ ENDDO ; ENDDO ; ENDDO
 !
     ! acc region   
     mxm(ZQL,PRHO)
-    ZQR =  ZCR* ZQL*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + ZFNEG*(0.5-SIGN(0.5,ZCR)) ) 
+    ZQR =  PCR* ZQL*( ZFPOS*(0.5+SIGN(0.5,PCR)) + ZFNEG*(0.5-SIGN(0.5,PCR)) ) 
     dxf(PR,ZQR)
 
-!!$   PR = DXF( ZCR*MXM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                             ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+!!$   PR = DXF( PCR*MXM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                             ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 
-    !$acc end region 
-! acc update host (pr)
+#ifdef JUAN_ACC_01_X
+!$acc end region 
 !$acc end data region 
+#endif
 
    CALL GET_HALO(PR,HDIR="01_X")   
 !
@@ -649,7 +675,7 @@ END FUNCTION PPM_01_X
 !-------------------------------------------------------------------------------
 !
 !     ########################################################################
-      FUNCTION PPM_01_Y(HLBCY, KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+      FUNCTION PPM_01_Y(HLBCY, KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 !     ########################################################################
 !!
@@ -669,7 +695,6 @@ USE MODI_SHUMAN
 USE MODI_GET_HALO
 !
 USE MODD_CONF
-!USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
 !BEG JUAN PPM_LL
 USE MODD_LUNIT
 !END JUAN PPM_LL
@@ -683,31 +708,38 @@ IMPLICIT NONE
 CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY  ! Y direction LBC type
 !
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
+REAL,                   INTENT(IN)  :: PTSTEP   ! Time step 
 !
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
-!
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
-REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC  &  ! variable at t
+                                     , PCR   &  ! Courant number
+                                     , PRHO     ! density
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR &
 !
 !*       0.2   Declarations of local variables :
 !
-INTEGER:: IIB,IJB    ! Begining useful area in x,y,z directions
-INTEGER:: IIE,IJE    ! End useful area in x,y,z directions
-!
 ! terms used in parabolic interpolation, dmq, qL, qR, dq, q6
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZQL,ZQR
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZDQ,ZQ6
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZDMQ
-!
+    , ZQL,ZQR , ZDQ,ZQ6 , ZDMQ & 
 ! extra variables for the initial guess of parabolae parameters
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZQL0,ZQR0,ZQ60
-!
+    , ZQL0,ZQR0,ZQ60 &
 ! advection fluxes
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFPOS, ZFNEG
+    , ZFPOS, ZFNEG
+
+
+INTEGER:: IIB,IJB    ! Begining useful area in x,y,z directions
+INTEGER:: IIE,IJE    ! End useful area in x,y,z directions
+!
+
+!!$REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZQL,ZQR
+!!$REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZDQ,ZQ6
+!!$REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZDMQ 
+!!$!
+!!$
+!!$REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZQL0,ZQR0,ZQ60
+!!$!
+!!$
+!!$REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZFPOS, ZFNEG
 !
 !BEG JUAN PPM_LL
 INTEGER                          :: ILUOUT,IRESP             ! for prints
@@ -750,14 +782,18 @@ IF(NHALO /= 1) THEN
     CALL ABORT
     STOP
 ENDIF
+!
 CALL GET_HALO(PSRC,HDIR="01_Y")
-
+!
 !
 !-------------------------------------------------------------------------------
 !
 !
-!$acc data region local (ZQL,ZQR,ZDQ,ZQ6,ZDMQ,ZQL0,ZQR0,ZQ60,ZFPOS,ZFNEG) copyin (psrc,zcr,prho) copyout(pr) !
-!$acc region ! local (ZQL,ZQR,ZDQ,ZQ6,ZDMQ,ZQL0,ZQR0,ZQ60,ZFPOS,ZFNEG) copyin (psrc,zcr,prho) copyout(pr)
+#define JUAN_ACC_01_Y
+#ifdef JUAN_ACC_01_Y
+!$acc data region local (ZQL,ZQR,ZDQ,ZQ6,ZDMQ,ZQL0,ZQR0,ZQ60,ZFPOS,ZFNEG) copyin (psrc,pcr,prho) copyout(pr)
+!$acc region
+#endif
 PR=PSRC
 ZQL=PSRC
 ZQR=PSRC
@@ -844,10 +880,10 @@ ZFNEG=PSRC
 !!$!
 !!$! and finally calculate fluxes for the advection
 !!$!
-!!$!  ZFPOS(j) = Fct[ ZQR(j-1),ZCR(i),ZDQ(j-1),ZQ6(j-1) ]
+!!$!  ZFPOS(j) = Fct[ ZQR(j-1),PCR(i),ZDQ(j-1),ZQ6(j-1) ]
 !!$!
-!!$   ZFPOS(IIW:IIA,IJB:IJE+1,:) = ZQR(IIW:IIA,IJB-1:IJE,:) - 0.5*ZCR(IIW:IIA,IJB:IJE+1,:) * &            
-!!$        (ZDQ(IIW:IIA,IJB-1:IJE,:) - (1.0 - 2.0*ZCR(IIW:IIA,IJB:IJE+1,:)/3.0)        &
+!!$   ZFPOS(IIW:IIA,IJB:IJE+1,:) = ZQR(IIW:IIA,IJB-1:IJE,:) - 0.5*PCR(IIW:IIA,IJB:IJE+1,:) * &            
+!!$        (ZDQ(IIW:IIA,IJB-1:IJE,:) - (1.0 - 2.0*PCR(IIW:IIA,IJB:IJE+1,:)/3.0)        &
 !!$        * ZQ6(IIW:IIA,IJB-1:IJE,:))
 !!$!
 !!$   CALL GET_HALO(ZFPOS,HDIR="01_Y")
@@ -857,15 +893,15 @@ ZFNEG=PSRC
 !!$! PPOSX(:,IJB-1,:) is not important for the calc of advection so 
 !!$! we set it to 0
 !!$!
-!!$   ZFNEG(IIW:IIA,:,:) = ZQL(IIW:IIA,:,:) - 0.5*ZCR(IIW:IIA,:,:) *      &            
-!!$        ( ZDQ(IIW:IIA,:,:) + (1.0 + 2.0*ZCR(IIW:IIA,:,:)/3.0) * ZQ6(IIW:IIA,:,:) )
+!!$   ZFNEG(IIW:IIA,:,:) = ZQL(IIW:IIA,:,:) - 0.5*PCR(IIW:IIA,:,:) *      &            
+!!$        ( ZDQ(IIW:IIA,:,:) + (1.0 + 2.0*PCR(IIW:IIA,:,:)/3.0) * ZQ6(IIW:IIA,:,:) )
 !!$!
 !!$   CALL GET_HALO(ZFNEG,HDIR="01_Y")
 !!$!
 !!$! advect the actual field in Y direction by V*dt
 !!$!
-!!$   PR = DYF( ZCR*MYM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                             ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+!!$   PR = DYF( PCR*MYM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                             ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 !!$   CALL GET_HALO(PR,HDIR="01_Y") 
 !!$   CALL MPPDB_CHECK3DM("PPM::PPM_01_Y CYCL ::PR",PRECISION,PR) 
 
@@ -925,7 +961,11 @@ ZFNEG=PSRC
 !
 !  update ZDMQ HALO before next/further  utilisation 
 !
-   !TEMPO_JUAN   CALL  GET_HALO(ZDMQ,HDIR="01_Y")  
+#ifdef TEMPO_JUAN   
+!$acc end region   
+CALL  GET_HALO_D(ZDMQ,HDIR="01_Y")  
+!$acc region   
+#endif
 !
 ! calculate qL and qR with the modified dmq
 !
@@ -941,7 +981,11 @@ ZFNEG=PSRC
 !!$   ENDDO ; ENDDO ; ENDDO
 !!$   !#acc end region
 !
-   !TEMPO_JUAN   CALL  GET_HALO(ZQL0,HDIR="01_Y")
+#ifdef TEMPO_JUAN  
+!$acc end region    
+CALL  GET_HALO_D(ZQL0,HDIR="01_Y")
+!$acc region   
+#endif
 !  
 !  SOUTH BOUND
 !
@@ -1013,21 +1057,25 @@ ZFNEG=PSRC
 !
    ! and finally calculate fluxes for the advection
    !#acc region
-   ZFPOS(IIW:IIA,IJB:IJE+1,:) = ZQR(IIW:IIA,IJB-1:IJE,:) - 0.5*ZCR(IIW:IIA,IJB:IJE+1,:) * &            
-        (ZDQ(IIW:IIA,IJB-1:IJE,:) - (1.0 - 2.0*ZCR(IIW:IIA,IJB:IJE+1,:)/3.0)        &
+   ZFPOS(IIW:IIA,IJB:IJE+1,:) = ZQR(IIW:IIA,IJB-1:IJE,:) - 0.5*PCR(IIW:IIA,IJB:IJE+1,:) * &            
+        (ZDQ(IIW:IIA,IJB-1:IJE,:) - (1.0 - 2.0*PCR(IIW:IIA,IJB:IJE+1,:)/3.0)        &
         * ZQ6(IIW:IIA,IJB-1:IJE,:))
    !#acc end region 
 
 !!$   !#acc region
 !!$   DO K=1,IKU ; DO J=IJB,IJE+1
 !!$         !#acc do independent
-!!$         DO I=IIW,IIA ; ZFPOS(I,J,K) = ZQR(I,J-1,K) - 0.5*ZCR(I,J,K) * &            
-!!$                                       (ZDQ(I,J-1,K) - (1.0 - 2.0*ZCR(I,J,K)/3.0)        &
+!!$         DO I=IIW,IIA ; ZFPOS(I,J,K) = ZQR(I,J-1,K) - 0.5*PCR(I,J,K) * &            
+!!$                                       (ZDQ(I,J-1,K) - (1.0 - 2.0*PCR(I,J,K)/3.0)        &
 !!$                                       * ZQ6(I,J-1,K))
 !!$   ENDDO ; ENDDO ; ENDDO
 !!$   !#acc end region 
 !
-   !TEMPO_JUAN   CALL GET_HALO(ZFPOS,HDIR="01_Y")
+#ifdef TEMPO_JUAN
+!$acc end region   
+CALL GET_HALO_D(ZFPOS,HDIR="01_Y")
+!$acc region   
+#endif
 !
 !
 ! advection flux at open boundary when u(IJB) > 0
@@ -1037,7 +1085,7 @@ ZFNEG=PSRC
 !!$   IF (LSOUTH_ll()) THEN
    IF (GSOUTH) THEN
    !#acc region
-    ZFPOS(IIW:IIA,IJB,:) = (PSRC(IIW:IIA,IJB-1,:) - ZQR(IIW:IIA,IJB-1,:))*ZCR(IIW:IIA,IJB,:) + &
+    ZFPOS(IIW:IIA,IJB,:) = (PSRC(IIW:IIA,IJB-1,:) - ZQR(IIW:IIA,IJB-1,:))*PCR(IIW:IIA,IJB,:) + &
                       ZQR(IIW:IIA,IJB-1,:)
    !#acc end region 
    ENDIF
@@ -1045,19 +1093,23 @@ ZFNEG=PSRC
 ! PPOSX(:,IJB-1,:) is not important for the calc of advection so 
 ! we set it to 0
    !#acc region
-   ZFNEG(IIW:IIA,:,:) = ZQL(IIW:IIA,:,:) - 0.5*ZCR(IIW:IIA,:,:) *      &            
-        ( ZDQ(IIW:IIA,:,:) + (1.0 + 2.0*ZCR(IIW:IIA,:,:)/3.0) * ZQ6(IIW:IIA,:,:) )
+   ZFNEG(IIW:IIA,:,:) = ZQL(IIW:IIA,:,:) - 0.5*PCR(IIW:IIA,:,:) *      &            
+        ( ZDQ(IIW:IIA,:,:) + (1.0 + 2.0*PCR(IIW:IIA,:,:)/3.0) * ZQ6(IIW:IIA,:,:) )
    !#acc end region
 
 !!$   !#acc region
 !!$   DO K=1,IKU ; DO J=1,IJU
 !!$         !#acc do independent
-!!$         DO I=IIW,IIA ; ZFNEG(I,J,K) = ZQL(I,J,K) - 0.5*ZCR(I,J,K) *      &            
-!!$        ( ZDQ(I,J,K) + (1.0 + 2.0*ZCR(I,J,K)/3.0) * ZQ6(I,J,K) )
+!!$         DO I=IIW,IIA ; ZFNEG(I,J,K) = ZQL(I,J,K) - 0.5*PCR(I,J,K) *      &            
+!!$        ( ZDQ(I,J,K) + (1.0 + 2.0*PCR(I,J,K)/3.0) * ZQ6(I,J,K) )
 !!$   ENDDO ; ENDDO ; ENDDO
 !!$   !#acc end region 
 !
-   !TEMPO_JUAN   CALL GET_HALO(ZFNEG,HDIR="01_Y")
+#ifdef TEMPO_JUAN
+!$acc end region   
+   CALL GET_HALO_D(ZFNEG,HDIR="01_Y")
+!$acc region   
+#endif
 !
 ! advection flux at open boundary when u(IJE+1) < 0
 !
@@ -1065,7 +1117,7 @@ ZFNEG=PSRC
 !
    IF (GNORTH) THEN 
    !#acc region     
-    ZFNEG(IIW:IIA,IJE+1,:) = (ZQR(IIW:IIA,IJE,:)-PSRC(IIW:IIA,IJE+1,:))*ZCR(IIW:IIA,IJE+1,:) + &
+    ZFNEG(IIW:IIA,IJE+1,:) = (ZQR(IIW:IIA,IJE,:)-PSRC(IIW:IIA,IJE+1,:))*PCR(IIW:IIA,IJE+1,:) + &
                         ZQR(IIW:IIA,IJE,:)
    !#acc end region 
    ENDIF
@@ -1074,19 +1126,20 @@ ZFNEG=PSRC
 !!$!
     !#acc region   
     mym(ZQL,PRHO)
-    ZQR =  ZCR* ZQL*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + ZFNEG*(0.5-SIGN(0.5,ZCR)) ) 
+    ZQR =  PCR* ZQL*( ZFPOS*(0.5+SIGN(0.5,PCR)) + ZFNEG*(0.5-SIGN(0.5,PCR)) ) 
     dyf(PR,ZQR)
     !#acc end region 
 
 !!$      PR = ZDMQ + ZQL0 + ZDQ +  ZQ60 +  ZQL + ZQR  + ZQ6 + ZQR0 + ZFPOS + ZFNEG
 
-!!$   PR = DYF( ZCR*MYM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                             ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+!!$   PR = DYF( PCR*MYM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                             ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 !!$!
 
+#ifdef JUAN_ACC_01_Y
 !$acc end region 
-! acc update host (pr)
 !$acc end data region 
+#endif
 
    CALL GET_HALO(PR,HDIR="01_Y")
 
@@ -1160,7 +1213,7 @@ END FUNCTION PPM_01_Y
 !-------------------------------------------------------------------------------
 !
 !     ########################################################################
-      FUNCTION PPM_01_Z(KGRID, PSRC, ZCR, PRHO, PTSTEP) RESULT(PR)
+      FUNCTION PPM_01_Z(KGRID, PSRC, PCR, PRHO, PTSTEP) RESULT(PR)
 !     ########################################################################
 !!
 !!****  PPM_01_Z - PPM_01 fully monotonic PPM advection scheme in Z direction
@@ -1179,7 +1232,6 @@ USE MODI_GET_HALO
 !
 USE MODD_CONF
 USE MODD_PARAMETERS
-!USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
 USE MODE_MPPDB
 !
 IMPLICIT NONE
@@ -1189,29 +1241,27 @@ IMPLICIT NONE
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR &
 !
 !*       0.2   Declarations of local variables :
 !
-INTEGER:: IKB    ! Begining useful area in x,y,z directions
-INTEGER:: IKE    ! End useful area in x,y,z directions
-!
 ! terms used in parabolic interpolation, dmq, qL, qR, dq, q6
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZQL,ZQR
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZDQ,ZQ6
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZDMQ
+                                                     , ZQL, ZQR, ZDQ, ZQ6, ZDMQ &
 !
 ! extra variables for the initial guess of parabolae parameters
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZQL0,ZQR0,ZQ60
+                                                     , ZQL0,ZQR0,ZQ60 &
 !
 ! advection fluxes
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFPOS, ZFNEG
+                                                     , ZFPOS, ZFNEG
+!
+INTEGER:: IKB    ! Begining useful area in x,y,z directions
+INTEGER:: IKE    ! End useful area in x,y,z directions
 !
 !JUAN ACC
 INTEGER                          :: I,J,K ,IIU,IJU,IKU
@@ -1236,9 +1286,11 @@ IIU=size(psrc,1)
 IJU=size(psrc,2)
 IKU=size(psrc,3)
 
-
-!$acc data region local (ZDMQ,ZQL0,ZQR0,ZDQ,ZQ60,ZQL,ZQR,ZQ6,ZFPOS,ZFNEG) copyin (psrc,zcr,prho) copyout(pr)
+#define JUAN_ACC_01_Z
+#ifdef JUAN_ACC_01_Z
+!$acc data region local (ZDMQ,ZQL0,ZQR0,ZDQ,ZQ60,ZQL,ZQR,ZQ6,ZFPOS,ZFNEG) copyin (psrc,pcr,prho) copyout(pr)
 !$acc region  
+#endif
 !
 !-------------------------------------------------------------------------------
 !
@@ -1312,40 +1364,40 @@ ZDQ = ZQR - ZQL
 !
 ! and finally calculate fluxes for the advection
 !
-ZFPOS(:,:,IKB+1:IKE+1) = ZQR(:,:,IKB:IKE) - 0.5*ZCR(:,:,IKB+1:IKE+1) * &            
-     (ZDQ(:,:,IKB:IKE) - (1.0 - 2.0*ZCR(:,:,IKB+1:IKE+1)/3.0)        &
+ZFPOS(:,:,IKB+1:IKE+1) = ZQR(:,:,IKB:IKE) - 0.5*PCR(:,:,IKB+1:IKE+1) * &            
+     (ZDQ(:,:,IKB:IKE) - (1.0 - 2.0*PCR(:,:,IKB+1:IKE+1)/3.0)        &
      * ZQ6(:,:,IKB:IKE))
 !
 ! advection flux at open boundary when u(IKB) > 0
-ZFPOS(:,:,IKB) = (PSRC(:,:,IKB-1) - ZQR(:,:,IKB-1))*ZCR(:,:,IKB) + &
+ZFPOS(:,:,IKB) = (PSRC(:,:,IKB-1) - ZQR(:,:,IKB-1))*PCR(:,:,IKB) + &
                  ZQR(:,:,IKB-1)
 !
 ! PPOSX(IKB-1) is not important for the calc of advection so 
 ! we set it to 0
 ZFPOS(:,:,IKB-1) = 0.0
 !
-ZFNEG(:,:,IKB-1:IKE) = ZQL(:,:,IKB-1:IKE) - 0.5*ZCR(:,:,IKB-1:IKE) *      &            
-     ( ZDQ(:,:,IKB-1:IKE) + (1.0 + 2.0*ZCR(:,:,IKB-1:IKE)/3.0) * &
+ZFNEG(:,:,IKB-1:IKE) = ZQL(:,:,IKB-1:IKE) - 0.5*PCR(:,:,IKB-1:IKE) *      &            
+     ( ZDQ(:,:,IKB-1:IKE) + (1.0 + 2.0*PCR(:,:,IKB-1:IKE)/3.0) * &
        ZQ6(:,:,IKB-1:IKE) )
 !
 ! advection flux at open boundary when u(IKE+1) < 0
-ZFNEG(:,:,IKE+1) = (ZQR(:,:,IKE)-PSRC(:,:,IKE+1))*ZCR(:,:,IKE+1) + &
+ZFNEG(:,:,IKE+1) = (ZQR(:,:,IKE)-PSRC(:,:,IKE+1))*PCR(:,:,IKE+1) + &
                    ZQR(:,:,IKE)
 !
 ! advect the actual field in Z direction by W*dt
 !
 
     mzm(ZQL,PRHO)
-    ZQR =  ZCR* ZQL*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + ZFNEG*(0.5-SIGN(0.5,ZCR)) ) 
+    ZQR =  PCR* ZQL*( ZFPOS*(0.5+SIGN(0.5,PCR)) + ZFNEG*(0.5-SIGN(0.5,PCR)) ) 
     dzf(PR,ZQR)
 
-!!$PR = DZF( ZCR*MZM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                          ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
-
+!!$PR = DZF( PCR*MZM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                          ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 
+#ifdef JUAN_ACC_01_Z
 !$acc end region
-! acc update host (pr)
 !$acc end data region  
+#endif
 
 CALL GET_HALO(PR)
 CALL MPPDB_CHECK3DM("PPM::PPM_01_Z ::PR",PRECISION,PR)
@@ -1418,7 +1470,7 @@ END FUNCTION PPM_01_Z
 !-------------------------------------------------------------------------------
 !
 !     ########################################################################
-      FUNCTION PPM_S0_X(HLBCX, KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+      FUNCTION PPM_S0_X(HLBCX, KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 !     ########################################################################
 !!
@@ -1454,13 +1506,13 @@ CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX  ! X direction LBC type
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 !*       0.2   Declarations of local variables :
 !
@@ -1468,10 +1520,10 @@ INTEGER:: IIB,IJB    ! Begining useful area in x,y,z directions
 INTEGER:: IIE,IJE    ! End useful area in x,y,z directions
 !
 ! advection fluxes
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFPOS, ZFNEG
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZFPOS, ZFNEG
 !
 ! variable at cell edges
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZPHAT
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZPHAT
 !
 !BEG JUAN PPM_LL
 TYPE(HALO2LIST_ll), POINTER      :: TZ_PSRC_HALO2_ll         ! halo2 for PSRC
@@ -1479,9 +1531,12 @@ INTEGER                          :: ILUOUT,IRESP             ! for prints
 INTEGER                          :: IJS,IJN
 !END JUAN PPM_LL
 !JUAN ACC
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZRHO_MXM , ZCR_MXM , ZCR_DXF
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZRHO_MXM , ZCR_MXM , ZCR_DXF
 INTEGER                          :: I,J,K ,IIU,IJU,IKU
 LOGICAL                          :: GWEST , GEAST
+!
+REAL, DIMENSION(SIZE(PCR,2),SIZE(PCR,3))             :: ZPSRC_HALO2_WEST
+!
 ! inline shuman with macro 
 #define dxf(PDXF,PA) PDXF(1:IIU-1,:,:) = PA(2:IIU,:,:) - PA(1:IIU-1,:,:) ; PDXF(IIU,:,:)    = PDXF(2*JPHEXT,:,:) ! DXF(PDXF,PA)
 #define mxm(PMXM,PA) PMXM(2:IIU,:,:) = 0.5*( PA(2:IIU,:,:)+PA(1:IIU-1,:,:) ) ; PMXM(1,:,:) = PMXM(IIU-2*JPHEXT+1,:,:) !  MXM(PMXM,PA)
@@ -1518,8 +1573,14 @@ IF(NHALO /= 1) THEN
 ENDIF
 !
 CALL GET_HALO2(PSRC,TZ_PSRC_HALO2_ll)
-!$acc data region local (PR,ZPHAT,ZFPOS,ZFNEG,ZRHO_MXM,ZCR_MXM,ZCR_DXF) copyin (psrc,zcr,prho)
+ZPSRC_HALO2_WEST(:,:) = TZ_PSRC_HALO2_ll%HALO2%WEST(:,:)
+!
+#define JUAN_ACC_S0_X
+#ifdef JUAN_ACC_S0_X
+!$acc data region local (ZPHAT,ZFPOS,ZFNEG,ZRHO_MXM,ZCR_MXM,ZCR_DXF) copyin (psrc,pcr,prho,ZPSRC_HALO2_WEST)  copyout(pr) 
 !$acc region
+#endif
+!
 ZPHAT=PSRC
 ZFPOS=PSRC
 ZFNEG=PSRC
@@ -1572,16 +1633,16 @@ ZPHAT(IIB+1:IIE,IJS:IJN,:) = ( 7.0 * &
 !!$CALL  GET_HALO(ZPHAT,HDIR="Z0_X")
 !!$
 !!$   ZFPOS(IIB:IIE,IJS:IJN,:) = ZPHAT(IIB:IIE,IJS:IJN,:) - & 
-!!$       ZCR(IIB:IIE,IJS:IJN,:)*(ZPHAT(IIB:IIE,IJS:IJN,:) - PSRC(IIB-1:IIE-1,IJS:IJN,:)) - &
-!!$       ZCR(IIB:IIE,IJS:IJN,:)*(1.0 - ZCR(IIB:IIE,IJS:IJN,:)) * &
+!!$       PCR(IIB:IIE,IJS:IJN,:)*(ZPHAT(IIB:IIE,IJS:IJN,:) - PSRC(IIB-1:IIE-1,IJS:IJN,:)) - &
+!!$       PCR(IIB:IIE,IJS:IJN,:)*(1.0 - PCR(IIB:IIE,IJS:IJN,:)) * &
 !!$       (ZPHAT(IIB-1:IIE-1,IJS:IJN,:) - 2.0*PSRC(IIB-1:IIE-1,IJS:IJN,:) + ZPHAT(IIB:IIE,IJS:IJN,:))
 !!$!
 !!$!
 !!$CALL GET_HALO(ZFPOS,HDIR="Z0_X") ! JUAN
 !!$!
 !!$   ZFNEG(IIB:IIE,IJS:IJN,:) = ZPHAT(IIB:IIE,IJS:IJN,:) + & 
-!!$        ZCR(IIB:IIE,IJS:IJN,:)*(ZPHAT(IIB:IIE,IJS:IJN,:) - PSRC(IIB:IIE,IJS:IJN,:)) + &
-!!$        ZCR(IIB:IIE,IJS:IJN,:)*(1.0 + ZCR(IIB:IIE,IJS:IJN,:)) * &
+!!$        PCR(IIB:IIE,IJS:IJN,:)*(ZPHAT(IIB:IIE,IJS:IJN,:) - PSRC(IIB:IIE,IJS:IJN,:)) + &
+!!$        PCR(IIB:IIE,IJS:IJN,:)*(1.0 + PCR(IIB:IIE,IJS:IJN,:)) * &
 !!$        (ZPHAT(IIB:IIE,IJS:IJN,:) - 2.0*PSRC(IIB:IIE,IJS:IJN,:) + ZPHAT(IIB+1:IIE+1,IJS:IJN,:))
 !!$!
 !!$! define fluxes for CYCL BC outside physical domain
@@ -1591,8 +1652,8 @@ ZPHAT(IIB+1:IIE,IJS:IJN,:) = ( 7.0 * &
 !!$! calculate the advection
 !!$!
 !!$   PR = PSRC * PRHO - &
-!!$        DXF( ZCR*MXM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                             ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+!!$        DXF( PCR*MXM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                             ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 !!$   CALL GET_HALO(PR,HDIR="S0_X") ! JUAN
 !!$   CALL MPPDB_CHECK3DM("PPM::PPM_S0_X CYCL ::PR",PRECISION,PR)
 !!$!
@@ -1601,12 +1662,16 @@ ZPHAT(IIB+1:IIE,IJS:IJN,:) = ( 7.0 * &
 !
 !  WEST BOUND 
 !
-!TEMPO_JUAN  IF (.NOT. GWEST) THEN
-!TEMPO_JUAN   ZPHAT(IIB  ,IJS:IJN,:) = ( 7.0 * &
-!TEMPO_JUAN                      ( PSRC(IIB  ,IJS:IJN,:) + PSRC(IIB-1,IJS:IJN,:)                  ) - &
-!TEMPO_JUAN                      ( PSRC(IIB+1,IJS:IJN,:) + TZ_PSRC_HALO2_ll%HALO2%WEST(IJS:IJN,:) ) ) / 12.0
-!TEMPO_JUAN! <=>  WEST BOUND     ( PSRC(IIB+1,IJS:IJN,:) + PSRC(IIB-2,IJS:IJN,:)                  ) ) / 12.0
-!TEMPO_JUAN  ENDIF
+#ifdef TEMPO_JUAN  
+IF (.NOT. GWEST) THEN
+   ZPHAT(IIB  ,IJS:IJN,:) = ( 7.0 * &
+                      ( PSRC(IIB  ,IJS:IJN,:) + PSRC(IIB-1,IJS:IJN,:)                  ) - &
+                      ( PSRC(IIB+1,IJS:IJN,:) + ZPSRC_HALO2_WEST(IJS:IJN,:)            ) ) / 12.0
+!                     ( PSRC(IIB+1,IJS:IJN,:) + TZ_PSRC_HALO2_ll%HALO2%WEST(IJS:IJN,:) ) ) / 12.0
+! <=>  WEST BOUND     ( PSRC(IIB+1,IJS:IJN,:) + PSRC(IIB-2,IJS:IJN,:)                  ) ) / 12.0
+ENDIF
+!TEMPO_JUAN 
+#endif 
 !
 !   update ZPHAT HALO before next/further  utilisation 
 !
@@ -1626,15 +1691,15 @@ ZPHAT(IIB+1:IIE,IJS:IJN,:) = ( 7.0 * &
 ! update ZPHAT HALO before next/further  utilisation
 !
    ZFPOS(IIB:IIE,IJS:IJN,:) = ZPHAT(IIB:IIE,IJS:IJN,:) - & 
-        ZCR(IIB:IIE,IJS:IJN,:)*(ZPHAT(IIB:IIE,IJS:IJN,:) - PSRC(IIB-1:IIE-1,IJS:IJN,:)) - &
-        ZCR(IIB:IIE,IJS:IJN,:)*(1.0 - ZCR(IIB:IIE,IJS:IJN,:)) * &
+        PCR(IIB:IIE,IJS:IJN,:)*(ZPHAT(IIB:IIE,IJS:IJN,:) - PSRC(IIB-1:IIE-1,IJS:IJN,:)) - &
+        PCR(IIB:IIE,IJS:IJN,:)*(1.0 - PCR(IIB:IIE,IJS:IJN,:)) * &
         (ZPHAT(IIB-1:IIE-1,IJS:IJN,:) - 2.0*PSRC(IIB-1:IIE-1,IJS:IJN,:) + ZPHAT(IIB:IIE,IJS:IJN,:))
 !
 !TEMPO_JUAN CALL GET_HALO(ZFPOS,HDIR="Z0_X") ! JUAN
 !
 ! positive flux on the WEST boundary
   IF (GWEST) THEN
-   ZFPOS(IIB,IJS:IJN,:) = (PSRC(IIB-1,IJS:IJN,:) - ZPHAT(IIB,IJS:IJN,:))*ZCR(IIB,IJS:IJN,:) + &
+   ZFPOS(IIB,IJS:IJN,:) = (PSRC(IIB-1,IJS:IJN,:) - ZPHAT(IIB,IJS:IJN,:))*PCR(IIB,IJS:IJN,:) + &
                      ZPHAT(IIB,IJS:IJN,:) 
 ! this is not used
    ZFPOS(IIB-1,IJS:IJN,:) = 0.0
@@ -1643,40 +1708,40 @@ ZPHAT(IIB+1:IIE,IJS:IJN,:) = ( 7.0 * &
 ! negative fluxes
 
    ZFNEG(IIB:IIE,IJS:IJN,:) = ZPHAT(IIB:IIE,IJS:IJN,:) + & 
-        ZCR(IIB:IIE,IJS:IJN,:)*(ZPHAT(IIB:IIE,IJS:IJN,:) - PSRC(IIB:IIE,IJS:IJN,:)) + &
-        ZCR(IIB:IIE,IJS:IJN,:)*(1.0 + ZCR(IIB:IIE,IJS:IJN,:)) * &
+        PCR(IIB:IIE,IJS:IJN,:)*(ZPHAT(IIB:IIE,IJS:IJN,:) - PSRC(IIB:IIE,IJS:IJN,:)) + &
+        PCR(IIB:IIE,IJS:IJN,:)*(1.0 + PCR(IIB:IIE,IJS:IJN,:)) * &
         (ZPHAT(IIB:IIE,IJS:IJN,:) - 2.0*PSRC(IIB:IIE,IJS:IJN,:) + ZPHAT(IIB+1:IIE+1,IJS:IJN,:))
 !
-!TEMPO_JUAN   CALL GET_HALO(ZFNEG,HDIR="Z0_X") ! JUAN
+!TEMPO_JUAN CALL GET_HALO(ZFNEG,HDIR="Z0_X") ! JUAN
 !
   IF (GEAST) THEN
 !
-! in OPEN case ZCR(IIB-1) is not used, so we also set ZFNEG(IIB-1) = 0
+! in OPEN case PCR(IIB-1) is not used, so we also set ZFNEG(IIB-1) = 0
 !
    ZFNEG(IIB-1,IJS:IJN,:) = 0.0
 !
 ! modified negative flux on EAST boundary. We use linear function instead of a
 ! parabola to represent the tracer field, so it simplifies the flux expresion
 !
-   ZFNEG(IIE+1,IJS:IJN,:) = (ZPHAT(IIE+1,IJS:IJN,:) - PSRC(IIE+1,IJS:IJN,:))*ZCR(IIE+1,IJS:IJN,:) + &
+   ZFNEG(IIE+1,IJS:IJN,:) = (ZPHAT(IIE+1,IJS:IJN,:) - PSRC(IIE+1,IJS:IJN,:))*PCR(IIE+1,IJS:IJN,:) + &
                        ZPHAT(IIE+1,IJS:IJN,:)
   ENDIF
 !
 ! calculate the advection
 !
    mxm(ZRHO_MXM,PRHO)
-   ZCR_MXM =  ZCR* ZRHO_MXM*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + ZFNEG*(0.5-SIGN(0.5,ZCR)) ) 
+   ZCR_MXM =  PCR* ZRHO_MXM*( ZFPOS*(0.5+SIGN(0.5,PCR)) + ZFNEG*(0.5-SIGN(0.5,PCR)) ) 
    dxf(ZCR_DXF,ZCR_MXM)
    PR = PSRC * PRHO - ZCR_DXF
 !!$   PR = PSRC * PRHO - &
-!!$        DXF( ZCR*MXM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                             ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+!!$        DXF( PCR*MXM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                             ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 !
 !
 ! in OPEN case fix boundary conditions
 !
   IF (GWEST) THEN
-   WHERE ( ZCR(IIB,IJS:IJN,:) <= 0. ) !  OUTFLOW condition
+   WHERE ( PCR(IIB,IJS:IJN,:) <= 0. ) !  OUTFLOW condition
       PR(IIB-1,IJS:IJN,:) = 2.*PR(IIB,IJS:IJN,:) - PR(IIB+1,IJS:IJN,:)
    ELSEWHERE
       PR(IIB-1,IJS:IJN,:) = PR(IIB,IJS:IJN,:)
@@ -1684,16 +1749,17 @@ ZPHAT(IIB+1:IIE,IJS:IJN,:) = ( 7.0 * &
   ENDIF
 !
   IF (GEAST) THEN 
-   WHERE ( ZCR(IIE,IJS:IJN,:) >= 0. ) !  OUTFLOW condition
+   WHERE ( PCR(IIE,IJS:IJN,:) >= 0. ) !  OUTFLOW condition
       PR(IIE+1,IJS:IJN,:) = 2.*PR(IIE,IJS:IJN,:) - PR(IIE-1,IJS:IJN,:)
    ELSEWHERE
       PR(IIE+1,IJS:IJN,:) = PR(IIE,IJS:IJN,:)
    END WHERE
   ENDIF
-
+!
+#ifdef JUAN_ACC_S0_X
 !$acc end region 
-!$acc update host (pr)
 !$acc end data region 
+#endif
 !
 CALL GET_HALO(PR,HDIR="S0_X") 
 CALL MPPDB_CHECK3DM("PPM::PPM_S0_X OPEN ::PR",PRECISION,PR)
@@ -1709,7 +1775,7 @@ END FUNCTION PPM_S0_X
 !-------------------------------------------------------------------------------
 !
 !     ########################################################################
-      FUNCTION PPM_S0_Y(HLBCY, KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+      FUNCTION PPM_S0_Y(HLBCY, KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 !     ########################################################################
 !!
@@ -1743,13 +1809,13 @@ CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY ! Y direction LBC type
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 !*       0.2   Declarations of local variables :
 !
@@ -1757,21 +1823,23 @@ INTEGER:: IIB,IJB    ! Begining useful area in x,y,z directions
 INTEGER:: IIE,IJE    ! End useful area in x,y,z directions
 !
 ! advection fluxes
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFPOS, ZFNEG
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZFPOS, ZFNEG
 !
 ! variable at cell edges
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZPHAT
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZPHAT
 !
 !BEG JUAN PPM_LL
 TYPE(HALO2LIST_ll), POINTER      :: TZ_PSRC_HALO2_ll         ! halo2 for PSRC
-TYPE(HALO2LIST_ll), POINTER      :: TZ_PHAT_HALO2_ll         ! halo2 for ZPHAT
 INTEGER                          :: ILUOUT,IRESP             ! for prints
 INTEGER                          :: IIW,IIA
 !END JUAN PPM_LL
 !JUAN ACC
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZRHO_MYM , ZCR_MYM , ZCR_DYF
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZRHO_MYM , ZCR_MYM , ZCR_DYF
 INTEGER                          :: I,J,K ,IIU,IJU,IKU
 LOGICAL                          :: GSOUTH , GNORTH
+!
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,3))             :: ZPSRC_HALO2_SOUTH
+!
 ! inline shuman with macro 
 #define  dyf(PDYF,PA) PDYF(:,1:IJU-1,:) = PA(:,2:IJU,:) - PA(:,1:IJU-1,:); PDYF(:,IJU,:) = PDYF(:,2*JPHEXT,:) !   DYF(PDYF,PA)
 #define mym(PMYM,PA) PMYM(:,2:IJU,:) = 0.5*( PA(:,2:IJU,:)+PA(:,1:IJU-1,:) ) ; PMYM(:,1,:) = PMYM(:,IJU-2*JPHEXT+1,:) !  MYM(PMYM,PA)
@@ -1802,12 +1870,16 @@ IF ( L2D ) THEN
 END IF
 !
 CALL GET_HALO2(PSRC,TZ_PSRC_HALO2_ll)
+ZPSRC_HALO2_SOUTH(:,:) = TZ_PSRC_HALO2_ll%HALO2%SOUTH(:,:)
 !
 ! Initialize with relalistic value all work array 
 !
-!$acc data region local(PR,ZPHAT,ZFPOS,ZFNEG,ZRHO_MYM,ZCR_MYM,ZCR_DYF) copyin (psrc,zcr,prho)
+#define JUAN_ACC_S0_Y
+#ifdef JUAN_ACC_S0_Y
+!$acc data region local(ZPHAT,ZFPOS,ZFNEG,ZRHO_MYM,ZCR_MYM,ZCR_DYF) copyin (psrc,pcr,prho,zpsrc_halo2_south)  copyout(pr) 
 !$acc region
-
+#endif
+!
 ZPHAT=PSRC
 ZFPOS=PSRC
 ZFNEG=PSRC
@@ -1848,16 +1920,16 @@ ZPHAT(IIW:IIA,IJB+1:IJE,:) = (7.0 * &
 !!$! calculate the fluxes:
 !!$!
 !!$   ZFPOS(IIW:IIA,IJB:IJE,:) = ZPHAT(IIW:IIA,IJB:IJE,:) - & 
-!!$        ZCR(IIW:IIA,IJB:IJE,:)*(ZPHAT(IIW:IIA,IJB:IJE,:) - PSRC(IIW:IIA,IJB-1:IJE-1,:)) - &
-!!$        ZCR(IIW:IIA,IJB:IJE,:)*(1.0 - ZCR(IIW:IIA,IJB:IJE,:)) * &
+!!$        PCR(IIW:IIA,IJB:IJE,:)*(ZPHAT(IIW:IIA,IJB:IJE,:) - PSRC(IIW:IIA,IJB-1:IJE-1,:)) - &
+!!$        PCR(IIW:IIA,IJB:IJE,:)*(1.0 - PCR(IIW:IIA,IJB:IJE,:)) * &
 !!$        (ZPHAT(IIW:IIA,IJB-1:IJE-1,:) - 2.0*PSRC(IIW:IIA,IJB-1:IJE-1,:) + ZPHAT(IIW:IIA,IJB:IJE,:))
 !!$!
 !!$!
 !!$CALL GET_HALO(ZFPOS,HDIR="Z0_Y") ! JUAN
 !!$!
 !!$   ZFPOS(IIW:IIA,IJB:IJE,:) = ZPHAT(IIW:IIA,IJB:IJE,:) - & 
-!!$        ZCR(IIW:IIA,IJB:IJE,:)*(ZPHAT(IIW:IIA,IJB:IJE,:) - PSRC(IIW:IIA,IJB-1:IJE-1,:)) - &
-!!$        ZCR(IIW:IIA,IJB:IJE,:)*(1.0 - ZCR(IIW:IIA,IJB:IJE,:)) * &
+!!$        PCR(IIW:IIA,IJB:IJE,:)*(ZPHAT(IIW:IIA,IJB:IJE,:) - PSRC(IIW:IIA,IJB-1:IJE-1,:)) - &
+!!$        PCR(IIW:IIA,IJB:IJE,:)*(1.0 - PCR(IIW:IIA,IJB:IJE,:)) * &
 !!$        (ZPHAT(IIW:IIA,IJB-1:IJE-1,:) - 2.0*PSRC(IIW:IIA,IJB-1:IJE-1,:) + ZPHAT(IIW:IIA,IJB:IJE,:))
 !!$!
 !!$!
@@ -1867,8 +1939,8 @@ ZPHAT(IIW:IIA,IJB+1:IJE,:) = (7.0 * &
 !!$! calculate the advection
 !!$!
 !!$   PR = PSRC * PRHO - &
-!!$        DYF( ZCR*MYM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                             ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+!!$        DYF( PCR*MYM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                             ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 !!$   CALL GET_HALO(PR,HDIR="S0_Y") ! JUAN
 !!$   CALL MPPDB_CHECK3DM("PPM::PPM_S0_Y CYCL ::PR",PRECISION,PR
 !!$!
@@ -1877,14 +1949,18 @@ ZPHAT(IIW:IIA,IJB+1:IJE,:) = (7.0 * &
 !
 !  SOUTH BOUND 
 !
-!TEMPO_JUAN  IF ( .NOT. GSOUTH) THEN
-!TEMPO_JUAN   ZPHAT(IIW:IIA,IJB  ,:) = (7.0 * &
-!TEMPO_JUAN                      (PSRC(IIW:IIA,IJB  ,:) + PSRC(IIW:IIA,IJB-1,:)) - &
-!TEMPO_JUAN                      (PSRC(IIW:IIA,IJB+1,:) + TZ_PSRC_HALO2_ll%HALO2%SOUTH(IIW:IIA,:) )) / 12.0
-!TEMPO_JUAN! <=> SOUTH BOUND     (PSRC(IIW:IIA,IJB+1,:) + PSRC(IIW:IIA,IJB-2,:)                   )) / 12.0
-!TEMPO_JUAN  ENDIF
+#ifdef TEMPO_JUAN  
+IF ( .NOT. GSOUTH) THEN
+   ZPHAT(IIW:IIA,IJB  ,:) = (7.0 * &
+                      (PSRC(IIW:IIA,IJB  ,:) + PSRC(IIW:IIA,IJB-1,:)) - &
+                      (PSRC(IIW:IIA,IJB+1,:) + ZPSRC_HALO2_SOUTH(IIW:IIA,:)            )) / 12.0
+!                     (PSRC(IIW:IIA,IJB+1,:) + TZ_PSRC_HALO2_ll%HALO2%SOUTH(IIW:IIA,:) )) / 12.0
+! <=> SOUTH BOUND     (PSRC(IIW:IIA,IJB+1,:) + PSRC(IIW:IIA,IJB-2,:)                   )) / 12.0
+ENDIF
+!TEMPO_JUAN  
+#endif
 !
-!JUAN_TEMP CALL  GET_HALO(ZPHAT,HDIR="Z0_Y")
+!TEMPO_JUAN CALL  GET_HALO(ZPHAT,HDIR="Z0_Y")
 !
   IF (GSOUTH) THEN
    ZPHAT(IIW:IIA,IJB  ,:) = 0.5*(PSRC(IIW:IIA,IJB-1,:) + PSRC(IIW:IIA,IJB,:))
@@ -1905,15 +1981,15 @@ ZPHAT(IIW:IIA,IJB+1:IJE,:) = (7.0 * &
 ! positive fluxes
 !
    ZFPOS(IIW:IIA,IJB:IJE,:) = ZPHAT(IIW:IIA,IJB:IJE,:) - & 
-        ZCR(IIW:IIA,IJB:IJE,:)*(ZPHAT(IIW:IIA,IJB:IJE,:) - PSRC(IIW:IIA,IJB-1:IJE-1,:)) - &
-        ZCR(IIW:IIA,IJB:IJE,:)*(1.0 - ZCR(IIW:IIA,IJB:IJE,:)) * &
+        PCR(IIW:IIA,IJB:IJE,:)*(ZPHAT(IIW:IIA,IJB:IJE,:) - PSRC(IIW:IIA,IJB-1:IJE-1,:)) - &
+        PCR(IIW:IIA,IJB:IJE,:)*(1.0 - PCR(IIW:IIA,IJB:IJE,:)) * &
         (ZPHAT(IIW:IIA,IJB-1:IJE-1,:) - 2.0*PSRC(IIW:IIA,IJB-1:IJE-1,:) + ZPHAT(IIW:IIA,IJB:IJE,:))
 !
-!JUAN_TEMP CALL GET_HALO(ZFPOS,HDIR="Z0_Y") ! JUAN
+!TEMPO_JUAN CALL GET_HALO(ZFPOS,HDIR="Z0_Y") ! JUAN
 !
 ! positive flux on the SOUTH boundary
   IF (GSOUTH) THEN
-   ZFPOS(IIW:IIA,IJB,:) = (PSRC(IIW:IIA,IJB-1,:) - ZPHAT(IIW:IIA,IJB,:))*ZCR(IIW:IIA,IJB,:) + &
+   ZFPOS(IIW:IIA,IJB,:) = (PSRC(IIW:IIA,IJB-1,:) - ZPHAT(IIW:IIA,IJB,:))*PCR(IIW:IIA,IJB,:) + &
                      ZPHAT(IIW:IIA,IJB,:)
 !
 ! this is not used
@@ -1923,36 +1999,36 @@ ZPHAT(IIW:IIA,IJB+1:IJE,:) = (7.0 * &
 ! negative fluxes
 !
    ZFNEG(IIW:IIA,IJB:IJE,:) = ZPHAT(IIW:IIA,IJB:IJE,:) + & 
-        ZCR(IIW:IIA,IJB:IJE,:)*(ZPHAT(IIW:IIA,IJB:IJE,:) - PSRC(IIW:IIA,IJB:IJE,:)) + &
-        ZCR(IIW:IIA,IJB:IJE,:)*(1.0 + ZCR(IIW:IIA,IJB:IJE,:)) * &
+        PCR(IIW:IIA,IJB:IJE,:)*(ZPHAT(IIW:IIA,IJB:IJE,:) - PSRC(IIW:IIA,IJB:IJE,:)) + &
+        PCR(IIW:IIA,IJB:IJE,:)*(1.0 + PCR(IIW:IIA,IJB:IJE,:)) * &
         (ZPHAT(IIW:IIA,IJB:IJE,:) - 2.0*PSRC(IIW:IIA,IJB:IJE,:) +ZPHAT(IIW:IIA,IJB+1:IJE+1,:))
 !
-   !JUAN_TEMP CALL GET_HALO(ZFNEG,HDIR="Z0_Y") ! JUAN
+   !TEMPO_JUAN CALL GET_HALO(ZFNEG,HDIR="Z0_Y") ! JUAN
 !
   IF (GNORTH) THEN
 ! this is not used
    ZFNEG(IIW:IIA,IJB-1,:) = 0.0
 !
 ! negative flux on the NORTH boundary
-   ZFNEG(IIW:IIA,IJE+1,:) = (ZPHAT(IIW:IIA,IJE+1,:) - PSRC(IIW:IIA,IJE+1,:))*ZCR(IIW:IIA,IJE+1,:) + &
+   ZFNEG(IIW:IIA,IJE+1,:) = (ZPHAT(IIW:IIA,IJE+1,:) - PSRC(IIW:IIA,IJE+1,:))*PCR(IIW:IIA,IJE+1,:) + &
                        ZPHAT(IIW:IIA,IJE+1,:)
   ENDIF
 !
 ! calculate the advection
 !
    mym(ZRHO_MYM,PRHO)
-   ZCR_MYM =  ZCR* ZRHO_MYM*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + ZFNEG*(0.5-SIGN(0.5,ZCR)) ) 
+   ZCR_MYM =  PCR* ZRHO_MYM*( ZFPOS*(0.5+SIGN(0.5,PCR)) + ZFNEG*(0.5-SIGN(0.5,PCR)) ) 
    dyf(ZCR_DYF,ZCR_MYM)
    PR = PSRC * PRHO - ZCR_DYF
 !!$   PR = PSRC * PRHO - &
-!!$        DYF( ZCR*MYM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-!!$                             ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+!!$        DYF( PCR*MYM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+!!$                             ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 !
 !
 ! in OPEN case fix boundary conditions
 !
   IF (GSOUTH) THEN
-   WHERE ( ZCR(IIW:IIA,IJB,:) <= 0. ) !  OUTFLOW condition
+   WHERE ( PCR(IIW:IIA,IJB,:) <= 0. ) !  OUTFLOW condition
       PR(IIW:IIA,IJB-1,:) = 1.0 * 2.*PR(IIW:IIA,IJB,:) - PR(IIW:IIA,IJB+1,:)
    ELSEWHERE
       PR(IIW:IIA,IJB-1,:) = PR(IIW:IIA,IJB,:) 
@@ -1960,16 +2036,17 @@ ZPHAT(IIW:IIA,IJB+1:IJE,:) = (7.0 * &
   ENDIF
 !
   IF (GNORTH) THEN
-   WHERE ( ZCR(IIW:IIA,IJE,:) >= 0. ) !  OUTFLOW condition
+   WHERE ( PCR(IIW:IIA,IJE,:) >= 0. ) !  OUTFLOW condition
       PR(IIW:IIA,IJE+1,:) = 1.0 * 2.*PR(IIW:IIA,IJE,:) - PR(IIW:IIA,IJE-1,:)
    ELSEWHERE
       PR(IIW:IIA,IJE+1,:) = PR(IIW:IIA,IJE,:) 
    END WHERE
   ENDIF
-
+!
+#ifdef JUAN_ACC_S0_Y
 !$acc end region 
-!$acc update host (pr)
 !$acc end data region 
+#endif
 !
    CALL GET_HALO(PR,HDIR="S0_Y") 
    CALL MPPDB_CHECK3DM("PPM::PPM_S0_Y OPEN ::PR",PRECISION,PR)  
@@ -1986,7 +2063,7 @@ END FUNCTION PPM_S0_Y
 !-------------------------------------------------------------------------------
 !
 !     ########################################################################
-      FUNCTION PPM_S0_Z(KGRID, PSRC, ZCR, PRHO, PTSTEP) &
+      FUNCTION PPM_S0_Z(KGRID, PSRC, PCR, PRHO, PTSTEP) &
                RESULT(PR)
 !     ########################################################################
 !!
@@ -2006,7 +2083,6 @@ USE MODI_GET_HALO
 !
 USE MODD_CONF
 USE MODD_PARAMETERS
-!USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
 USE MODE_MPPDB
 !
 IMPLICIT NONE
@@ -2016,13 +2092,13 @@ IMPLICIT NONE
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 !*       0.2   Declarations of local variables :
 !
@@ -2030,10 +2106,10 @@ INTEGER:: IKB    ! Begining useful area in x,y,z directions
 INTEGER:: IKE    ! End useful area in x,y,z directions
 !
 ! advection fluxes
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFPOS, ZFNEG
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZFPOS, ZFNEG
 !
 ! interpolated variable at cell edges
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZPHAT
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZPHAT
 !
 !-------------------------------------------------------------------------------
 !
@@ -2064,14 +2140,14 @@ ZPHAT(:,:,IKE+1) = 0.5*(PSRC(:,:,IKE) + PSRC(:,:,IKE+1))
 ! (for inflow or outflow situation)
 !
 ZFPOS(:,:,IKB+1:IKE+1) = ZPHAT(:,:,IKB+1:IKE+1) - & 
-     ZCR(:,:,IKB+1:IKE+1)*(ZPHAT(:,:,IKB+1:IKE+1) - PSRC(:,:,IKB:IKE)) - &
-     ZCR(:,:,IKB+1:IKE+1)*(1.0 - ZCR(:,:,IKB+1:IKE+1)) * &
+     PCR(:,:,IKB+1:IKE+1)*(ZPHAT(:,:,IKB+1:IKE+1) - PSRC(:,:,IKB:IKE)) - &
+     PCR(:,:,IKB+1:IKE+1)*(1.0 - PCR(:,:,IKB+1:IKE+1)) * &
      (ZPHAT(:,:,IKB:IKE) - 2.0*PSRC(:,:,IKB:IKE) + ZPHAT(:,:,IKB+1:IKE+1))
 !
 !!$CALL GET_HALO(ZFPOS) ! JUAN
 !
 ! positive flux on the BOTTOM boundary
-ZFPOS(:,:,IKB) = (PSRC(:,:,IKB-1) - ZPHAT(:,:,IKB))*ZCR(:,:,IKB) + &
+ZFPOS(:,:,IKB) = (PSRC(:,:,IKB-1) - ZPHAT(:,:,IKB))*PCR(:,:,IKB) + &
                   ZPHAT(:,:,IKB)
 !
 ! below bottom flux - not used
@@ -2080,8 +2156,8 @@ ZFPOS(:,:,IKB-1) = 0.0
 ! negative fluxes:
 !
 ZFNEG(:,:,IKB:IKE) = ZPHAT(:,:,IKB:IKE) + & 
-     ZCR(:,:,IKB:IKE)*(ZPHAT(:,:,IKB:IKE) - PSRC(:,:,IKB:IKE)) + &
-     ZCR(:,:,IKB:IKE)*(1.0 + ZCR(:,:,IKB:IKE)) * &
+     PCR(:,:,IKB:IKE)*(ZPHAT(:,:,IKB:IKE) - PSRC(:,:,IKB:IKE)) + &
+     PCR(:,:,IKB:IKE)*(1.0 + PCR(:,:,IKB:IKE)) * &
      (ZPHAT(:,:,IKB:IKE) - 2.0*PSRC(:,:,IKB:IKE) +ZPHAT(:,:,IKB+1:IKE+1))
 !
 !!$   CALL GET_HALO(ZFNEG) ! JUAN
@@ -2090,14 +2166,14 @@ ZFNEG(:,:,IKB:IKE) = ZPHAT(:,:,IKB:IKE) + &
 ZFNEG(:,:,IKB-1) = 0.0 
 !
 ! negative flux at the TOP
-ZFNEG(:,:,IKE+1) = (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE+1))*ZCR(:,:,IKE+1) + &
+ZFNEG(:,:,IKE+1) = (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE+1))*PCR(:,:,IKE+1) + &
                     ZPHAT(:,:,IKE+1) 
 !
 ! calculate the advection
 !
 PR = PSRC * PRHO - &
-     DZF( ZCR*MZM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,ZCR)) + & 
-                          ZFNEG*(0.5-SIGN(0.5,ZCR)) ) )
+     DZF( PCR*MZM(PRHO)*( ZFPOS*(0.5+SIGN(0.5,PCR)) + & 
+                          ZFNEG*(0.5-SIGN(0.5,PCR)) ) )
 !
 ! in OPEN case fix boundary conditions
 !
@@ -2113,7 +2189,7 @@ END FUNCTION PPM_S0_Z
 !-------------------------------------------------------------------------------
 !
 !     ########################################################################
-      FUNCTION PPM_S1_X(HLBCX, KGRID, PSRC, ZCR, PRHO, PRHOT, &
+      FUNCTION PPM_S1_X(HLBCX, KGRID, PSRC, PCR, PRHO, PRHOT, &
                         PTSTEP) RESULT(PR)
 !     ########################################################################
 !!
@@ -2134,7 +2210,6 @@ USE MODI_SHUMAN
 USE MODD_CONF
 USE MODD_LUNIT
 USE MODD_PARAMETERS
-!USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
 !
 IMPLICIT NONE
 !
@@ -2145,14 +2220,14 @@ CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX  ! X direction LBC type
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHOT ! density at t+dt
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 !*       0.2   Declarations of local variables :
 !
@@ -2160,13 +2235,13 @@ INTEGER:: IIB,IJB,IKB    ! Begining useful area in x,y,z directions
 INTEGER:: IIE,IJE,IKE    ! End useful area in x,y,z directions
 !
 ! variable at cell edges
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZPHAT, ZRUT
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZPHAT, ZRUT
 !
 ! advection fluxes, upwind and correction
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFUP, ZFCOR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZFUP, ZFCOR
 !
 ! ratios for limiting the correction flux
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZRPOS, ZRNEG
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZRPOS, ZRNEG
 !
 ! variables for limiting the correction flux
 REAL :: ZSRCMAX, ZSRCMIN, ZFIN, ZFOUT
@@ -2189,7 +2264,7 @@ IKE = SIZE(PSRC,3) - JPVEXT
 !
 ! Calculate contravariant component rho*u/dx
 !
-ZRUT = ZCR/PTSTEP * MXM(PRHO)
+ZRUT = PCR/PTSTEP * MXM(PRHO)
 !
 ! calculate 4th order fluxes at cell edges in the inner domain
 !
@@ -2221,47 +2296,47 @@ END SELECT
 ! that makes it equivalent to the PPM flux
 ! flux_ppm = flux_up + flux_corr
 !
-WHERE ( ZCR(IIB:IIE,:,:) .GT. 0.0 )
+WHERE ( PCR(IIB:IIE,:,:) .GT. 0.0 )
    ZFUP(IIB:IIE,:,:) = ZRUT(IIB:IIE,:,:) * PSRC(IIB-1:IIE-1,:,:)
    ZFCOR(IIB:IIE,:,:) = ZRUT(IIB:IIE,:,:) * &
-        (1.0 - ZCR(IIB:IIE,:,:)) * &
-        (ZPHAT(IIB:IIE,:,:) - PSRC(IIB-1:IIE-1,:,:) - ZCR(IIB:IIE,:,:) * &
+        (1.0 - PCR(IIB:IIE,:,:)) * &
+        (ZPHAT(IIB:IIE,:,:) - PSRC(IIB-1:IIE-1,:,:) - PCR(IIB:IIE,:,:) * &
         (ZPHAT(IIB-1:IIE-1,:,:) - 2.0*PSRC(IIB-1:IIE-1,:,:)+ZPHAT(IIB:IIE,:,:)))
 ELSEWHERE
    ZFUP(IIB:IIE,:,:) = ZRUT(IIB:IIE,:,:) * PSRC(IIB:IIE,:,:)
    ZFCOR(IIB:IIE,:,:) = ZRUT(IIB:IIE,:,:) * &
-        (1.0 + ZCR(IIB:IIE,:,:)) * &
-        (ZPHAT(IIB:IIE,:,:) - PSRC(IIB:IIE,:,:) + ZCR(IIB:IIE,:,:) * &
+        (1.0 + PCR(IIB:IIE,:,:)) * &
+        (ZPHAT(IIB:IIE,:,:) - PSRC(IIB:IIE,:,:) + PCR(IIB:IIE,:,:) * &
         (ZPHAT(IIB:IIE,:,:) - 2.0*PSRC(IIB:IIE,:,:) + ZPHAT(IIB+1:IIE+1,:,:)))
 END WHERE
 !
 ! set boundaries to CYCL
 !
-WHERE ( ZCR(IIB-1,:,:) .GT. 0.0 )
+WHERE ( PCR(IIB-1,:,:) .GT. 0.0 )
    ZFUP(IIB-1,:,:) = ZRUT(IIB-1,:,:) * PSRC(IIE-1,:,:)
    ZFCOR(IIB-1,:,:) =  ZRUT(IIB-1,:,:) * &
-        (1.0 - ZCR(IIB-1,:,:)) * &
-        (ZPHAT(IIB-1,:,:) - PSRC(IIE-1,:,:) - ZCR(IIB-1,:,:) * &
+        (1.0 - PCR(IIB-1,:,:)) * &
+        (ZPHAT(IIB-1,:,:) - PSRC(IIE-1,:,:) - PCR(IIB-1,:,:) * &
         (ZPHAT(IIE-1,:,:) - 2.0*PSRC(IIE-1,:,:) + ZPHAT(IIB-1,:,:)))
 ELSEWHERE
    ZFUP(IIB-1,:,:) = ZRUT(IIB-1,:,:) * PSRC(IIB-1,:,:)
    ZFCOR(IIB-1,:,:) =  ZRUT(IIB-1,:,:) * &
-        (1.0 + ZCR(IIB-1,:,:)) * &
-        (ZPHAT(IIB-1,:,:) - PSRC(IIB-1,:,:) + ZCR(IIB-1,:,:) * &
+        (1.0 + PCR(IIB-1,:,:)) * &
+        (ZPHAT(IIB-1,:,:) - PSRC(IIB-1,:,:) + PCR(IIB-1,:,:) * &
         (ZPHAT(IIB-1,:,:) - 2.0*PSRC(IIB-1,:,:) + ZPHAT(IIB,:,:)))
 END WHERE
 !
-WHERE ( ZCR(IIE+1,:,:) .GT. 0.0 )
+WHERE ( PCR(IIE+1,:,:) .GT. 0.0 )
    ZFUP(IIE+1,:,:) = ZRUT(IIE+1,:,:) * PSRC(IIE,:,:)
    ZFCOR(IIE+1,:,:) =  ZRUT(IIE+1,:,:) * &
-        (1.0 - ZCR(IIE+1,:,:)) * &
-        (ZPHAT(IIE+1,:,:) - PSRC(IIE,:,:) - ZCR(IIE+1,:,:) * &
+        (1.0 - PCR(IIE+1,:,:)) * &
+        (ZPHAT(IIE+1,:,:) - PSRC(IIE,:,:) - PCR(IIE+1,:,:) * &
         (ZPHAT(IIE,:,:) - 2.0*PSRC(IIE,:,:) + ZPHAT(IIE+1,:,:)))
 ELSEWHERE
    ZFUP(IIE+1,:,:) = ZRUT(IIE+1,:,:) * PSRC(IIE+1,:,:)
    ZFCOR(IIE+1,:,:) =  ZRUT(IIE+1,:,:) * &
-        (1.0 + ZCR(IIE+1,:,:)) * &
-        (ZPHAT(IIE+1,:,:) - PSRC(IIE+1,:,:) + ZCR(IIE+1,:,:) * &
+        (1.0 + PCR(IIE+1,:,:)) * &
+        (ZPHAT(IIE+1,:,:) - PSRC(IIE+1,:,:) + PCR(IIE+1,:,:) * &
         (ZPHAT(IIE+1,:,:) - 2.0*PSRC(IIE+1,:,:) + ZPHAT(IIB+1,:,:)))
 END WHERE
 !
@@ -2367,7 +2442,7 @@ END FUNCTION PPM_S1_X
 !-------------------------------------------------------------------------------
 !
 !     ########################################################################
-      FUNCTION PPM_S1_Y(HLBCY, KGRID, PSRC, ZCR, PRHO, PRHOT, &
+      FUNCTION PPM_S1_Y(HLBCY, KGRID, PSRC, PCR, PRHO, PRHOT, &
                         PTSTEP) RESULT(PR)
 !     ########################################################################
 !!
@@ -2388,7 +2463,6 @@ USE MODI_SHUMAN
 USE MODD_LUNIT
 USE MODD_CONF
 USE MODD_PARAMETERS
-!USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
 !
 IMPLICIT NONE
 !
@@ -2399,14 +2473,14 @@ CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY  ! X direction LBC type
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHOT ! density at t+dt
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 !*       0.2   Declarations of local variables :
 !
@@ -2414,13 +2488,13 @@ INTEGER:: IIB,IJB,IKB   ! Begining useful area in x,y,z directions
 INTEGER:: IIE,IJE,IKE   ! End useful area in x,y,z directions
 !
 ! variable at cell edges
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZPHAT, ZRVT
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZPHAT, ZRVT
 !
 ! advection fluxes, upwind and correction
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFUP, ZFCOR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZFUP, ZFCOR
 !
 ! ratios for limiting the correction flux
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZRPOS, ZRNEG
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZRPOS, ZRNEG
 !
 ! variables for limiting the correction flux
 REAL :: ZSRCMAX, ZSRCMIN, ZFIN, ZFOUT
@@ -2449,7 +2523,7 @@ IKE = SIZE(PSRC,3) - JPVEXT
 !
 !-------------------------------------------------------------------------------
 !
-ZRVT = ZCR/PTSTEP * MYM(PRHO)
+ZRVT = PCR/PTSTEP * MYM(PRHO)
 !
 ! calculate 4th order fluxes at cell edges in the inner domain !
 ZPHAT(:,IJB+1:IJE,:) = (7.0 * &
@@ -2480,47 +2554,47 @@ END SELECT
 ! that makes it equivalent to the PPM flux
 ! flux_ppm = flux_up + flux_corr
 !
-WHERE ( ZCR(:,IJB:IJE,:) .GT. 0.0 )
+WHERE ( PCR(:,IJB:IJE,:) .GT. 0.0 )
    ZFUP(:,IJB:IJE,:) = ZRVT(:,IJB:IJE,:) * PSRC(:,IJB-1:IJE-1,:)
    ZFCOR(:,IJB:IJE,:) = ZRVT(:,IJB:IJE,:) * &
-        (1.0 - ZCR(:,IJB:IJE,:)) * &
-        (ZPHAT(:,IJB:IJE,:) - PSRC(:,IJB-1:IJE-1,:) - ZCR(:,IJB:IJE,:) * &
+        (1.0 - PCR(:,IJB:IJE,:)) * &
+        (ZPHAT(:,IJB:IJE,:) - PSRC(:,IJB-1:IJE-1,:) - PCR(:,IJB:IJE,:) * &
         (ZPHAT(:,IJB-1:IJE-1,:) - 2.0*PSRC(:,IJB-1:IJE-1,:)+ZPHAT(:,IJB:IJE,:)))
 ELSEWHERE
    ZFUP(:,IJB:IJE,:) = ZRVT(:,IJB:IJE,:) * PSRC(:,IJB:IJE,:)
    ZFCOR(:,IJB:IJE,:) = ZRVT(:,IJB:IJE,:) * &
-        (1.0 + ZCR(:,IJB:IJE,:)) * &
-        (ZPHAT(:,IJB:IJE,:) - PSRC(:,IJB:IJE,:) + ZCR(:,IJB:IJE,:) * &
+        (1.0 + PCR(:,IJB:IJE,:)) * &
+        (ZPHAT(:,IJB:IJE,:) - PSRC(:,IJB:IJE,:) + PCR(:,IJB:IJE,:) * &
         (ZPHAT(:,IJB:IJE,:) - 2.0*PSRC(:,IJB:IJE,:) + ZPHAT(:,IJB+1:IJE+1,:)))
 END WHERE
 !
 ! set boundaries to CYCL
 !
-WHERE ( ZCR(:,IJB-1,:) .GT. 0.0 )
+WHERE ( PCR(:,IJB-1,:) .GT. 0.0 )
    ZFUP(:,IJB-1,:) = ZRVT(:,IJB-1,:) * PSRC(:,IJE-1,:)
    ZFCOR(:,IJB-1,:) =  ZRVT(:,IJB-1,:) * &
-        (1.0 - ZCR(:,IJB-1,:)) * &
-        (ZPHAT(:,IJB-1,:) - PSRC(:,IJE-1,:) - ZCR(:,IJB-1,:) * &
+        (1.0 - PCR(:,IJB-1,:)) * &
+        (ZPHAT(:,IJB-1,:) - PSRC(:,IJE-1,:) - PCR(:,IJB-1,:) * &
         (ZPHAT(:,IJE-1,:) - 2.0*PSRC(:,IJE-1,:) + ZPHAT(:,IJB-1,:)))
 ELSEWHERE
    ZFUP(:,IJB-1,:) = ZRVT(:,IJB-1,:) * PSRC(:,IJB-1,:)
    ZFCOR(:,IJB-1,:) =  ZRVT(:,IJB-1,:) * &
-        (1.0 + ZCR(:,IJB-1,:)) * &
-        (ZPHAT(:,IJB-1,:) - PSRC(:,IJB-1,:) + ZCR(:,IJB-1,:) * &
+        (1.0 + PCR(:,IJB-1,:)) * &
+        (ZPHAT(:,IJB-1,:) - PSRC(:,IJB-1,:) + PCR(:,IJB-1,:) * &
         (ZPHAT(:,IJB-1,:) - 2.0*PSRC(:,IJB-1,:) + ZPHAT(:,IJB,:)))
 END WHERE
 !
-WHERE ( ZCR(:,IJE+1,:) .GT. 0.0 )
+WHERE ( PCR(:,IJE+1,:) .GT. 0.0 )
    ZFUP(:,IJE+1,:) = ZRVT(:,IJE+1,:) * PSRC(:,IJE,:)
    ZFCOR(:,IJE+1,:) =  ZRVT(:,IJE+1,:) * &
-        (1.0 - ZCR(:,IJE+1,:)) * &
-        (ZPHAT(:,IJE+1,:) - PSRC(:,IJE,:) - ZCR(:,IJE+1,:) * &
+        (1.0 - PCR(:,IJE+1,:)) * &
+        (ZPHAT(:,IJE+1,:) - PSRC(:,IJE,:) - PCR(:,IJE+1,:) * &
         (ZPHAT(:,IJE,:) - 2.0*PSRC(:,IJE,:) + ZPHAT(:,IJE+1,:)))
 ELSEWHERE
    ZFUP(:,IJE+1,:) = ZRVT(:,IJE+1,:) * PSRC(:,IJE+1,:)
    ZFCOR(:,IJE+1,:) =  ZRVT(:,IJE+1,:) * &
-        (1.0 + ZCR(:,IJE+1,:)) * &
-        (ZPHAT(:,IJE+1,:) - PSRC(:,IJE+1,:) + ZCR(:,IJE+1,:) * &
+        (1.0 + PCR(:,IJE+1,:)) * &
+        (ZPHAT(:,IJE+1,:) - PSRC(:,IJE+1,:) + PCR(:,IJE+1,:) * &
         (ZPHAT(:,IJE+1,:) - 2.0*PSRC(:,IJE+1,:) + ZPHAT(:,IJB+1,:)))
 END WHERE
 !
@@ -2624,7 +2698,7 @@ END FUNCTION PPM_S1_Y
 !-------------------------------------------------------------------------------
 !
 !     ########################################################################
-      FUNCTION PPM_S1_Z(KGRID, PSRC, ZCR, PRHO, PRHOT, PTSTEP) &
+      FUNCTION PPM_S1_Z(KGRID, PSRC, PCR, PRHO, PRHOT, PTSTEP) &
                         RESULT(PR)
 !     ########################################################################
 !!
@@ -2643,7 +2717,6 @@ USE MODI_SHUMAN
 !
 USE MODD_CONF
 USE MODD_PARAMETERS
-!USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
 !
 IMPLICIT NONE
 !
@@ -2652,14 +2725,14 @@ IMPLICIT NONE
 INTEGER,                INTENT(IN)  :: KGRID   ! C grid localisation
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PSRC    ! variable at t
-REAL, DIMENSION(:,:,:), INTENT(IN)  :: ZCR     ! Courant number
+REAL, DIMENSION(:,:,:), INTENT(IN)  :: PCR     ! Courant number
 !
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHO  ! density
 REAL, DIMENSION(:,:,:), INTENT(IN)  :: PRHOT ! density at t+dt
 REAL,                   INTENT(IN)  :: PTSTEP  ! Time step 
 !
 ! output source term
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: PR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: PR
 !
 !*       0.2   Declarations of local variables :
 !
@@ -2667,13 +2740,13 @@ INTEGER:: IIB,IJB,IKB   ! Begining useful area in x,y,z directions
 INTEGER:: IIE,IJE,IKE   ! End useful area in x,y,z directions
 !
 ! variable at cell edges
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZPHAT, ZRVT
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZPHAT, ZRVT
 !
 ! advection fluxes, upwind and correction
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZFUP, ZFCOR
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZFUP, ZFCOR
 !
 ! ratios for limiting the correction flux
-REAL, DIMENSION(SIZE(ZCR,1),SIZE(ZCR,2),SIZE(ZCR,3)) :: ZRPOS, ZRNEG
+REAL, DIMENSION(SIZE(PCR,1),SIZE(PCR,2),SIZE(PCR,3)) :: ZRPOS, ZRNEG
 !
 ! variables for limiting the correction flux
 REAL :: ZSRCMAX, ZSRCMIN, ZFIN, ZFOUT
@@ -2693,7 +2766,7 @@ IKE = SIZE(PSRC,3) - JPVEXT
 !
 !-------------------------------------------------------------------------------
 !
-ZRVT = ZCR/PTSTEP * MZM(PRHO)
+ZRVT = PCR/PTSTEP * MZM(PRHO)
 !
 ! calculate 4th order fluxes at cell edges in the inner domain !
 ZPHAT(:,:,IKB+1:IKE) = (7.0 * &
@@ -2721,78 +2794,78 @@ ZPHAT(:,:,IKE+1) = (7.0 * &
 ! that makes it equivalent to the PPM flux
 ! flux_ppm = flux_up + flux_corr
 !
-WHERE ( ZCR(:,:,IKB:IKE) .GT. 0.0 )
+WHERE ( PCR(:,:,IKB:IKE) .GT. 0.0 )
    ZFUP(:,:,IKB:IKE) = ZRVT(:,:,IKB:IKE) * PSRC(:,:,IKB-1:IKE-1)
    ZFCOR(:,:,IKB:IKE) = ZRVT(:,:,IKB:IKE) * &
-        (1.0 - ZCR(:,:,IKB:IKE)) * &
-        (ZPHAT(:,:,IKB:IKE) - PSRC(:,:,IKB-1:IKE-1) - ZCR(:,:,IKB:IKE) * &
+        (1.0 - PCR(:,:,IKB:IKE)) * &
+        (ZPHAT(:,:,IKB:IKE) - PSRC(:,:,IKB-1:IKE-1) - PCR(:,:,IKB:IKE) * &
         (ZPHAT(:,:,IKB-1:IKE-1) - 2.0*PSRC(:,:,IKB-1:IKE-1)+ZPHAT(:,:,IKB:IKE)))
 ELSEWHERE
    ZFUP(:,:,IKB:IKE) = ZRVT(:,:,IKB:IKE) * PSRC(:,:,IKB:IKE)
    ZFCOR(:,:,IKB:IKE) = ZRVT(:,:,IKB:IKE) * &
-        (1.0 + ZCR(:,:,IKB:IKE)) * &
-        (ZPHAT(:,:,IKB:IKE) - PSRC(:,:,IKB:IKE) + ZCR(:,:,IKB:IKE) * &
+        (1.0 + PCR(:,:,IKB:IKE)) * &
+        (ZPHAT(:,:,IKB:IKE) - PSRC(:,:,IKB:IKE) + PCR(:,:,IKB:IKE) * &
         (ZPHAT(:,:,IKB:IKE) - 2.0*PSRC(:,:,IKB:IKE) + ZPHAT(:,:,IKB+1:IKE+1)))
 END WHERE
 !
 ! set BC to WALL
 !
-WHERE ( ZCR(:,:,IKB-1) .GT. 0.0 )
+WHERE ( PCR(:,:,IKB-1) .GT. 0.0 )
    ZFUP(:,:,IKB-1) = ZRVT(:,:,IKB-1) * PSRC(:,:,IKB+2)
    ZFCOR(:,:,IKB-1) =  ZRVT(:,:,IKB-1) * &
-        (1.0 - ZCR(:,:,IKB-1)) * &
-        (ZPHAT(:,:,IKB+1) - PSRC(:,:,IKB+2) - ZCR(:,:,IKB+1) * &
+        (1.0 - PCR(:,:,IKB-1)) * &
+        (ZPHAT(:,:,IKB+1) - PSRC(:,:,IKB+2) - PCR(:,:,IKB+1) * &
         (ZPHAT(:,:,IKB+2) - 2.0*PSRC(:,:,IKB+2) + ZPHAT(:,:,IKB+1)))
 ELSEWHERE
    ZFUP(:,:,IKB-1) = ZRVT(:,:,IKB-1) * PSRC(:,:,IKB+1)
    ZFCOR(:,:,IKB-1) =  ZRVT(:,:,IKB-1) * &
-        (1.0 + ZCR(:,:,IKB-1)) * &
-        (ZPHAT(:,:,IKB+1) - PSRC(:,:,IKB+1) + ZCR(:,:,IKB+1) * &
+        (1.0 + PCR(:,:,IKB-1)) * &
+        (ZPHAT(:,:,IKB+1) - PSRC(:,:,IKB+1) + PCR(:,:,IKB+1) * &
         (ZPHAT(:,:,IKB+1) - 2.0*PSRC(:,:,IKB+1) + ZPHAT(:,:,IKB)))
 END WHERE
 !
-WHERE ( ZCR(:,:,IKE+1) .GT. 0.0 )
+WHERE ( PCR(:,:,IKE+1) .GT. 0.0 )
    ZFUP(:,:,IKE+1) = ZRVT(:,:,IKE+1) * PSRC(:,:,IKE)
    ZFCOR(:,:,IKE+1) =  ZRVT(:,:,IKE+1) * &
-        (1.0 - ZCR(:,:,IKE+1)) * &
-        (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE) - ZCR(:,:,IKE+1) * &
+        (1.0 - PCR(:,:,IKE+1)) * &
+        (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE) - PCR(:,:,IKE+1) * &
         (ZPHAT(:,:,IKE) - 2.0*PSRC(:,:,IKE) + ZPHAT(:,:,IKE+1)))
 ELSEWHERE
    ZFUP(:,:,IKE+1) = ZRVT(:,:,IKE+1) * PSRC(:,:,IKE+1)
    ZFCOR(:,:,IKE+1) =  ZRVT(:,:,IKE+1) * &
-        (1.0 + ZCR(:,:,IKE+1)) * &
-        (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE+1) + ZCR(:,:,IKE+1) * &
+        (1.0 + PCR(:,:,IKE+1)) * &
+        (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE+1) + PCR(:,:,IKE+1) * &
         (ZPHAT(:,:,IKE+1) - 2.0*PSRC(:,:,IKE+1) + ZPHAT(:,:,IKE)))
 END WHERE
 !
 !
 !!$! set boundaries to CYCL
 !!$!
-!!$WHERE ( ZCR(:,:,IKB-1) .GT. 0.0 )
+!!$WHERE ( PCR(:,:,IKB-1) .GT. 0.0 )
 !!$   ZFUP(:,:,IKB-1) = ZRVT(:,:,IKB-1) * PSRC(:,:,IKE-1)
 !!$   ZFCOR(:,:,IKB-1) =  ZRVT(:,:,IKB-1) * &
-!!$        (1.0 - ZCR(:,:,IKB-1)) * &
-!!$        (ZPHAT(:,:,IKB-1) - PSRC(:,:,IKE-1) - ZCR(:,:,IKB-1) * &
+!!$        (1.0 - PCR(:,:,IKB-1)) * &
+!!$        (ZPHAT(:,:,IKB-1) - PSRC(:,:,IKE-1) - PCR(:,:,IKB-1) * &
 !!$        (ZPHAT(:,:,IKE-1) - 2.0*PSRC(:,:,IKE-1) + ZPHAT(:,:,IKB-1)))
 !!$ELSEWHERE
 !!$   ZFUP(:,:,IKB-1) = ZRVT(:,:,IKB-1) * PSRC(:,:,IKB-1)
 !!$   ZFCOR(:,:,IKB-1) =  ZRVT(:,:,IKB-1) * &
-!!$        (1.0 + ZCR(:,:,IKB-1)) * &
-!!$        (ZPHAT(:,:,IKB-1) - PSRC(:,:,IKB-1) + ZCR(:,:,IKB-1) * &
+!!$        (1.0 + PCR(:,:,IKB-1)) * &
+!!$        (ZPHAT(:,:,IKB-1) - PSRC(:,:,IKB-1) + PCR(:,:,IKB-1) * &
 !!$        (ZPHAT(:,:,IKB-1) - 2.0*PSRC(:,:,IKB-1) + ZPHAT(:,:,IKB)))
 !!$END WHERE
 !!$!
-!!$WHERE ( ZCR(:,:,IKE+1) .GT. 0.0 )
+!!$WHERE ( PCR(:,:,IKE+1) .GT. 0.0 )
 !!$   ZFUP(:,:,IKE+1) = ZRVT(:,:,IKE+1) * PSRC(:,:,IKE)
 !!$   ZFCOR(:,:,IKE+1) =  ZRVT(:,:,IKE+1) * &
-!!$        (1.0 - ZCR(:,:,IKE+1)) * &
-!!$        (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE) - ZCR(:,:,IKE+1) * &
+!!$        (1.0 - PCR(:,:,IKE+1)) * &
+!!$        (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE) - PCR(:,:,IKE+1) * &
 !!$        (ZPHAT(:,:,IKE) - 2.0*PSRC(:,:,IKE) + ZPHAT(:,:,IKE+1)))
 !!$ELSEWHERE
 !!$   ZFUP(:,:,IKE+1) = ZRVT(:,:,IKE+1) * PSRC(:,:,IKE+1)
 !!$   ZFCOR(:,:,IKE+1) =  ZRVT(:,:,IKE+1) * &
-!!$        (1.0 + ZCR(:,:,IKE+1)) * &
-!!$        (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE+1) + ZCR(:,:,IKE+1) * &
+!!$        (1.0 + PCR(:,:,IKE+1)) * &
+!!$        (ZPHAT(:,:,IKE+1) - PSRC(:,:,IKE+1) + PCR(:,:,IKE+1) * &
 !!$        (ZPHAT(:,:,IKE+1) - 2.0*PSRC(:,:,IKE+1) + ZPHAT(:,:,IKB+1)))
 !!$END WHERE
 !