C $Header: /home/ubuntu/mnt/e9_copy/MITgcm_contrib/nesting_sannino/nest_driver/main.F,v 1.1 2010/11/28 03:27:56 jmc Exp $ C $Name: $ #include "CPP_OPTIONS.h" C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| 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,nSyC :: Domain decomposition used for CHILD C nSxP,nSyP :: Domain decomposition used for PARENT C OLX,OLY :: Domain dec. overlapping (same for both models) C NrP,NyP,NxP :: Dimension PARENT model C NrC,NyC,NxC :: 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 n3dC :: number of 3-D fields sent from CHILD 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) INTEGER :: P2C_linU(NyC) !Linear interp. INTEGER :: WO3_linU(NyC) !Linear interp. !Which Of 3 INTEGER :: P2C_linV(NyC) !Linear interp. INTEGER :: WO3_linV(NyC) !Linear interp. !Which Of 3 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 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) REAL*8,allocatable :: VAR_C1(:,:,:,:) 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) REAL*8,allocatable :: localP3D_a(:,:,:), localP2D_a(:,:) 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) INTEGER :: index C---------------------------------------------------- C Define Global Variables to Exchange C---------------------------------------------------- INTEGER :: n3dC PARAMETER ( n3dC = 15 ) REAL*8 :: globalC3D(NxC,NyC,NrC,n3dC) C |___________ 15 fields REAL*8 globalC2D(NxC,NyC,4) C |___________ 4 fields REAL*8,allocatable :: globalC3D_a(:,:,:,:),globalC2D_a(:,:,:) INTEGER :: indexF,index2F,index3F REAL*4 :: AREA_VOL INTEGER :: vstart,vstop,VCONT,VCONTP(0:3) C- log-file IO-unit and name: STDlog.xxxx INTEGER iUnit PARAMETER ( iUnit = 35 ) CHARACTER*11 fNam INTEGER mLoop 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) CALL MPI_Comm_size(NEST_COMM,isize,ierr) CALL MPI_Comm_rank(NEST_COMM,irank,ierr) C-------------------------------------------------------- C- change local dir to rank_N and open log file CALL SETDIR( rank ) WRITE(fNam,'(A,I4.4)') 'STDlog.', rank OPEN( iUnit, FILE=fNam, STATUS='unknown') mLoop = 0 C-------------------------------------------------------- C COMPUTE MASTER VALUES C-------------------------------------------------------- MSTR_DRV(1) = 0 MSTR_PRNT(1) = npd MSTR_CHLD(1) = NCPUs_PRNT + npd DO Count_Lev = 2, NST_LEV_TOT MSTR_DRV(Count_Lev) = MSTR_CHLD(Count_Lev-1) + & NCPUs_CHLD(Count_Lev - 1) MSTR_CHLD(Count_Lev) = MSTR_DRV(Count_Lev) + 1 MSTR_PRNT(Count_Lev) = MSTR_CHLD(Count_Lev-1) ENDDO vstart = 1+rank*(nSxP/MSTR_PRNT(1)) vstop = (1+rank)*(nSxP/MSTR_PRNT(1)) VCONT = (nSxP/npd)*(nSyP)*rank-1 VCONTP(rank) = VCONT C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| C-- Print out nesting parameter: WRITE(iUnit,'(A)') '// ===================================' WRITE(iUnit,'(A)') '// NEST_DRIVER parameters :' WRITE(iUnit,'(A)') '// ===================================' WRITE(iUnit,*) 'NEST_DRIVER: rank =', rank, ' ; color =',color WRITE(iUnit,*) 'NEST_DRIVER: size =', size, ' ; npd =', npd WRITE(iUnit,*) 'NEST_DRIVER: irank =', irank, ' ; isize =', isize WRITE(iUnit,*) 'NEST_DRIVER: vstart =', vstart WRITE(iUnit,*) 'NEST_DRIVER: vstop =', vstop WRITE(iUnit,*) 'NEST_DRIVER: VCONTP =', VCONTP(rank) WRITE(iUnit,*) 'NEST_DRIVER: NST_LEV_TOT =', NST_LEV_TOT WRITE(iUnit,*) 'NEST_DRIVER: NST_LEV =', NST_LEV WRITE(iUnit,*) 'NEST_DRIVER: NCPUs_PRNT =', NCPUs_PRNT WRITE(iUnit,*) 'NEST_DRIVER: NCPUs_CHLD =', NCPUs_CHLD WRITE(iUnit,*) 'NEST_DRIVER: MSTR_DRV =', MSTR_DRV WRITE(iUnit,*) 'NEST_DRIVER: MSTR_PRNT =', MSTR_PRNT WRITE(iUnit,*) 'NEST_DRIVER: MSTR_CHLD =', MSTR_CHLD C-------------------------------------------------------- C COMPUTE DOMAIN DECOMPOSITION C-------------------------------------------------------- c if(rank.eq.0) then IM_C = int(NxC/nSxC) JM_C = int(NyC/nSyC) IM_P = int(NxP/nSxP) JM_P = int(NyP/nSyP) 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 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 index = 6*JM_P*NrP*5 index_var3D = IM_P*JM_P*NrP index_var2D = IM_P*JM_P indexF = (JM_C+OLY+OLY)*NrC*2*5 index3F = IM_C*JM_C*NrC*n3dC index2F = IM_C*JM_C*4 ALLOCATE( globalPA(6,JM_P,NrP,5) ) ALLOCATE( localP3D_a(IM_P,JM_P,NrP) ) ALLOCATE( localP2D_a(IM_P,JM_P) ) ALLOCATE( VAR_C1((JM_C+OLY+OLY),NrC,2,5) ) ALLOCATE( globalC3D_a(IM_C,JM_C,NrC,n3dC) ) ALLOCATE( globalC2D_a(IM_C,JM_C,4)) IF ( rank.EQ.0 ) THEN C-------------------------------------------------------- C WARNING C-------------------------------------------------------- write(iUnit,*) '*************************************' write(iUnit,*) ' have you update geometric files?' write(iUnit,*) ' in ./CHILD e ./PARENT' write(iUnit,*) '*************************************' C-------------------------------------------------------- C PARENT MODEL C-------------------------------------------------------- write(iUnit,*) ' [1] Read PARENT model geometry' C---------------------------------------------------- C XC & YC C---------------------------------------------------- MSIZE = NxP*NyP*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/XC.data', & form='unformatted') read (1,REC=1) Xo_P(:,:) close(1) open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/YC.data', & form='unformatted') read (1,REC=1) Yo_P(:,:) close(1) C---------------------------------------------------- C XG & YG C---------------------------------------------------- MSIZE = NxP*NyP*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/XG.data', & form='unformatted') read (1,REC=1) Xg_P(:,:) close(1) open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/YG.data', & form='unformatted') 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*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/hFacC.data', & form='unformatted') read (1,REC=1) hFacC_P(:,:,:) close(1) C---------------------------------------------------- C hFacW C---------------------------------------------------- MSIZE = NxP*NyP*NrP*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/hFacW.data', & form='unformatted') read (1,REC=1) hFacW_P(:,:,:) close(1) C---------------------------------------------------- C hFacS C---------------------------------------------------- MSIZE = NxP*NyP*NrP*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/hFacS.data', & form='unformatted') read (1,REC=1) hFacS_P(:,:,:) close(1) C---------------------------------------------------- C RAC C---------------------------------------------------- MSIZE = NxP*NyP*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/RAC.data', & form='unformatted') read (1,REC=1) RAC_P(:,:) close(1) C---------------------------------------------------- C RAW C---------------------------------------------------- MSIZE = NxP*NyP*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/PARENT/RAW.data', & form='unformatted') read (1,REC=1) RAW_P(:,:) close(1) C---------------------------------------------------- C RAS C---------------------------------------------------- MSIZE = NxP*NyP*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & 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-------------------------------------------------------- write(iUnit,*) ' [2] Read CHILD model geometry' C---------------------------------------------------- C XC & YC C---------------------------------------------------- MSIZE = NxC*NyC*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/CHILD/XC.data', & form='unformatted') read (1,REC=1) Xo_C(:,:) close(1) open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/CHILD/YC.data', & form='unformatted') read (1,REC=1) Yo_C(:,:) close(1) C---------------------------------------------------- C XG & YG C---------------------------------------------------- MSIZE = NxC*NyC*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/CHILD/XG.data', & form='unformatted') read (1,REC=1) Xg_C(:,:) close(1) open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & 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*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/CHILD/hFacC.data', & form='unformatted') read (1,REC=1) hFacC_C(:,:,:) close(1) C---------------------------------------------------- C hFacW C---------------------------------------------------- MSIZE = NxC*NyC*NrC*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/CHILD/hFacW.data', & form='unformatted') read (1,REC=1) hFacW_C(:,:,:) close(1) C---------------------------------------------------- C hFacC C---------------------------------------------------- MSIZE = NxC*NyC*NrC*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/CHILD/hFacS.data', & form='unformatted') read (1,REC=1) hFacS_C(:,:,:) close(1) C---------------------------------------------------- C RAC C---------------------------------------------------- MSIZE = NxC*NyC*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/CHILD/RAC.data', & form='unformatted') read (1,REC=1) RAC_C(:,:) close(1) C---------------------------------------------------- C RAW C---------------------------------------------------- MSIZE = NxC*NyC*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & file=trim(dirNEST)//'/CHILD/RAW.data', & form='unformatted') read (1,REC=1) RAW_C(:,:) close(1) C---------------------------------------------------- C RAS C---------------------------------------------------- MSIZE = NxC*NyC*WORDLENGTH open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD', & 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 > o > | | | C |__/__|_____| C | | | C > o > | C |_____|_____|_____|_____| C C---------------------------------------------------- write(iUnit,*) ' [3] Compute J index P-->C' C-------------------------------------------------------- C Compute J indicies for mapping P->C (I) C-------------------------------------------------------- I = 1 II = WesternB 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 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 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 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 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-------------------------------------------------------- write(iUnit,*) ' [5] Compute J index P-->C for (o)' I = 1 II = WesternB 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-------------------------------------------------------- write(iUnit,*) ' [6] Compute J index P-->C for (v bilinear)' I = 1 II = WesternB 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) 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-------------------------------------------------------- write(iUnit,*) ' [8] Compute I J index C-->P for (o)' DO J = 1,NyP DO I = 1,NxP I_C2P(:,I,J) = 0 J_C2P(:,I,J) = 0 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 ENDDO ENDDO C-- end if rank=0 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) 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) 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) 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 C- WesternB side DO II = 1,9 IF (I_C2P(II,I,J).eq.0) cycle IF (J_C2P(II,I,J).eq.0) cycle 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 mLoop = mLoop + 1 c DO mLoop=1,24 WRITE(iUnit,'(A,I6)') '== Main Loop , iter=', mLoop 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) 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 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 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. J1 = 1+IndJ-OLY J2 = JM_C+IndJ+OLY JJ1 = 1 JJ2 = JM_C+OLY+OLY IF(1 +IndJ-OLY.LT.0) THEN J1 = 1 JJ1 = 4 ENDIF IF(JM_C+IndJ+OLY.GT.NyC) THEN J2 = NyC JJ2 = JM_C ENDIF 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 write(*,*) 'VIC: MANDO SEGNALE DI OK AL CHILD PER INIZIALIZZARE' ONOFF=1 ENDIF c 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) VAR_C1(:,:,:,:) = 0. J1 = 1+IndJ-OLY J2 = JM_C+IndJ+OLY JJ1 = 1 JJ2 = JM_C+OLY+OLY IF(1 +IndJ-OLY.LT.0) THEN J1 = 1 JJ1 = 4 ENDIF IF(JM_C+IndJ+OLY.GT.NyC) THEN J2 = NyC JJ2 = JM_C ENDIF 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--------------------------------------------------------------------- goto 8888 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. J1 = 1+IndJ-OLY J2 = JM_C+IndJ+OLY JJ1 = 1 JJ2 = JM_C+OLY+OLY IF(1 +IndJ-OLY.LT.0) THEN J1 = 1 JJ1 = 4 ENDIF IF(JM_C+IndJ+OLY.GT.NyC) THEN J2 = NyC JJ2 = JM_C ENDIF 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 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. J1 = 1+IndJ-OLY J2 = JM_C+IndJ+OLY JJ1 = 1 JJ2 = JM_C+OLY+OLY IF(1 +IndJ-OLY.LT.0) THEN J1 = 1 JJ1 = 4 ENDIF IF(JM_C+IndJ+OLY.GT.NyC) THEN J2 = NyC JJ2 = JM_C ENDIF 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) write(iUnit,*) & 'CALL MPI_RECV 3-D var from CHILD, index3F=', index3F from= MPI_ANY_SOURCE CALL MPI_RECV (globalC3D_a,index3F, 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,index2F, 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 C-- end if rank=0 ENDIF CALL MPI_BCAST(globalC3D,NxC*NyC*NrC*n3dC,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 Indx = I_C2P(II,I,J) Jndx = J_C2P(II,I,J) 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 Indx = I_C2P(II,I,J) Jndx = J_C2P(II,I,J) 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 write(iUnit,*) ' (*) 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 c CALL MPI_SEND (VAR3D_P(1+IndI:IM_P+IndI c & ,1+IndJ:JM_P+IndJ,:,iVar), c & index_var3D,MPI_REAL8,MSTR_PRNT(NST_LEV)+VCONT, c & 3000,MPI_Comm_World,ierr) localP3D_a(:,:,:) = & VAR3D_P(1+IndI:IM_P+IndI,1+IndJ:JM_P+IndJ,:,iVar) CALL MPI_SEND ( localP3D_a, index_var3D, MPI_REAL8, & MSTR_PRNT(NST_LEV)+VCONT, & 3000, MPI_Comm_World, ierr ) ENDDO DO iVar=1,4 c CALL MPI_SEND (VAR2D_P(1+IndI:IM_P+IndI c & ,1+IndJ:JM_P+IndJ,iVar), c & index_var2D,MPI_REAL8,MSTR_PRNT(NST_LEV)+VCONT, c & 3000,MPI_Comm_World,ierr) localP2D_a(:,:) = & VAR2D_P(1+IndI:IM_P+IndI,1+IndJ:JM_P+IndJ,iVar) CALL MPI_SEND ( localP2D_a, index_var2D, MPI_REAL8, & MSTR_PRNT(NST_LEV)+VCONT, & 3000, MPI_Comm_World, ierr ) ENDDO C----------------------------------------------------------- ENDDO ENDDO CALL MPI_BARRIER(NEST_COMM,ierr) ENDDO C--------------------------------------------------------- C======================================================= C END MAIN LOOP C======================================================= CLOSE( iUnit ) CALL MPI_FINALIZE(ierr) C--------------------------------------------------------- 101 FORMAT (I1) STOP END