Philippe 27/03/2023: pressurez: merge ZSOLVER version into MNH/

src/MNH/pressurez.f90

@@ -288,7 +288,11 @@ USE MODI_SHUMAN
 #ifdef MNH_OPENACC
 USE MODI_SHUMAN_DEVICE
 #endif
-!
+#ifdef MNH_MGSOLVER
+USE MODI_ZCONJGRAD
+USE MODI_ZSOLVER
+#endif
+
 IMPLICIT NONE
 !
 !*       0.1   declarations of arguments
@@ -454,21 +458,21 @@ IF (MPPDB_INITIALIZED) THEN
   CALL MPPDB_CHECK(PDZZ,"PRESSUREZ beg:PDZZ")
   CALL MPPDB_CHECK(PDZX,"PRESSUREZ beg:PDZX")
   CALL MPPDB_CHECK(PRHOT,"PRESSUREZ beg:PRHOT")
-  CALL MPPDB_CHECK(PAF,"PRESSUREZ beg:PAF")
-  CALL MPPDB_CHECK(PCF,"PRESSUREZ beg:PCF")
-  CALL MPPDB_CHECK(PBF,"PRESSUREZ beg:PBF")
-  !CALL MPPDB_CHECK(PTRIGSX,"PRESSUREZ beg:PTRIGSX")
-  !CALL MPPDB_CHECK(PTRIGSY,"PRESSUREZ beg:PTRIGSY")
-  !CALL MPPDB_CHECK(KIFAXX,"PRESSUREZ beg:KIFAXX")
-  !CALL MPPDB_CHECK(KIFAXY,"PRESSUREZ beg:KIFAXY")
+!!$  CALL MPPDB_CHECK(PAF,"PRESSUREZ beg:PAF")
+!!$  CALL MPPDB_CHECK(PCF,"PRESSUREZ beg:PCF")
+!!$  CALL MPPDB_CHECK(PBF,"PRESSUREZ beg:PBF")
+!!$  CALL MPPDB_CHECK(PTRIGSX,"PRESSUREZ beg:PTRIGSX")
+!!$  CALL MPPDB_CHECK(PTRIGSY,"PRESSUREZ beg:PTRIGSY")
+!!$  CALL MPPDB_CHECK(KIFAXX,"PRESSUREZ beg:KIFAXX")
+!!$  CALL MPPDB_CHECK(KIFAXY,"PRESSUREZ beg:KIFAXY")
   CALL MPPDB_CHECK(PTHT,"PRESSUREZ beg:PTHT")
   CALL MPPDB_CHECK(PRT,"PRESSUREZ beg:PRT")
   CALL MPPDB_CHECK(PRHODREF,"PRESSUREZ beg:PRHODREF")
   CALL MPPDB_CHECK(PTHVREF,"PRESSUREZ beg:PTHVREF")
   CALL MPPDB_CHECK(PRVREF,"PRESSUREZ beg:PRVREF")
   CALL MPPDB_CHECK(PEXNREF,"PRESSUREZ beg:PEXNREF")
-  !CALL MPPDB_CHECK(PBFB,"PRESSUREZ beg:PBFB")
-  CALL MPPDB_CHECK(PBF_SXP2_YP1_Z,"PRESSUREZ beg:PBF_SXP2_YP1_Z")
+!!$  CALL MPPDB_CHECK(PBFB,"PRESSUREZ beg:PBFB")
+  !/!\ don't work not of size NX/NY/NZ CALL MPPDB_CHECK(PBF_SXP2_YP1_Z,"PRESSUREZ beg:PBF_SXP2_YP1_Z")
   !Check all INOUT arrays
   CALL MPPDB_CHECK(PRUS,"PRESSUREZ beg:PRUS")
   CALL MPPDB_CHECK(PRVS,"PRESSUREZ beg:PRVS")
@@ -487,8 +491,13 @@ END IF
 ILUOUT = TLUOUT%NLU
 !
 CALL GET_GLOBALDIMS_ll (IIMAX_ll,IJMAX_ll)
+#ifndef MNH_MGSOLVER
 IF ( ( MIN(IIMAX_ll,IJMAX_ll) < NPROC  ) .AND. ( HPRESOPT /= 'ZRESI' ) ) THEN
-   WRITE(UNIT=ILUOUT,FMT=*) 'ERROR IN PRESSUREZ:: YOU WANT TO USE TO MANY PROCESSES WITHOUT CPRESOPT="ZRESI" '
+   WRITE(UNIT=ILUOUT,FMT=*) 'ERROR IN PRESSUREZ:: YOU WANT TO USE TOO MANY PROCESSES WITHOUT CPRESOPT="ZRESI" '
+#else
+IF ( ( MIN(IIMAX_ll,IJMAX_ll) < NPROC  ) .AND. ( HPRESOPT /= 'ZRESI' ) .AND. ( HPRESOPT /= 'ZSOLV' )) THEN
+   WRITE(UNIT=ILUOUT,FMT=*) 'ERROR IN PRESSUREZ:: YOU WANT TO USE TOO MANY PROCESSES WITHOUT CPRESOPT="ZRESI/ZSOLV" '
+#endif
    WRITE(UNIT=ILUOUT,FMT=*) 'MIN(IIMAX_ll,IJMAX_ll)=',MIN(IIMAX_ll,IJMAX_ll),' < NPROC =', NPROC
    WRITE(UNIT=ILUOUT,FMT=*) 'YOU HAVE TO SET CPRESOPT="ZRESI => JOB ABORTED '
    CALL PRINT_MSG(NVERB_FATAL,'GEN','PRESSUREZ','')
@@ -735,12 +744,28 @@ ELSE
      CALL CONJGRAD(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,ZTHETAV,       &
      PDXHATM,PDYHATM,PRHOT,PAF,PBF,PCF,PTRIGSX,PTRIGSY,                       &
      KIFAXX,KIFAXY,KITR,ZDV_SOURCE,ZPHIT)
+#ifdef MNH_MGSOLVER
+  CASE('ZGRAD')     ! Conjugate Gradient method
+     CALL ZCONJGRAD(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,ZTHETAV,       &
+     PDXHATM,PDYHATM,PRHOT,PAF,PBF,PCF,PTRIGSX,PTRIGSY,                       &
+     KIFAXX,KIFAXY,KITR,ZDV_SOURCE,ZPHIT)
+#endif
 !
   CASE('CRESI')     ! Conjugate Residual method
     CALL CONRESOL(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,ZTHETAV,       &
     PDXHATM,PDYHATM,PRHOT,PAF,PBF,PCF,PTRIGSX,PTRIGSY,                       &
     KIFAXX,KIFAXY,KITR,ZDV_SOURCE,ZPHIT)
 !
+#ifdef MNH_MGSOLVER
+  CASE('ZSOLV')     ! Conjugate Residual method
+    CALL ZSOLVER(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,ZTHETAV,       &
+    PDXHATM,PDYHATM,PRHOT,PAF,PBF,PCF,PTRIGSX,PTRIGSY,                      &
+    KIFAXX,KIFAXY,KITR,ZDV_SOURCE,ZPHIT,                                    &
+    PAF_ZS,PBF_ZS,PCF_ZS,                                                   &
+    PDXATH_ZS,PDYATH_ZS,PRHO_ZS,PBFB,                                       &
+    A_K,B_K,C_K,D_K)
+!
+#endif
   CASE('ZRESI')     ! Conjugate Residual method
     CALL CONRESOLZ(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,ZTHETAV,      &
     PDXHATM,PDYHATM,PRHOT,PAF,PBF,PCF,PTRIGSX,PTRIGSY,                       &
@@ -757,26 +782,51 @@ END IF
 !
 IF ( GWEST ) THEN
    !$acc kernels async
+#if 0
+!PW: original version
     ZPHIT(IIB-1,:,IKB-1) = ZPHIT(IIB,:,IKB)
     ZPHIT(IIB-1,:,IKE+1) = ZPHIT(IIB,:,IKE)
+#else
+!PW: version from ZSOLVER
+   ZPHIT(IIB-1,:,IKB-1) = ZPHIT(IIB,:,IKB-1)
+#endif
    !$acc end kernels
 END IF
 IF ( GEAST ) THEN
    !$acc kernels async
+#if 0
+!PW: original version
     ZPHIT(IIE+1,:,IKB-1) = ZPHIT(IIE,:,IKB)
     ZPHIT(IIE+1,:,IKE+1) = ZPHIT(IIE,:,IKE)
+#else
+!PW: version from ZSOLVER
+   ZPHIT(IIE+1,:,IKB-1) = ZPHIT(IIE,:,IKB-1)
+#endif
    !$acc end kernels
 END IF
 !
 IF ( GSOUTH ) THEN
    !$acc kernels async
+#if 0
+!PW: original version
     ZPHIT(:,IJB-1,IKB-1) = ZPHIT(:,IJB,IKB)
     ZPHIT(:,IJB-1,IKE+1) = ZPHIT(:,IJB,IKE)    
+#else
+!PW: version from ZSOLVER
+   ZPHIT(:,IJB-1,IKB-1) = ZPHIT(:,IJB,IKB-1)
+#endif
    !$acc end kernels
 END IF
 IF ( GNORTH ) THEN
    !$acc kernels async
+#if 0
+!PW: original version
+    ZPHIT(:,IJE+1,IKB-1) = ZPHIT(:,IJE,IKB)
+    ZPHIT(:,IJE+1,IKE+1) = ZPHIT(:,IJE,IKE)  
+#else
+!PW: version from ZSOLVER
    ZPHIT(:,IJE+1,IKB-1) = ZPHIT(:,IJE,IKB-1)
+#endif
    !$acc end kernels
 END IF
 IF ( GSOUTH2D ) THEN