C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| C PROGRAM NEST_DRIVER C !DESCRIPTION: C Routine that manages the MPI communication between the CHILD C and PARENT models. It performs also the necessary C interpolations from PARENT2CHILD and CHILD2PARENT. C C ver 1.0 by G. Sannino, V. Ruggiero, A. Carillo, P. Heimbach C C First application described in: C Sannino G.,Herrmann, Carillo, Rupolo, Ruggiero, Artale, Heimbach, 2009: C An eddy-permitting model of the Mediterranean Sea with a two-way grid C refinement at Gibraltar. C Ocean Modelling, 30(1), 56-72, doi: 10.1016/j.ocemod.2009.06.002 C C !LOCAL INPUT VARIABLES: c --------------------------------------------------------------------------------- c NST_LEV_TOT :: Total nesting levels (1 for only one nesting) c NST_LEV :: Number of the actual nesting c NCPUs_CHLD :: Number of CPUs used for the CHILD model c at NST_LEV nesting level c NCPUs_PRNT :: Number of CPUs used for the PARENT model c at NST_LEV nesting level c nSxC,nSyF :: Domain decomposition used for CHILD c nSxP,nSyP :: Domain decomposition used for PARENT c OLX,OLY :: Domain dec. overlapping (same for both models) c NrC,NyC,NxC :: Dimension PARENT model c NrF,NyF,NxF :: Dimension CHILD model c WesterB :: Western (i) PARENT index where start the refinement c EasterB :: Eastern (i) PARENT index where finish the refinement c dirNEST :: Directory where are stored the geometry data of both models c --------------------------------------------------------------------------------- CEOP implicit none c-------------------------------------------------------- c INPUT VARIABLE DEFINITION c-------------------------------------------------------- INTEGER :: NST_LEV_TOT, NST_LEV, NCPUs_PRNT INTEGER :: Count_Lev PARAMETER (NST_LEV_TOT = 1) !Number of Total Nesting Levels PARAMETER (NST_LEV = 1) !Which level am I? INTEGER :: NCPUs_CHLD(NST_LEV_TOT) INTEGER :: MSTR_DRV(NST_LEV_TOT) INTEGER :: MSTR_PRNT(NST_LEV_TOT) INTEGER :: MSTR_CHLD(NST_LEV_TOT) PARAMETER (NCPUs_PRNT = 40) DATA NCPUs_CHLD / 20 / c-------------------------------------------------------- INTEGER :: nSxC,nSyC !Domain decomposition CHILD INTEGER :: nSxP,nSyP !Domain decomposition PARENT PARAMETER (nSxC = 4 , nSyC = 5) PARAMETER (nSxP = 8 , nSyP = 5) c-------------------------------------------------------- INTEGER :: OLY,OLX !Domain decomposition overlapping c !(the same for both models) PARAMETER (OLX = 3, OLY = 3) c-------------------------------------------------------- INTEGER :: NrP,NxP,NyP INTEGER :: NrC,NxC,NyC INTEGER :: IM_C,JM_C INTEGER :: IM_P,JM_P INTEGER :: IndI,IndJ INTEGER :: IndI_P(nSxP*nSyP),IndJ_P(nSxP*nSyP) INTEGER :: IndI_C(nSxC*nSyC),IndJ_C(nSxC*nSyC) INTEGER :: WesternB,EasternB c-------------------------------------------------------- PARAMETER (NrP=42, NyP=120,NxP = 400) !PARENT model PARAMETER (NrC=42, NyC=105,NxC = 140) !CHILD model c-------------------------------------------------------- PARAMETER (WesternB = 43,EasternB=90) c-------------------------------------------------------- CHARACTER :: dirNEST*80 c-------------------------------------------------------- PARAMETER (dirNEST ="/home/sannino/NESTING/") c-------------------------------------------------------- INCLUDE 'mpif.h' INTEGER :: ierr,rank,size,npd INTEGER :: irank,isize INTEGER :: color INTEGER :: istatus,NEST_comm INTEGER :: from,whm,status(MPI_STATUS_SIZE),st_count INTEGER :: I,J,K,II,JJ,Irec,III,JJJ,KK,ICONT INTEGER :: iVar,Indx,Jndx INTEGER :: J1,J2,JJ1,JJ2 INTEGER :: I_START,I_END,I_STEP REAL*4 :: XF,YF,XP1,YP1,XP2,YP2,YP3 REAL*8 :: TRANSPORT_WEST,TRANSPORT_EAST CHARACTER*10 :: c2i(30) c---------------------------------------------------- c Define PARENT Model Geometry c---------------------------------------------------- c REAL*4 Xu_P(NxP,NyP) REAL*4 :: Yu_P(NxP,NyP) REAL*4 :: Xv_P(NxP,NyP) REAL*4 :: Yv_P(NxP,NyP) REAL*4 :: Xo_P(NxP,NyP) REAL*4 :: Yo_P(NxP,NyP) REAL*4 :: Xg_P(NxP,NyP) REAL*4 :: Yg_P(NxP,NyP) REAL*4 :: hFacW_P(NxP,NyP,NrP) REAL*4 :: hFacS_P(NxP,NyP,NrP) REAL*4 :: RAC_P(NxP,NyP) REAL*4 :: RAW_P(NxP,NyP) REAL*4 :: RAS_P(NxP,NyP) REAL*4 :: hFacC_P(NxP,NyP,NrP) REAL*4 :: DEEP_P(NxP,NyP,NrP) REAL*4 :: INV_VOL_C_P(NxP,NyP,NrP) REAL*4 :: INV_VOL_S_P(NxP,NyP,NrP) REAL*4 :: INV_VOL_W_P(NxP,NyP,NrP) c---------------------------------------------------- c Define CHILD Model Geometry c---------------------------------------------------- REAL*4 :: Xu_C(NxC,NyC) REAL*4 :: Yu_C(NxC,NyC) REAL*4 :: Xv_C(NxC,NyC) REAL*4 :: Yv_C(NxC,NyC) REAL*4 :: Xo_C(NxC,NyC) REAL*4 :: Yo_C(NxC,NyC) REAL*4 :: Xg_C(NxC,NyC) REAL*4 :: Yg_C(NxC,NyC) REAL*4 :: hFacW_C(NxC,NyC,NrC) REAL*4 :: hFacS_C(NxC,NyC,NrC) REAL*4 :: RAC_C(NxC,NyC) REAL*4 :: RAW_C(NxC,NyC) REAL*4 :: RAS_C(NxC,NyC) REAL*4 :: hFacC_C(NxC,NyC,NrC) REAL*4 :: DEEP_C(NxC,NyC,NrC) c---------------------------------------------------- c Define relative (PARENT-->CHILD) indicies c---------------------------------------------------- INTEGER :: P2C_U(NyC) c INTEGER :: P2C_linU(NyC) !Linear interp. INTEGER :: WO3_linU(NyC) !Linear interp. !Which Of 3 c INTEGER :: P2C_linV(NyC) !Linear interp. INTEGER :: WO3_linV(NyC) !Linear interp. !Which Of 3 c INTEGER :: P2C_V(NyC) !Linear interp. INTEGER :: P2C_o(NyC) !Linear interp. INTEGER :: P2C1_V(NyC) !BiLinear interp. INTEGER :: P2C2_V(NyC) !BiLinear interp. INTEGER :: P2C1_o(NyC) !BiLinear interp. INTEGER :: P2C2_o(NyC) !BiLinear interp. c---------------------------------------------------- c Define relative (CHILD-->PARENT) indicies c---------------------------------------------------- INTEGER I_C2P(9,NxP,NyP) INTEGER J_C2P(9,NxP,NyP) c---------------------------------------------------- c Define CHILD model variable c---------------------------------------------------- c _____________ (1) WesternB (2) EasternB c | REAL*8 :: U_C1(NyC,NrC,2) REAL*8 :: V_C1(NyC,NrC,2) REAL*8 :: T_C1(NyC,NrC,2) REAL*8 :: S_C1(NyC,NrC,2) REAL*8 :: ETA_C1(NyC,NrC,2) INTEGER :: MSIZE c REAL*8 :: U_C2(NyC,NrC,2) REAL*8 :: V_C2(NyC,NrC,2) REAL*8 :: T_C2(NyC,NrC,2) REAL*8 :: S_C2(NyC,NrC,2) REAL*8 :: ETA_C2(NyC,NrC,2) c REAL*8,allocatable :: VAR_C1(:,:,:,:) c REAL*8 :: DIFF_U(NyC,NrC,2) REAL*8 :: DIFF_V(NyC,NrC,2) REAL*8 :: DIFF_T(NyC,NrC,2) REAL*8 :: DIFF_S(NyC,NrC,2) REAL*8 :: DIFF_ETA(NyC,NrC,2) c---------------------------------------------------- c Define PARENT model variable c---------------------------------------------------- REAL*8 :: VAR3D_P(NxP,NyP,NrP,15) REAL*8 :: VAR2D_P(NxP,NyP,4) INTEGER :: ONOFF INTEGER :: index_var3D,index_var2D c---------------------------------------------------------------| c (1) U || (2) V || (3) T || (4) S | c---------------------------------------------------------------| c (5) gU || (6) gV || (7) gT || (8) gS | c---------------------------------------------------------------| c (9) gUNm1 || (10) gVNm1 || (11) gTNm1 || (12) gSNm1 | c---------------------------------------------------------------| c (13) totPhiHyd || (14) IVDConvCount || | c---------------------------------------------------------------| c (15) etaN || (16) etaH || (17) phiHydLow | c---------------------------------------------------------------| c (18) etaNm1 || (19) etaHm1|| | c---------------------------------------------------------------| c---------------------------------------------------------------| c Define Global Variables to Exchange | c---------------------------------------------------------------| REAL*8,allocatable :: globalPA (:,:,:,:) !(6,NyP,NrP,5) REAL*8 :: globalP1(6,NyP,NrP) REAL*8 :: globalP2(6,NyP,NrP) REAL*8 :: globalP3(6,NyP,NrP) REAL*8 :: globalP4(6,NyP,NrP) REAL*8 :: globalP5(6,NyP,NrP) REAL*8 :: globalP6(6,NyP,NrP) REAL*8 :: globalP7(6,NyP,NrP) REAL*8 :: globalP8(6,NyP,NrP) REAL*8 :: globalP9(6,NyP,NrP) REAL*8 :: globalP10(6,NyP,NrP) REAL*8 :: globalP11(6,NyP,NrP) REAL*8 :: globalP12(6,NyP,NrP) REAL*8 :: globalP13(6,NyP,NrP) REAL*8 :: globalP14(6,NyP,NrP) c INTEGER :: index c---------------------------------------------------- c Define Global Variables to Exchange c---------------------------------------------------- REAL*8 :: globalC3D(NxC,NyC,NrC,15) c |___________ 15 fields c REAL*8 globalC2D(NxC,NyC,4) c |___________ 4 fields c c REAL*8,allocatable :: globalC3D_a(:,:,:,:),globalC2D_a(:,:,:) c INTEGER :: indexF,index1F REAL*4 :: AREA_VOL INTEGER :: vstart,vstop,VCONT,VCONTP(0:3) c---------------------------------------------------- c MPI starts here c---------------------------------------------------- call MPI_Init(ierr) call MPI_Comm_size(MPI_COMM_WORLD,size,ierr) call MPI_Comm_rank(MPI_COMM_WORLD,rank,ierr) npd=size-(NCPUs_PRNT+NCPUs_CHLD(1)) if(rank.lt.npd) color=0 call MPI_COMM_SPLIT (MPI_COMM_WORLD, color,0, & NEST_COMM,ierr) c call MPI_Comm_size(NEST_COMM,isize,ierr) call MPI_Comm_rank(NEST_COMM,irank,ierr) c-------------------------------------------------------- c COMPUTE MASTER VALUES c-------------------------------------------------------- MSTR_DRV(1) = 0 c MSTR_PRNT(1) = npd MSTR_CHLD(1) = NCPUs_PRNT + npd c DO Count_Lev = 2, NST_LEV_TOT MSTR_DRV(Count_Lev) = MSTR_CHLD(Count_Lev-1) + & NCPUs_CHLD(Count_Lev - 1) c MSTR_CHLD(Count_Lev) = MSTR_DRV(Count_Lev) + 1 MSTR_PRNT(Count_Lev) = MSTR_CHLD(Count_Lev-1) ENDDO c vstart = 1+rank*(nSxP/MSTR_PRNT(1)) vstop = (1+rank)*(nSxP/MSTR_PRNT(1)) VCONT = (nSxP/npd)*(nSyP)*rank-1 VCONTP(rank) = VCONT c-------------------------------------------------------- c COMPUTE DOMAIN DECOMPOSITION c-------------------------------------------------------- if(rank.eq.0) then c IM_C = int(NxC/nSxC) JM_C = int(NyC/nSyC) c IM_P = int(NxP/nSxP) JM_P = int(NyP/nSyP) c indexF = IM_C*JM_C*NrC*15 index1F = IM_C*JM_C*4 c indexF = (JM_C+OLY+OLY)*NrC*2*5 c index_var3D = IM_P*JM_P*NrP index_var2D = IM_P*JM_P c ICONT = 0 DO I = 1,nSxP DO J = 1,nSyP ICONT = ICONT + 1 IndI_P(ICONT) = IM_P*(I-1) IndJ_P(ICONT) = JM_P*(J-1) END DO END DO c ICONT = 0 DO I = 1,nSxC DO J = 1,nSyC ICONT = ICONT + 1 IndI_C(ICONT) = IM_C*(I-1) IndJ_C(ICONT) = JM_C*(J-1) END DO END DO c ALLOCATE (globalPA (6,JM_P,NrP,5)) index=6*JM_P*NrP*5 c ALLOCATE(globalC3D_a(IM_C,JM_C,NrC,15)) ALLOCATE(globalC2D_a(IM_C,JM_C,4)) c ALLOCATE (VAR_C1((JM_C+OLY+OLY),NrC,2,5)) c-------------------------------------------------------- c WARNING c-------------------------------------------------------- print*,'*************************************' print*,' have you update geometric files?' print*,' in ./CHILD e ./PARENT' print*,'*************************************' c-------------------------------------------------------- c PARENT MODEL c-------------------------------------------------------- print*,' [1] Read PARENT model geometry' c---------------------------------------------------- c XC & YC c---------------------------------------------------- MSIZE = NxP*NyP c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/XC.data', & form='unformatted') c read (1,REC=1) Xo_P(:,:) close(1) c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/YC.data', & form='unformatted') c read (1,REC=1) Yo_P(:,:) close(1) c---------------------------------------------------- c XG & YG c---------------------------------------------------- MSIZE = NxP*NyP c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/XG.data', & form='unformatted') read (1,REC=1) Xg_P(:,:) close(1) c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/YG.data', & form='unformatted') c read (1,REC=1) Yg_P(:,:) close(1) c---------------------------------------------------- c Yu c---------------------------------------------------- DO J = 1,NyP DO I = 1,NxP c Xu_P(I,J) = Xg_P(I,J) Yu_P(I,J) = Yo_P(I,J) ENDDO ENDDO c---------------------------------------------------- c Xv & Yv c---------------------------------------------------- DO J = 1,NyP DO I = 1,NxP Xv_P(I,J) = Xo_P(I,J) Yv_P(I,J) = Yg_P(I,J) ENDDO ENDDO c---------------------------------------------------- c hFacC c---------------------------------------------------- MSIZE = NxP*NyP*NrP c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/hFacC.data', & form='unformatted') read (1,REC=1) hFacC_P(:,:,:) close(1) c---------------------------------------------------- c hFacW c---------------------------------------------------- MSIZE = NxP*NyP*NrP c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/hFacW.data', & form='unformatted') read (1,REC=1) hFacW_P(:,:,:) close(1) c---------------------------------------------------- c hFacS c---------------------------------------------------- MSIZE = NxP*NyP*NrP c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/hFacS.data', & form='unformatted') read (1,REC=1) hFacS_P(:,:,:) close(1) c---------------------------------------------------- c RAC c---------------------------------------------------- MSIZE = NxP*NyP c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/RAC.data', & form='unformatted') read (1,REC=1) RAC_P(:,:) close(1) c---------------------------------------------------- c RAW c---------------------------------------------------- MSIZE = NxP*NyP c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/RAW.data', & form='unformatted') read (1,REC=1) RAW_P(:,:) close(1) c---------------------------------------------------- c RAS c---------------------------------------------------- MSIZE = NxP*NyP c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/PARENT/RAS.data', & form='unformatted') read (1,REC=1) RAS_P(:,:) close(1) c---------------------------------------------------- c MASK x PARENT c---------------------------------------------------- DO K = 1,NrP DO J = 1,NyP DO I = 1,NxP DEEP_P(i,j,k) = 0. IF (hFacC_P(i,j,k).ne.0) then DEEP_P(I,J,K) = 1. ENDIF ENDDO ENDDO ENDDO c---------------------------------------------------- c 1/Volume (C) c---------------------------------------------------- DO K = 1,NrP DO J = 1,NyP DO I = 1,NxP INV_VOL_C_P(I,J,K) = 1. IF ((RAC_P(I,J)*hFacC_P(I,J,K)).ne.0.) THEN INV_VOL_C_P(I,J,K) = 1./(RAC_P(I,J)*hFacC_P(I,J,K)) ENDIF ENDDO ENDDO ENDDO c---------------------------------------------------- c 1/Volume (W) c---------------------------------------------------- DO K = 1,NrP DO J = 1,NyP DO I = 1,NxP INV_VOL_W_P(I,J,K) = 1. IF ((RAW_P(I,J)*hFacW_P(I,J,K)).ne.0.) THEN INV_VOL_W_P(I,J,K) = 1./(RAW_P(I,J)*hFacW_P(I,J,K)) ENDIF ENDDO ENDDO ENDDO c---------------------------------------------------- c 1/Volume (S) c---------------------------------------------------- DO K = 1,NrP DO J = 1,NyP DO I = 1,NxP INV_VOL_S_P(I,J,K) = 1. IF ((RAS_P(I,J)*hFacS_P(I,J,K)).ne.0.) THEN INV_VOL_S_P(I,J,K) = 1./(RAS_P(I,J)*hFacS_P(I,J,K)) ENDIF ENDDO ENDDO ENDDO c-------------------------------------------------------- c CHILD MODEL c-------------------------------------------------------- print*,' [2] Read CHILD model geometry' c---------------------------------------------------- c XC & YC c---------------------------------------------------- MSIZE = NxC*NyC c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/XC.data', & form='unformatted') read (1,REC=1) Xo_C(:,:) close(1) c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/YC.data', & form='unformatted') c read (1,REC=1) Yo_C(:,:) close(1) c---------------------------------------------------- c XG & YG c---------------------------------------------------- MSIZE = NxC*NyC c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/XG.data', & form='unformatted') read (1,REC=1) Xg_C(:,:) close(1) c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/YG.data', & form='unformatted') read (1,REC=1) Yg_C(:,:) close(1) c---------------------------------------------------- c Xu & Yu c---------------------------------------------------- DO J = 1,NyC DO I = 1,NxC Xu_C(I,J) = Xg_C(I,J) Yu_C(I,J) = Yo_C(I,J) ENDDO ENDDO c---------------------------------------------------- c Xv & Yv c---------------------------------------------------- DO J = 1,NyC DO I = 1,NxC Xv_C(I,J) = Xo_C(I,J) Yv_C(I,J) = Yg_C(I,J) ENDDO ENDDO c---------------------------------------------------- c hFacC c---------------------------------------------------- MSIZE = NxC*NyC*NrC c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/hFacC.data', & form='unformatted') read (1,REC=1) hFacC_C(:,:,:) close(1) c---------------------------------------------------- c hFacW c---------------------------------------------------- MSIZE = NxC*NyC*NrC c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/hFacW.data', & form='unformatted') read (1,REC=1) hFacW_C(:,:,:) close(1) c---------------------------------------------------- c hFacC c---------------------------------------------------- MSIZE = NxC*NyC*NrC c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/hFacS.data', & form='unformatted') read (1,REC=1) hFacS_C(:,:,:) close(1) c---------------------------------------------------- c RAC c---------------------------------------------------- MSIZE = NxC*NyC c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/RAC.data', & form='unformatted') read (1,REC=1) RAC_C(:,:) close(1) c---------------------------------------------------- c RAW c---------------------------------------------------- MSIZE = NxC*NyC c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/RAW.data', & form='unformatted') read (1,REC=1) RAW_C(:,:) close(1) c---------------------------------------------------- c RAS c---------------------------------------------------- MSIZE = NxC*NyC c open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & convert='big_endian', & file=trim(dirNEST)//'/CHILD/RAS.data', & form='unformatted') read (1,REC=1) RAS_C(:,:) close(1) c---------------------------------------------------- c MASK x CHILD c---------------------------------------------------- DO K = 1,NrC DO J = 1,NyC DO I = 1,NxC DEEP_C(i,j,k) = 0. IF (hFacC_C(i,j,k).ne.0) then DEEP_C(I,J,K) = 1. ENDIF ENDDO ENDDO ENDDO c---------------------------------------------------- c c __/________ ___________ c | | | | || c > o > | | | c |__/__|_____| c | | | c > o > | c |_____|_____|_____|_____| c c c---------------------------------------------------- print*,' [3] Compute J index P-->C' c-------------------------------------------------------- c Compute J indicies for mapping P->C (I) c-------------------------------------------------------- I = 1 II = WesternB c DO J = 1,NyC P2C_U(J) = 0. DO JJ = 1,NyP-1 YF = Yg_C(I,J) YP1 = Yg_P(II,JJ) YP2 = Yg_P(II,JJ+1) IF (YF.ge.YP1.and.YF.lt.YP2) THEN P2C_U(J) = JJ ENDIF ENDDO ENDDO c-------------------------------------------------------- c Compute J indicies for mapping P->C (II) c-------------------------------------------------------- I = 1 II = WesternB c DO J = 1,NyC P2C_linU(J) = 0. DO JJ = 1,NyP-1 YF = Yu_C(I,J) YP1 = Yu_P(II,JJ) YP2 = Yu_P(II,JJ+1) IF (YF.ge.YP1.and.YF.lt.YP2) THEN P2C_linU(J) = JJ ENDIF ENDDO ENDDO c-------------------------------------------------------- c Compute J indicies for mapping P->C (III) c-------------------------------------------------------- I = 1 II = WesternB c DO J = 1,NyC DO JJ = 1,NyP-1 YF = Yu_C(I,J) YP1 = Yu_P(II,JJ) IF (YF.eq.YP1) THEN WO3_linU(J) = 0 if (J+1.le.NyC) WO3_linU(J+1) = 1 if (J+2.le.NyC) WO3_linU(J+2) = 2 ENDIF ENDDO ENDDO c--------------------Lower bound DO J = 1,NyC DO JJ = 1,NyP-1 YF = Yu_C(I,J) YP1 = Yu_P(II,JJ) IF (YF.eq.YP1) THEN WO3_linU(J) = 0 if (J-1.gt.0) WO3_linU(J-1) = 2 if (J-2.gt.0) WO3_linU(J-2) = 1 GOTO 2345 ENDIF ENDDO ENDDO 2345 CONTINUE c--------------------Upper bound DO J = NyC,1,-1 DO JJ = 1,NyP-1 YF = Yu_C(I,J) YP1 = Yu_P(II,JJ) IF (YF.eq.YP1) THEN WO3_linU(J) = 0 if (J+1.le.NyC) WO3_linU(J+1) = 1 if (J+2.le.NyC) WO3_linU(J+2) = 2 GOTO 2346 ENDIF ENDDO ENDDO 2346 CONTINUE c-------------------------------------------------------- c Compute J indicies for mapping P->C (IV) c-------------------------------------------------------- I = 1 II = WesternB c DO J = 1,NyC P2C_linV(J) = 0. DO JJ = 1,NyP-1 YF = Yv_C(I,J) YP1 = Yv_P(II,JJ) YP2 = Yv_P(II,JJ+1) IF (YF.ge.YP1.and.YF.lt.YP2) THEN P2C_linV(J) = JJ ENDIF ENDDO ENDDO c-------------------------------------------------------- c Compute J indicies for mapping P->C (V) c-------------------------------------------------------- I = 1 II = WesternB c DO J = 1,NyC DO JJ = 1,NyP-1 YF = Yv_C(I,J) YP1 = Yv_P(II,JJ) IF (YF.eq.YP1) THEN WO3_linV(J) = 0 if (J+1.le.NyC) WO3_linV(J+1) = 1 if (J+2.le.NyC) WO3_linV(J+2) = 2 ENDIF ENDDO ENDDO c--------------------Lower bound DO J = 1,NyC DO JJ = 1,NyP-1 YF = Yv_C(I,J) YP1 = Yv_P(II,JJ) IF (YF.eq.YP1) THEN WO3_linV(J) = 0 if (J-1.gt.0) WO3_linV(J-1) = 2 if (J-2.gt.0) WO3_linV(J-2) = 1 GOTO 23451 ENDIF ENDDO ENDDO 23451 CONTINUE c--------------------Upper bound DO J = NyC,1,-1 DO JJ = 1,NyP-1 YF = Yv_C(I,J) YP1 = Yv_P(II,JJ) IF (YF.eq.YP1) THEN WO3_linV(J) = 0 if (J+1.le.NyC) WO3_linV(J+1) = 1 if (J+2.le.NyC) WO3_linV(J+2) = 2 GOTO 23461 ENDIF ENDDO ENDDO 23461 CONTINUE c-------------------------------------------------------- c Compute J indicies for mapping P->C (V) c-------------------------------------------------------- print*,' [5] Compute J index P-->C for (o)' I = 1 II = WesternB c DO J = 1,NyC P2C_o(J) = 0. DO JJ = 1,NyP-1 YF = Yo_C(I,J) YP1 = Yg_P(II,JJ) YP2 = Yg_P(II,JJ+1) IF (YF.gt.YP1.and.YF.lt.YP2) THEN P2C_o(J) = JJ ENDIF ENDDO ENDDO c-------------------------------------------------------- c Compute J indicies for mapping P->C (VI) c-------------------------------------------------------- print*,' [6] Compute J index P-->C for (v bilinear)' I = 1 II = WesternB c DO J = 1,NyC DO JJ = 2,NyP-1 YF = Yv_C(I,J) YP1 = Yv_P(II,JJ) YP2 = Yv_P(II,JJ+1) YP3 = Yv_P(II,JJ-1) c IF (YF.ge.YP1.and.YF.lt.YP2) THEN P2C1_V(J) = JJ P2C2_V(J) = JJ+1 ENDIF ENDDO ENDDO c-------------------------------------------------------- c Look for the 9 CHILD indicies in PARENT grid cell c-------------------------------------------------------- print*,' [8] Compute I J index C-->P for (o)' c DO J = 1,NyP DO I = 1,NxP I_C2P(:,I,J) = 0 J_C2P(:,I,J) = 0 c DO JJ = 1,NyC DO II = 1,NxC IF (Xo_C(II,JJ).eq.Xo_P(I,J).and. & Yo_C(II,JJ).eq.Yo_P(I,J)) then KK = 0 DO JJJ = JJ-1,JJ+1 DO III = II-1,II+1 KK = kk +1 if (III.lt.1.or.III.gt.NxC) cycle if (JJJ.lt.1.or.JJJ.gt.NyC) cycle I_C2P(KK,I,J) = III J_C2P(KK,I,J) = JJJ ENDDO ENDDO ENDIF ENDDO ENDDO c ENDDO ENDDO ENDIF c-------------------------------------------------------- c Broadcast all the above variables c-------------------------------------------------------- call MPI_BCAST(I_C2P,9*NxP*NyP,MPI_INTEGER, & 0,NEST_COMM,ierr) call MPI_BCAST(J_C2P,9*NxP*NyP,MPI_INTEGER, & 0,NEST_COMM,ierr) call MPI_BCAST(RAC_C,NxC*NyC,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(hFacC_C,NxC*NyC*NrC,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(INV_VOL_C_P,NxP*NyP*NrP,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(RAW_C,NxC*NyC,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(hFacW_C,NxC*NyC*NrC,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(INV_VOL_W_P,NxP*NyP*NrP,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(RAS_C,NxC*NyC,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(hFacS_C,NxC*NyC*NrC,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(INV_VOL_S_P,NxP*NyP*NrP,MPI_REAL, & 0,NEST_COMM,ierr) c call MPI_BCAST(DEEP_C,NxC*NyC*NrC,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(RAC_P,NxP*NyP,MPI_REAL, & 0,NEST_COMM,ierr) c call MPI_BCAST(IM_P,1,MPI_INTEGER, & 0,NEST_COMM,ierr) call MPI_BCAST(JM_P,1,MPI_INTEGER, & 0,NEST_COMM,ierr) call MPI_BCAST(index_var3D,1,MPI_INTEGER, & 0,NEST_COMM,ierr) call MPI_BCAST(index_var2D,1,MPI_INTEGER, & 0,NEST_COMM,ierr) c call MPI_BCAST(DEEP_P,NxP*NyP*NrP,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(hFacS_P,NxP*NyP*NrP,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(hFacC_P,NxP*NyP*NrP,MPI_REAL, & 0,NEST_COMM,ierr) call MPI_BCAST(hFacW_P,NxP*NyP*NrP,MPI_REAL, & 0,NEST_COMM,ierr) c-------------------------------------------------------- if(rank.eq.0) then c-------------------------------------------------------- DO K = 1,NrP DO J = 1,NyP DO I = WesternB+1,EasternB-1 cc WesternB side DO II = 1,9 IF (I_C2P(II,I,J).eq.0) cycle IF (J_C2P(II,I,J).eq.0) cycle c Indx = I_C2P(II,I,J) Jndx = J_C2P(II,I,J) ENDDO ENDDO ENDDO ENDDO c--------------------------------------------------------- ONOFF=0 endif c-------------------------------------------------------- c BEGIN MAIN LOOP c-------------------------------------------------------- do if(rank.eq.0) then c-------------------------------------------------------- c (1) READ FROM PARENT MODEL c-------------------------------------------------------- ICONT=1 DO WHILE(ICONT.le.nSxP*nSyP) from= MPI_ANY_SOURCE call MPI_RECV (globalPA, index, MPI_REAL8, & FROM, 3000, & MPI_COMM_World, status,ierr) c ICONT=ICONT+1 c whm=status(MPI_SOURCE)-MSTR_PRNT(NST_LEV)+1 c call MPI_GET_COUNT(status,MPI_REAL8,st_count,ierr) c DO II = 1,6 IF (globalPA(II,1,1,1).ne.-999.) THEN globalP1(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,1) globalP2(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,2) globalP3(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,3) globalP4(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,4) globalP5(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,5) ENDIF ENDDO ENDDO c-------------------------------------------------------- c Start interpolation for CHILD c-------------------------------------------------------- CALL INTERPOLATION_P2C ( & globalP1,globalP2,globalP3,globalP4,globalP5, & NxP,NyP,NrP, & NxC,NyC,NrC, $ WesternB,EasternB, $ P2C_U,P2C_V,P2C_o,P2C1_V,P2C2_V,P2C1_o,P2C2_o, $ P2C_linU,WO3_linU,P2C_linV,WO3_linV, $ Xv_C,Yv_C,Xv_P,Yv_P, $ T_C1,S_C1,U_C1,V_C1,ETA_C1, $ DEEP_C,DEEP_P & ) c============================================================== c Open Files from PARENT MODEL c============================================================== ICONT=1 do while(ICONT.le.nSxP*nSyP) from= MPI_ANY_SOURCE call MPI_RECV (globalPA, index, MPI_REAL8, & FROM, 3000, & MPI_COMM_World, status,ierr) ICONT=ICONT+1 whm=status(MPI_SOURCE)-MSTR_PRNT(NST_LEV)+1 call MPI_GET_COUNT(status,MPI_REAL8,st_count,ierr) DO II = 1,6 IF (globalPA(II,1,1,1).ne.-999.) THEN globalP1(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,1) globalP2(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,2) globalP3(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,3) globalP4(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,4) globalP5(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) = & globalPA(II,1:JM_P,:,5) ENDIF ENDDO end do c-------------------------------------------------------- c Start inteprolation for CHILD c-------------------------------------------------------- CALL INTERPOLATION_P2C ( & globalP1,globalP2,globalP3,globalP4,globalP5, & NxP,NyP,NrP, & NxC,NyC,NrC, $ WesternB,EasternB, $ P2C_U,P2C_V,P2C_o,P2C1_V,P2C2_V,P2C1_o,P2C2_o, $ P2C_linU,WO3_linU,P2C_linV,WO3_linV, $ Xv_C,Yv_C,Xv_P,Yv_P, $ T_C2,S_C2,U_C2,V_C2,ETA_C2, $ DEEP_C,DEEP_P & ) c============================================================== c Temporal Interpolation OBCs x CHILD MODEL c============================================================== c 0 1200 c ---+--.--.--+---- Parent c c |--|--|-- c 0 800 c 400 c------------------------------------------------------------ DO I = 1,2 ! WesternB & EasternB DIFF_T(:,:,I) = (T_C2(:,:,I) - T_C1(:,:,I))/3. DIFF_S(:,:,I) = (S_C2(:,:,I) - S_C1(:,:,I))/3. DIFF_U(:,:,I) = (U_C2(:,:,I) - U_C1(:,:,I))/3. DIFF_V(:,:,I) = (V_C2(:,:,I) - V_C1(:,:,I))/3. DIFF_eta(:,:,I) = (eta_C2(:,:,I) - eta_C1(:,:,I))/3. ENDDO c------------------------------------------------------------- c c Step 0 (Rec = 1 ==> WesternB) c------- (Rec = 2 ==> EasternB) DO I = 1,2 !WesternB & EasternB T_C1(:,:,I) = T_C2(:,:,I) !+ DIFF_T(:,:,I) S_C1(:,:,I) = S_C2(:,:,I) !+ DIFF_S(:,:,I) U_C1(:,:,I) = U_C2(:,:,I) !+ DIFF_U(:,:,I) V_C1(:,:,I) = V_C2(:,:,I) !+ DIFF_V(:,:,I) ETA_C1(:,:,I) = ETA_C2(:,:,I) !+ DIFF_ETA(:,:,I) ENDDO if(ONOFF.eq.0) then c--------------------------------------------------------------------- ICONT = -1 DO I = 1,nSxC DO J = 1,nSyC ICONT = ICONT + 1 IndI = IM_C*(I-1) IndJ = JM_C*(J-1) VAR_C1(:,:,:,:) = 0. c J1 = 1+IndJ-OLY J2 = JM_C+IndJ+OLY c JJ1 = 1 JJ2 = JM_C+OLY+OLY c IF(1 +IndJ-OLY.LT.0) THEN J1 = 1 JJ1 = 4 ENDIF c IF(JM_C+IndJ+OLY.GT.NyC) THEN J2 = NyC JJ2 = JM_C ENDIF c VAR_C1(JJ1:JJ2,:,:,1) = U_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,2) = V_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,3) = T_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,4) = S_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,5) = ETA_C1(J1:J2,:,:) c call MPI_SEND (VAR_C1, indexF, MPI_REAL8, & MSTR_CHLD(NST_LEV)+ICONT, 3000, & MPI_Comm_World,ierr) c ENDDO ENDDO c---------------------------------------------------------------------- cc write(*,*) 'VIC: MANDO SEGNALE DI OK AL CHILD PER INIZIALIZZARE' ONOFF=1 ENDIF cc write(*,*) 'VIC: MANDO SEGNALE DI OK AL CHILD PER IL PASSO 1' c----------------------------------------------------------------------- ICONT = -1 DO I = 1,nSxC DO J = 1,nSyC ICONT = ICONT + 1 IndI = IM_C*(I-1) IndJ = JM_C*(J-1) c VAR_C1(:,:,:,:) = 0. c J1 = 1+IndJ-OLY J2 = JM_C+IndJ+OLY c JJ1 = 1 JJ2 = JM_C+OLY+OLY c IF(1 +IndJ-OLY.LT.0) THEN J1 = 1 JJ1 = 4 ENDIF c IF(JM_C+IndJ+OLY.GT.NyC) THEN J2 = NyC JJ2 = JM_C ENDIF c VAR_C1(JJ1:JJ2,:,:,1) = U_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,2) = V_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,3) = T_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,4) = S_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,5) = ETA_C1(J1:J2,:,:) c call MPI_SEND (VAR_C1, indexF, MPI_REAL8, & MSTR_CHLD(NST_LEV)+ICONT, 3000, & MPI_Comm_World,ierr) ENDDO ENDDO c--------------------------------------------------------------------- goto 8888 c c Step 1 (Rec = 3 ==> WesternB) c------- (Rec = 4 ==> EasternB) DO I = 1,2 !WesternB & EasternB T_C1(:,:,I) = T_C2(:,:,I) !+ DIFF_T(:,:,I) S_C1(:,:,I) = S_C2(:,:,I) !+ DIFF_S(:,:,I) U_C1(:,:,I) = U_C2(:,:,I) !+ DIFF_U(:,:,I) V_C1(:,:,I) = V_C2(:,:,I) !+ DIFF_V(:,:,I) ETA_C1(:,:,I) = ETA_C2(:,:,I) !+ DIFF_ETA(:,:,I) ENDDO c---------------------------------------------------------- ICONT = -1 DO I = 1,nSxC DO J = 1,nSyC ICONT = ICONT + 1 IndI = IM_C*(I-1) IndJ = JM_C*(J-1) VAR_C1(:,:,:,:) = 0. c J1 = 1+IndJ-OLY J2 = JM_C+IndJ+OLY c JJ1 = 1 JJ2 = JM_C+OLY+OLY c IF(1 +IndJ-OLY.LT.0) THEN J1 = 1 JJ1 = 4 ENDIF c IF(JM_C+IndJ+OLY.GT.NyC) THEN J2 = NyC JJ2 = JM_C ENDIF c VAR_C1(JJ1:JJ2,:,:,1) = U_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,2) = V_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,3) = T_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,4) = S_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,5) = ETA_C1(J1:J2,:,:) call MPI_SEND (VAR_C1, indexF, MPI_REAL8, & MSTR_CHLD(NST_LEV)+ICONT, 3000, & MPI_Comm_World,ierr) ENDDO ENDDO c----------------------------------------------------------- c c Step 2 (Rec = 5 ==> WesternB) c------- (Rec = 6 ==> EasternB) DO I = 1,2 !WesternB & EasternB T_C1(:,:,I) = T_C2(:,:,I) !+ DIFF_T(:,:,I) S_C1(:,:,I) = S_C2(:,:,I) !+ DIFF_S(:,:,I) U_C1(:,:,I) = U_C2(:,:,I) !+ DIFF_U(:,:,I) V_C1(:,:,I) = V_C2(:,:,I) !+ DIFF_V(:,:,I) ETA_C1(:,:,I) = ETA_C2(:,:,I) !+ DIFF_ETA(:,:,I) ENDDO c---------------------------------------------------------- ICONT = -1 DO I = 1,nSxC DO J = 1,nSyC ICONT = ICONT + 1 IndI = IM_C*(I-1) IndJ = JM_C*(J-1) VAR_C1(:,:,:,:) = 0. c J1 = 1+IndJ-OLY J2 = JM_C+IndJ+OLY c JJ1 = 1 JJ2 = JM_C+OLY+OLY c IF(1 +IndJ-OLY.LT.0) THEN J1 = 1 JJ1 = 4 ENDIF c IF(JM_C+IndJ+OLY.GT.NyC) THEN J2 = NyC JJ2 = JM_C ENDIF c VAR_C1(JJ1:JJ2,:,:,1) = U_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,2) = V_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,3) = T_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,4) = S_C1(J1:J2,:,:) VAR_C1(JJ1:JJ2,:,:,5) = ETA_C1(J1:J2,:,:) call MPI_SEND (VAR_C1, indexF, MPI_REAL8, & MSTR_CHLD(NST_LEV)+ICONT, 3000, & MPI_Comm_World,ierr) ENDDO ENDDO c--------------------------------------------------------------- 8888 CONTINUE c-------------------------------------------------------- c Close OBCs Files x CHILD MODEL c-------------------------------------------------------- c---------------------------------------------------- c------------- MANDO SEGNALE DI OK AL CHILD c---------------------------------------------------- c-------------------------------------------------------- c (1) READ FROM CHILD MODEL c-------------------------------------------------------- call MPI_RECV (TRANSPORT_WEST, 1, MPI_REAL8, & MSTR_CHLD(NST_LEV), 3000, & MPI_COMM_World, status,ierr) call MPI_RECV (TRANSPORT_EAST, 1, MPI_REAL8, & MSTR_CHLD(NST_LEV), 3000, & MPI_COMM_World, status,ierr) c--------------------------------------------------------- c--------------------------------------------------------- ICONT=1 DO WHILE(ICONT.le.nSxC*nSyC) from= MPI_ANY_SOURCE call MPI_RECV (globalC3D_a,indexF, MPI_REAL8, & from, 3000, MPI_COMM_World, status,ierr) ICONT=ICONT+1 whm=status(MPI_SOURCE)-MSTR_CHLD(NST_LEV)+1 call MPI_GET_COUNT(status,MPI_REAL8,st_count,ierr) globalC3D(1+IndI_C(whm):IM_C+IndI_C(whm), & 1+IndJ_C(whm):JM_C+IndJ_C(whm),:,:)= & globalC3D_a(:,:,:,:) END DO c----------------------------- ICONT=1 DO WHILE(ICONT.le.nSxC*nSyC) from= MPI_ANY_SOURCE call MPI_RECV (globalC2D_a,index1F, MPI_REAL8, & from, 3000, MPI_COMM_World, status,ierr) ICONT=ICONT+1 whm=status(MPI_SOURCE)-MSTR_CHLD(NST_LEV)+1 call MPI_GET_COUNT(status,MPI_REAL8,st_count,ierr) globalC2D(1+IndI_C(whm):IM_C+IndI_C(whm), & 1+IndJ_C(whm):JM_C+IndJ_C(whm),:)= & globalC2D_a(:,:,:) END DO ENDIF call MPI_BCAST(globalC3D,NxC*NyC*NrC*15,MPI_REAL8, & 0,NEST_COMM,ierr) call MPI_BCAST(globalC2D,NxC*NyC*4,MPI_REAL8, & 0,NEST_COMM,ierr) 2323 CONTINUE c======================================================= c (1) READ FROM CHILD MODEL c======================================================= c======================================================= c (2) INTERPOLATIONS c======================================================= c 3D VAR c-------- DO iVar = 1,15 ! tipo di variabile DO K = 1,NrP DO J = 1,NyP DO I = WesternB+1,EasternB-1 VAR3D_P(I,J,K,iVar) = 0. ! inizializzo c WesternB side AREA_VOL = 0. !can be area or volume depend on the variable SELECT CASE(iVar) CASE(1,5,9) I_START = 1 I_END = 9 I_STEP = 1 !3 CASE(2,6,10) I_START = 1 I_END = 9 !3 I_STEP = 1 CASE DEFAULT I_START = 1 I_END = 9 I_STEP = 1 END SELECT DO II = I_START,I_END,I_STEP IF (I_C2P(II,I,J).eq.0) cycle IF (J_C2P(II,I,J).eq.0) cycle c Indx = I_C2P(II,I,J) Jndx = J_C2P(II,I,J) c c SELECT CASE(iVar) CASE (1,5,9) VAR3D_P(I,J,K,ivar) = VAR3D_P(I,J,K,iVar) + & globalC3D(Indx,Jndx,K,iVar)* $ RAW_C(Indx,Jndx)* & hFacW_C(Indx,Jndx,K) CASE (2,6,10) VAR3D_P(I,J,K,ivar) = VAR3D_P(I,J,K,iVar) + & globalC3D(Indx,Jndx,K,iVar)* $ RAS_C(Indx,Jndx)* & hFacS_C(Indx,Jndx,K) CASE DEFAULT VAR3D_P(I,J,K,ivar) = VAR3D_P(I,J,K,iVar) + & globalC3D(Indx,Jndx,K,iVar)* $ RAC_C(Indx,Jndx)* & hFacC_C(Indx,Jndx,K) AREA_VOL = AREA_VOL + & RAC_C(Indx,Jndx)* hFacC_C(Indx,Jndx,K) END SELECT ENDDO c----------------------------------------------- c Make a volume average c---------------------------------------------- SELECT CASE(iVar) CASE (1,5,9) VAR3D_P(I,J,K,ivar) = & VAR3D_P(I,J,K,iVar)* & INV_VOL_W_P(I,J,K) if (hFacW_P(I,J,K).eq.0.) VAR3D_P(I,J,K,ivar)=0. CASE (2,6,10) VAR3D_P(I,J,K,ivar) = & VAR3D_P(I,J,K,iVar)* & INV_VOL_S_P(I,J,K) if (hFacS_P(I,J,K).eq.0.) VAR3D_P(I,J,K,ivar)=0. CASE DEFAULT IF (AREA_VOL.ne.0.) then VAR3D_P(I,J,K,ivar) = & VAR3D_P(I,J,K,iVar)/AREA_VOL ENDIF if (hFacC_P(I,J,K).eq.0.) VAR3D_P(I,J,K,ivar)=0. END SELECT ENDDO ENDDO ENDDO ENDDO c 2D VAR c-------- DO iVar = 1,4 DO J = 1,NyP DO I = WesternB+1,EasternB-1 VAR2D_P(I,J,iVar) = 0. AREA_VOL = 0. DO II = 1,9 IF (I_C2P(II,I,J).eq.0) cycle IF (J_C2P(II,I,J).eq.0) cycle c Indx = I_C2P(II,I,J) Jndx = J_C2P(II,I,J) c VAR2D_P(I,J,ivar) = VAR2D_P(I,J,iVar) + & globalC2D(Indx,Jndx,iVar)* $ RAC_C(Indx,Jndx)* & DEEP_C(Indx,Jndx,1) AREA_VOL = AREA_VOL + & RAC_C(Indx,Jndx)* DEEP_C(Indx,Jndx,1) ENDDO c----------------------------- IF ((RAC_P(I,J)*DEEP_P(I,J,1)).ne.0.) then c IF (AREA_VOL.ne.0.) then VAR2D_P(I,J,ivar) = & VAR2D_P(I,J,iVar)/ & RAC_P(I,J) ENDIF c---------------------------- ENDDO ENDDO ENDDO if(rank.eq.0) then c-------------------------------------------------------- c Write Files for PARENT MODEL c-------------------------------------------------------- c print*,' (*) Open Files for PARENT MODEL' 7575 CONTINUE c---------------------------------------------------- c------------- MANDO SEGNALE DI OK AL PARENT c---------------------------------------------------- call MPI_SEND (TRANSPORT_WEST, 1, MPI_REAL8, & MSTR_PRNT(NST_LEV), 3000, & MPI_Comm_World,ierr) call MPI_SEND (TRANSPORT_EAST, 1, MPI_REAL8, & MSTR_PRNT(NST_LEV), 3000, & MPI_Comm_World,ierr) ENDIF c--------------------------------------------------------- !--------------------------------------------------------- VCONT=VCONTP(rank) DO I = vstart,vstop DO J = 1,nSyP VCONT = VCONT + 1 IndI = IM_P*(I-1) IndJ = JM_P*(J-1) c----------------------------------------------------------- DO iVar=1,15 CALL MPI_SEND (VAR3D_P(1+IndI:IM_P+IndI & ,1+IndJ:JM_P+IndJ,:,iVar), & index_var3D,MPI_REAL8,MSTR_PRNT(NST_LEV)+VCONT, & 3000,MPI_Comm_World,ierr) ENDDO DO iVar=1,4 call MPI_SEND (VAR2D_P(1+IndI:IM_P+IndI & ,1+IndJ:JM_P+IndJ,iVar), & index_var2D,MPI_REAL8,MSTR_PRNT(NST_LEV)+VCONT, & 3000,MPI_Comm_World,ierr) ENDDO c----------------------------------------------------------- END DO END DO call MPI_BARRIER(NEST_COMM,ierr) end do c--------------------------------------------------------- c======================================================= c END MAIN LOOP c======================================================= call MPI_FINALIZE(ierr) c--------------------------------------------------------- !--------------------------------------------------------- 101 FORMAT (I1) STOP END SUBROUTINE INTERPOLATION_P2C ( & globalP1,globalP2,globalP3,globalP4,globalP5, & NxP,NyP,NrP, & NxC,NyC,NrC, $ WesternB,EasternB, $ P2C_U,P2C_V,P2C_o,P2C1_V,P2C2_V,P2C1_o,P2C2_o, $ P2C_linU,WO3_linU,P2C_linV,WO3_linV, $ Xv_F,Yv_F,Xv_P,Yv_P, $ T_F,S_F,U_F,V_F,ETA_F, $ DEEP_F,DEEP_P & ) c---------------------------------------------------- implicit none c---------------------------------------------------- INTEGER :: I,J,K,II,JJ INTEGER :: WesternB,EasternB INTEGER :: NrP,NxP,NyP INTEGER :: NrC,NxC,NyC c---------------------------------------------------- REAL*8 :: Fp,Fm,Fo,VEL_MEMO INTEGER :: INDC c---------------------------------------------------- c Define Global Variables to Exchange c---------------------------------------------------- REAL*8 :: globalP1(6,NyP,NrP) REAL*8 :: globalP2(6,NyP,NrP) REAL*8 :: globalP3(6,NyP,NrP) REAL*8 :: globalP4(6,NyP,NrP) REAL*8 :: globalP5(6,NyP,NrP) c---------------------------------------------------- c Define CHILD Model Geometry c---------------------------------------------------- REAL*4 :: Xu_F(NxC,NyC) REAL*4 :: Yu_F(NxC,NyC) REAL*4 :: Xv_F(NxC,NyC) REAL*4 :: Yv_F(NxC,NyC) REAL*4 :: Xo_F(NxC,NyC) REAL*4 :: Yo_F(NxC,NyC) REAL*4 :: Xg_F(NxC,NyC) REAL*4 :: Yg_F(NxC,NyC) REAL*4 :: DEEP_F(NxC,NyC,NrC) c---------------------------------------------------- c Define PARENT Model Geometry c---------------------------------------------------- c REAL*4 :: Xu_P(NxP,NyP) REAL*4 :: Yu_P(NxP,NyP) REAL*4 :: Xv_P(NxP,NyP) REAL*4 :: Yv_P(NxP,NyP) REAL*4 :: Xo_P(NxP,NyP) REAL*4 :: Yo_P(NxP,NyP) REAL*4 :: Xg_P(NxP,NyP) REAL*4 :: Yg_P(NxP,NyP) REAL*4 :: DEEP_P(NxP,NyP,NrP) REAL*4 :: DEPDEP c----------------------------------------------------- REAL*4 :: X1,X2,X3,X4,Y1,Y2,Y3,Y4 REAL*4 :: f1,f2,f3,f4,f,x,y REAL*4 :: gammaT,gammaS,terzo,dueterzi REAL*4 :: gammaEta REAL*4 :: gammaV c---------------------------------------------------- c Define INDICIES MATRIX c---------------------------------------------------- INTEGER :: P2C_U(NyC) !x imposizione NETTA INTEGER :: P2C_linU(NyC) !x Lineare INTEGER :: WO3_linU(NyC) !x Lineare !Which Of 3 INTEGER :: P2C_linV(NyC) !x Lineare INTEGER :: WO3_linV(NyC) !x Lineare !Which Of 3 INTEGER :: P2C_V(NyC) !x Lineare INTEGER :: P2C_o(NyC) !x Lineare INTEGER :: P2C1_V(NyC) !x BiLineare INTEGER :: P2C2_V(NyC) !x BiLineare INTEGER :: P2C1_o(NyC) !x BiLineare INTEGER :: P2C2_o(NyC) !x BiLineare REAL*8 :: diff(NrC) REAL*8 :: DEPDEP_F_WesternB(NrC) REAL*8 :: DEPDEP_F_EasternB(NrC) c---------------------------------------------------- c Define CHILD model variable c---------------------------------------------------- c _____________ (1) WesternB (2) EasternB c | REAL*8 :: U_F(NyC,NrC,2) REAL*8 :: V_F(NyC,NrC,2) REAL*8 :: T_F(NyC,NrC,2) REAL*8 :: S_F(NyC,NrC,2) REAL*8 :: ETA_F(NyC,NrC,2) c---------------------------------------------------- PARAMETER ( terzo = 1./3.) PARAMETER (dueterzi = 2./3.) c======================================================= c (2) INTERPOLATIONS c======================================================= c (2.1) Linear for normal velocity component (u in this case) c------------------------------------------------------- DO K = 1,NrP DO J = 1,NyP-1 IF (globalP4(2,J,K).eq.0.) THEN !uso la salinità come discriminante globalP1 (2,J,K) = globalP1 (2,J+1,K) ENDIF ENDDO DO J = NyP,2,-1 IF (globalP4(2,J,K).eq.0.) THEN !uso la salinità come discriminante globalP1 (2,J,K) = globalP1 (2,J-1,K) ENDIF ENDDO ENDDO c======================================================= c (2.1) NOT Linear but simply imposed c------------------------------------------------------- DO J = 1,3 U_F(J,:,1) = globalP1(2,P2C_U(J),:) U_F(J,:,2) = globalP1(6,P2C_U(J),:) ENDDO DO J = NyC-2,NyC U_F(J,:,1) = globalP1(2,P2C_U(J),:) U_F(J,:,2) = globalP1(6,P2C_U(J),:) ENDDO c================================================= DO J = 4,NyC-3,3 INDC = P2C_U(J) DO K = 1,NrC c-------- WesternB ---------------------- Fp = globalP1(2,INDC+1,K) Fo = globalP1(2,INDC,K) Fm = globalP1(2,INDC-1,K) c VEL_MEMO = 0. DO I = -1,1 U_F(J+1+i,K,1) = ((Fp-2.*Fo+Fm)/24.)* & ((12.*float(i)**2+1.)/9.)+ & ((float(i)/6.)*(Fp-Fm))+(26.*Fo-Fp-Fm)/24. VEL_MEMO = VEL_MEMO + U_F(J+1+i,K,1) ENDDO VEL_MEMO = ((3.*Fo) - VEL_MEMO)/3. DO I = -1,1 U_F(J+1+i,K,1) = U_F(J+1+i,K,1) + VEL_MEMO ENDDO c-------- EasternB ---------------------- Fp = globalP1(6,INDC+1,K) Fo = globalP1(6,INDC,K) Fm = globalP1(6,INDC-1,K) VEL_MEMO = 0. DO I = -1,1 U_F(J+1+i,K,2) = ((Fp-2.*Fo+Fm)/24.)* & ((12.*float(i)**2+1.)/9.)+ & ((float(i)/6.)*(Fp-Fm))+(26.*Fo-Fp-Fm)/24. VEL_MEMO = VEL_MEMO + U_F(J+1+i,K,2) ENDDO VEL_MEMO = ((3.*Fo) - VEL_MEMO)/3. DO I = -1,1 U_F(J+1+i,K,2) = U_F(J+1+i,K,2) + VEL_MEMO ENDDO c------------------------------------------------------ ENDDO ENDDO c------------------------------------------------------- c (2.2) BiLinear for tangent velocity component (v in this case) c------------------------------------------------------- I = 1 II = WesternB V_F(:,:,:) = 0. DO K = 1,NrC DO J = 1,NyC x1 = Xv_P(II,P2C1_V(J)) x2 = Xv_P(II,P2C2_V(J)) x3 = Xv_P(II+1,P2C1_V(J)) x4 = Xv_P(II+1,P2C2_V(J)) y1 = Yv_P(II,P2C1_V(J)) y2 = Yv_P(II,P2C2_V(J)) y3 = Yv_P(II+1,P2C1_V(J)) y4 = Yv_P(II+1,P2C2_V(J)) x = Xv_F(I,J) y = Yv_F(I,J) f1 = globalP2(1,P2C1_V(J),K) f2 = globalP2(1,P2C2_V(J),K) f3 = globalP2(2,P2C1_V(J),K) f4 = globalP2(2,P2C2_V(J),K) call blint(x1,x2,x3,x4,y1,y2,y3,y4,f1,f2,f3,f4,x,y,f) V_F(J,K,1) = f ENDDO ENDDO c.............................................. I = NxC II = EasternB DO K = 1,NrC DO J = 1,NyC x1 = Xv_P(II,P2C1_V(J)) x2 = Xv_P(II,P2C2_V(J)) x3 = Xv_P(II-1,P2C1_V(J)) x4 = Xv_P(II-1,P2C2_V(J)) y1 = Yv_P(II,P2C1_V(J)) y2 = Yv_P(II,P2C2_V(J)) y3 = Yv_P(II-1,P2C1_V(J)) y4 = Yv_P(II-1,P2C2_V(J)) x = Xv_F(I,J) y = Yv_F(I,J) f1 = globalP2(6,P2C1_V(J),K) f2 = globalP2(6,P2C2_V(J),K) f3 = globalP2(5,P2C1_V(J),K) f4 = globalP2(5,P2C2_V(J),K) call blint(x1,x2,x3,x4,y1,y2,y3,y4,f1,f2,f3,f4,x,y,f) V_F(J,K,2) = f ENDDO ENDDO c------------------------------------------------------- c (2.3.1) Linear c------------------------------------------------------- c c WesternB c DO K = 1,NrP DO J = 1,NyP c DEPDEP = DEEP_P(WesternB,J,K) * DEEP_P(WesternB+1,J,K) c gammaT =(globalP3(2,J,K)-globalP3(1,J,K))*DEPDEP gammaS =(globalP4(2,J,K)-globalP4(1,J,K))*DEPDEP gammaeta =(globalP5(2,J,K)-globalP5(1,J,K))*DEPDEP cgm------- c gammaV =(globalP2(2,J,K)-globalP2(1,J,K))*DEPDEP cgm--------- globalP3(1,J,K) = globalP3(1,J,K) + terzo* gammaT globalP4(1,J,K) = globalP4(1,J,K) + terzo* gammaS globalP5(1,J,K) = globalP5(1,J,K) + terzo* gammaeta cgm------- c globalP2(1,J,K) = globalP2(1,J,K) + terzo* gammaV cgm--------- ENDDO c-------------------------------------------------------------- DO J = 1,NyP-1 IF (globalP4(1,J,K).eq.0.) THEN !uso la salinità come discriminante globalP3(1,J,K) = globalP3(1,J+1,K) globalP4(1,J,K) = globalP4(1,J+1,K) globalP5(1,J,K) = globalP5(1,J+1,K) ENDIF ENDDO DO J = NyP,2,-1 IF (globalP4(1,J,K).eq.0.) THEN !uso la salinità come discriminante globalP3(1,J,K) = globalP3(1,J-1,K) globalP4(1,J,K) = globalP4(1,J-1,K) globalP5(1,J,K) = globalP5(1,J-1,K) cgm--------- c globalP2(1,J,K) = globalP2(1,J-1,K) cgm------------- ENDIF ENDDO c--------------------------------------------------------------- ENDDO c c EasternB c DO K = 1,NrP DO J = 1,NyP c DEPDEP = DEEP_P(EasternB,J,K) * DEEP_P(EasternB-1,J,K) c gammaT =(globalP3(5,J,K)-globalP3(6,J,K))*DEPDEP gammaS =(globalP4(5,J,K)-globalP4(6,J,K))*DEPDEP gammaeta =(globalP5(5,J,K)-globalP5(6,J,K))*DEPDEP cgm------------- c gammaV =(globalP2(5,J,K)-globalP2(6,J,K))*DEPDEP cgm---------------- globalP3(6,J,K) = globalP3(6,J,K) + terzo* gammaT globalP4(6,J,K) = globalP4(6,J,K) + terzo* gammaS globalP5(6,J,K) = globalP5(6,J,K) + terzo* gammaeta cgm---------------- c globalP2(6,J,K) = globalP2(6,J,K) + terzo* gammaV cgm---------------- ENDDO c-------------------------------------------------------------- DO J = 1,NyP-1 IF (globalP4(6,J,K).eq.0.) THEN !uso la salinità come discriminante globalP3(6,J,K) = globalP3(6,J+1,K) globalP4(6,J,K) = globalP4(6,J+1,K) globalP5(6,J,K) = globalP5(6,J+1,K) cgm---------------- c globalP2(6,J,K) = globalP2(6,J+1,K) cgm---------------- ENDIF ENDDO DO J = NyP,2,-1 IF (globalP4(6,J,K).eq.0.) THEN !uso la salinità come discriminante globalP3(6,J,K) = globalP3(6,J-1,K) globalP4(6,J,K) = globalP4(6,J-1,K) globalP5(6,J,K) = globalP5(6,J-1,K) cgm---------------- c globalP2(6,J,K) = globalP2(6,J-1,K) cgm---------------- ENDIF ENDDO c--------------------------------------------------------------- ENDDO c------------------------------------------------------- c (2.3.2) Linear c------------------------------------------------------- DO J = 1,NyC c....MASK DEEP_F IF (J.EQ.NyC) THEN !EVITO ERRORI DI CHECK BOUND x Vittorio DEPDEP_F_WesternB (:) = DEEP_F(1 ,J,:)*DEEP_F(1 ,J,:) DEPDEP_F_EasternB(:) = DEEP_F(NxC,J,:)*DEEP_F(NxC,J,:) ELSE DEPDEP_F_WesternB (:) = DEEP_F(1 ,J,:)*DEEP_F(1 ,J+1,:) DEPDEP_F_EasternB(:) = DEEP_F(NxC,J,:)*DEEP_F(NxC,J+1,:) ENDIF c....WesternB...T..................... diff(:) = globalP3(1,P2C_linU(J)+1,:)- & globalP3(1,P2C_linU(J) ,:) T_F(J,:,1) = globalP3(1,P2C_linU(J) ,:)+ & (diff(:)/3.)*float((WO3_linU(J))) & *DEPDEP_F_WesternB(:) c.....EasternB..T..................... diff(:) = globalP3(6,P2C_linU(J)+1,:)- & globalP3(6,P2C_linU(J) ,:) T_F(J,:,2) = globalP3(6,P2C_linU(J) ,:)+ & (diff(:)/3.)*float((WO3_linU(J))) & *DEPDEP_F_EasternB(:) c.....WesternB...S..................... diff(:) = globalP4(1,P2C_linU(J)+1,:)- & globalP4 (1,P2C_linU(J) ,:) S_F(J,:,1) = globalP4(1,P2C_linU(J) ,:)+ & (diff(:)/3.)*float((WO3_linU(J))) & *DEPDEP_F_WesternB(:) c.... EasternB..S..................... diff(:) = (globalP4(6,P2C_linU(J)+1,:)- & globalP4 (6,P2C_linU(J) ,:)) S_F(J,:,2) = globalP4(6,P2C_linU(J) ,:)+ & (diff(:)/3.)*float((WO3_linU(J))) & *DEPDEP_F_EasternB(:) c.... WesternB...Eta..................... diff(:) = globalP5(1,P2C_linU(J)+1,:)- & globalP5(1,P2C_linU(J) ,:) eta_F(J,:,1) = globalP5(1,P2C_linU(J) ,:)+ & (diff(:)/3.)*float((WO3_linU(J))) & *DEPDEP_F_WesternB(:) c.... EasternB..Eta..................... diff(:) = globalP5(6,P2C_linU(J)+1,:)- & globalP5(6,P2C_linU(J) ,:) eta_F(J,:,2) = globalP5(6,P2C_linU(J) ,:)+ & (diff(:)/3.)*float(WO3_linU(J)) & *DEPDEP_F_EasternB(:) cgm-------------------------- c.... WesternB...V..................... c diff(:) = globalP2(1,P2C_linU(J)+1,:)- c & globalP2(1,P2C_linU(J) ,:) c c V_F(J,:,1) = globalP2(1,P2C_linU(J) ,:)+ c & (diff(:)/3.)*float((WO3_linU(J))) c & *DEPDEP_F_WesternB(:) c.... EasternB..V..................... c diff(:) = globalP2(6,P2C_linU(J)+1,:)- c & globalP2(6,P2C_linU(J) ,:) c c V_F(J,:,2) = globalP2(6,P2C_linU(J) ,:)+ c & (diff(:)/3.)*float((WO3_linU(J))) c & *DEPDEP_F_EasternB(:) cgm---------------------------- ENDDO 8765 CONTINUE c=============== END INTERPOLAZIONI x CHILD ==================== RETURN END c================================================================ SUBROUTINE BLINT(x1,x2,x3,x4,y1,y2,y3,y4,f1,f2,f3,f4,x,y,f) c================================================================ C C Bilinear interpolation subroutine. C (Xi,Yi,fi) = data grid & values surounding model point (x,y) C f = interpolated value at the model grid point. IMPLICIT NONE real*4 a1,a2,a3,a4 real*4 b1,b2,b3,b4 real*4 x1,x2,x3,x4 real*4 y1,y2,y3,y4 real*4 f1,f2,f3,f4 real*4 x,y,A,B,C,t,f,s C a1=x1-x2+x3-x4 a2=-x1+x4 a3=-x1+x2 a4=x1-x b1=y1-y2+y3-y4 b2=-y1+y4 b3=-y1+y2 b4=y1-y A=a3*b1-a1*b3 B=b2*a3+b1*a4-a1*b4-a2*b3 C=-a2*b4+a4*b2 if(ABS(A*C).gt.0.002*B**2) then t=(-B-sqrt(B*B-4.*A*C))/(2.*A) else t=C/ABS(B) endif 10 CONTINUE A=a2*b1-a1*b2 B=b3*a2+b1*a4-a1*b4-a3*b2 C=-a3*b4+a4*b3 if(ABS(A*C).gt.0.002*B**2) then s=(-B+sqrt(B*B-4.*A*C))/(2.*A) else s=-C/ABS(B) endif 20 CONTINUE f=f1*(1.-t)*(1.-s)+f2*t*(1.-s)+f3*s*t+f4*(1.-t)*s return end