MODULE W3PARALL !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Parallel routines for implicit solver ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 ! !/PDLIB INTEGER :: PDLIB_NSEAL, PDLIB_NSEALM !/PDLIB INTEGER, ALLOCATABLE :: JX_TO_JSEA(:), ISEA_TO_JSEA(:) REAL, ALLOCATABLE :: AC_tot(:,:) INTEGER, ALLOCATABLE :: ListISPnextDir(:), ListISPprevDir(:) INTEGER, ALLOCATABLE :: ListISPnextFreq(:), ListISPprevFreq(:) INTEGER, PARAMETER :: IMEM = 2 REAL, PARAMETER :: ONESIXTH = 1.0d0/6.0d0 REAL, PARAMETER :: ONETHIRD = 1.0d0/3.0d0 REAL, PARAMETER :: ZERO = 0.0d0 REAL, PARAMETER :: THR8 = TINY(1.0) REAL, PARAMETER :: THR = TINY(1.0) !!/S CALL STRACE (IENT, 'W3XXXX') CONTAINS !/ ------------------------------------------------------------------- / SUBROUTINE WAV_MY_WTIME(eTime) !/ ------------------------------------------------------------------- / !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ IMPLICIT NONE !/S INTEGER, SAVE :: IENT = 0 INTEGER mpimode REAL(8), intent(out) :: eTime !/MPI REAL(8) mpi_wtime mpimode=0 !/MPI mpimode=1 !/MPI eTime=mpi_wtime() !/S CALL STRACE (IENT, 'WAV_MY_WTIME') IF (mpimode .eq. 0) THEN CALL CPU_TIME(eTime) END IF !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / SUBROUTINE PRINT_MY_TIME(string) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Print timings ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE USE W3ODATMD, ONLY : IAPROC IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 !/ !/ ------------------------------------------------------------------- / ! character(*), intent(in) :: string REAL(8) :: eTime !/S CALL STRACE (IENT, 'PRINT_MY_TIME') CALL WAV_MY_WTIME(eTime) WRITE(740+IAPROC,*) 'TIMING time=', eTime, ' at step ', string !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / SUBROUTINE PROP_REFRACTION_PR1(ISEA,DTG, CAD) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Compute refraction part in matrix ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE USE W3GDATMD, ONLY: NK, NK2, NTH, NSPEC, SIG, DSIP, ECOS, ESIN, & EC2, ESC, ES2, FACHFA, MAPWN, FLCTH, FLCK, & CTMAX, DMIN, DTH, CTHG0S, MAPSF USE W3ADATMD, ONLY: CG, WN, DCXDX, DCXDY, DCYDX, DCYDY, DDDX, & DDDY, DW !/REFRX USE W3ADATMD, ONLY: DCDX, DCDY USE W3IDATMD, ONLY: FLCUR USE W3ODATMD, only : IAPROC IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 !/ !/ ------------------------------------------------------------------- / !/ !/ REAL, intent(out) :: CAD(NSPEC) INTEGER, intent(in) :: ISEA REAL, intent(in) :: DTG INTEGER :: ISP, IK, ITH, IX, IY REAL :: FRK(NK), FRG(NK), DSDD(0:NK+1) REAL :: FACTH, DCXY, DCYX, DCXXYY, DTTST REAL :: eDCXDX, eDCXDY, eDCYDX, eDCYDY, eDDDX, eDDDY, eCTHG0 REAL :: VCFLT(NSPEC), DEPTH, FDG REAL :: FDDMAX !/S CALL STRACE (IENT, 'PROP_REFRACTION_PR1') IX=MAPSF(ISEA,1) IY=MAPSF(ISEA,2) eDDDX=DDDX(IY,IX) eDDDY=DDDY(IY,IX) eCTHG0 = CTHG0S(ISEA) FACTH = DTG / DTH ! FDG = FACTH * eCTHG0 DEPTH = MAX ( DMIN , DW(ISEA) ) DO IK=0, NK+1 IF ( DEPTH*WN(IK,ISEA) .LT. 5. ) THEN DSDD(IK) = MAX ( 0. , CG(IK,ISEA)*WN(IK,ISEA)-0.5*SIG(IK) ) / DEPTH ELSE DSDD(IK) = 0. END IF END DO FDDMAX=0 DO ITH=1, NTH FDDMAX = MAX ( FDDMAX , ABS(ESIN(ITH)*eDDDX - ECOS(ITH)*eDDDY ) ) END DO !/DEBUG WRITE(740+IAPROC,*) 'eDDDX=', eDDDX, ' Y=', eDDDY !/DEBUG WRITE(740+IAPROC,*) 'FDDMAX=', FDDMAX !/DEBUG FLUSH(740+IAPROC) DO IK=1, NK FRK(IK) = FACTH * DSDD(IK) / WN(IK,ISEA) !FRK(IK) = FRK(IK) / MAX ( 1. , FRK(IK)*FDDMAX/CTMAX ) FRG(IK) = FDG * CG(IK,ISEA) END DO DO ISP=1, NSPEC VCFLT(ISP) = FRG(MAPWN(ISP)) * ECOS(ISP) + & FRK(MAPWN(ISP)) * ( ESIN(ISP)*eDDDX - ECOS(ISP)*eDDDY ) END DO !/DEBUG WRITE(740+IAPROC,*) 'pdlib: FACTH=', FACTH !/DEBUG WRITE(740+IAPROC,*) 'pdlib: CTHG0=', eCTHG0 !/DEBUG WRITE(740+IAPROC,*) 'pdlib: FDG=', FDG !/DEBUG WRITE(740+IAPROC,*) 'pdlib: FDDMAX=', FDDMAX !/DEBUG WRITE(740+IAPROC,*) 'pdlib: sum(FRK)=', sum(FRK) !/DEBUG WRITE(740+IAPROC,*) 'pdlib: sum(FRG)=', sum(FRG) !/DEBUG WRITE(740+IAPROC,*) 'pdlib: sum(DSDD)=', sum(DSDD) !/DEBUG WRITE(740+IAPROC,*) 'ISEA=', ISEA, ' sum(VCTH)=', sum(VCFLT) !/DEBUG FLUSH(740+IAPROC) ! !/REFRX! 3.c @C/@x refraction and great-circle propagation !/REFRX VCFLT = 0. !/REFRX FRK = 0. !/REFRX DO IK=1, NK !/REFRX FRK(IK) = FACTH * CG(IK,ISEA) * WN(IK,ISEA) / SIG(IK) !/REFRX END DO !/REFRX DO ISP=1, NSPEC !/REFRX VCFLT(ISP) = FRG(MAPWN(ISP)) * ECOS(ISP) & !/REFRX + FRK(MAPWN(ISP)) * ( ESIN(ISP)*DCDX(ISP,1,MAPWN(ISP)) & !/REFRX - ECOS(ISP)*DCDY(ISP,1,MAPWN(ISP)) ) !/REFRX END DO ! IF ( FLCUR ) THEN eDCXDX=DCXDX(IY,IX) eDCXDY=DCXDY(IY,IX) eDCYDX=DCYDX(IY,IX) eDCYDY=DCYDY(IY,IX) DCYX = FACTH * eDCYDX DCXXYY = FACTH * ( eDCXDX - eDCYDY ) DCXY = FACTH * eDCXDY DO ISP=1, NSPEC VCFLT(ISP) = VCFLT(ISP) + ES2(ISP)*DCYX + ESC(ISP)*DCXXYY - EC2(ISP)*DCXY END DO END IF DO ISP=1,NSPEC CAD(ISP)=DBLE(VCFLT(ISP)) END DO !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / ! SUBROUTINE PROP_REFRACTION_PR3(IP, ISEA, DTG, CAD, DoLimiter) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Compute refraction part in matrix alternative approach ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE USE CONSTANTS, ONLY : LPDLIB USE W3GDATMD, ONLY: NK, NK2, NTH, NSPEC, SIG, DSIP, ECOS, ESIN, & EC2, ESC, ES2, FACHFA, MAPWN, FLCTH, FLCK, & CTMAX, DMIN, DTH, CTHG0S, MAPSF USE W3ADATMD, ONLY: CG, WN, DCXDX, DCXDY, DCYDX, DCYDY, DDDX, & DDDY, DW USE W3IDATMD, ONLY: FLCUR USE W3ODATMD, only : IAPROC IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 REAL, intent(out) :: CAD(NSPEC) INTEGER, intent(in) :: ISEA, IP REAL, intent(in) :: DTG logical, intent(in) :: DoLimiter INTEGER :: ISP, IK, ITH, IX, IY REAL :: FRK(NK), FRG(NK), DSDD(0:NK+1) REAL :: FACTH, DCXY, DCYX, DCXXYY, DTTST REAL :: eDCXDX, eDCXDY, eDCYDX, eDCYDY, eDDDX, eDDDY, eCTHG0 REAL :: VCFLT(NSPEC), DEPTH, FDG REAL :: FDDMAX, CFLTHMAX, VELNOFILT, CTMAX_eff !/S CALL STRACE (IENT, 'PROP_REFRACTION_PR3') IX=MAPSF(ISEA,1) IY=MAPSF(ISEA,2) IF (LPDLIB) THEN eDDDX=DDDX(1,IP) eDDDY=DDDY(1,IP) ELSE eDDDX=DDDX(IY,IX) eDDDY=DDDY(IY,IX) ENDIF eCTHG0 = CTHG0S(ISEA) FACTH = 1.0 / DTH ! FDG = FACTH * eCTHG0 DEPTH = MAX ( DMIN , DW(ISEA) ) DO IK=0, NK+1 IF ( DEPTH*WN(IK,ISEA) .LT. 5. ) THEN DSDD(IK) = MAX ( 0. , CG(IK,ISEA)*WN(IK,ISEA)-0.5*SIG(IK) ) / DEPTH ELSE DSDD(IK) = 0. END IF END DO DO IK=1, NK FRK(IK) = FACTH * DSDD(IK) / WN(IK,ISEA) FRG(IK) = FDG * CG(IK,ISEA) END DO IF (FLCUR) THEN IF (LPDLIB) THEN eDCXDX = DCXDX(1,IP) eDCXDY = DCXDY(1,IP) eDCYDX = DCYDX(1,IP) eDCYDY = DCYDY(1,IP) ELSE eDCXDX = DCXDX(IY,IX) eDCXDY = DCXDY(IY,IX) eDCYDX = DCYDX(IY,IX) eDCYDY = DCYDY(IY,IX) ENDIF DCYX = FACTH * eDCYDX DCXXYY = FACTH * ( eDCXDX - eDCYDY ) DCXY = FACTH * eDCXDY DO ISP=1, NSPEC VCFLT(ISP) = ES2(ISP)*DCYX + ESC(ISP)*DCXXYY - EC2(ISP)*DCXY END DO ELSE VCFLT=0 END IF ! !/REFRX! 3.c @C/@x refraction and great-circle propagation !/REFRX DO IK=1, NK !/REFRX FRK(IK) = FACTH * CG(IK,ISEA) * WN(IK,ISEA) / SIG(IK) !/REFRX END DO ! CTMAX_eff=CTMAX/DTG DO ISP=1, NSPEC VELNOFILT = VCFLT(ISP) & + FRG(MAPWN(ISP)) * ECOS(ISP) & + FRK(MAPWN(ISP)) * (ESIN(ISP)*eDDDX - ECOS(ISP)*eDDDY) ! ! Puts filtering on total velocity (including currents and great circle effects) ! the filtering limits VCFLT to be less than CTMAX ! this modification was proposed by F. Ardhuin 2011/03/06 ! IF (DoLimiter .eqv. .TRUE.) THEN VCFLT(ISP)=SIGN(MIN(ABS(VELNOFILT),CTMAX_eff),VELNOFILT) ELSE VCFLT(ISP)=VELNOFILT END IF END DO DO ISP=1,NSPEC CAD(ISP)=DBLE(VCFLT(ISP)) END DO !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / SUBROUTINE PROP_FREQ_SHIFT(IP, ISEA, CAS, DMM, DTG) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Compute freq. shift in matrix ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE USE CONSTANTS, ONLY : LPDLIB USE W3GDATMD, ONLY: NK, NK2, NTH, NSPEC, SIG, DSIP, ECOS, ESIN, & EC2, ESC, ES2, FACHFA, MAPWN, FLCTH, FLCK, & CTMAX, DMIN, DTH, MAPSF USE W3ADATMD, ONLY: CG, WN, DCXDX, DCXDY, DCYDX, DCYDY, CX, CY, DDDX, DDDY, DW USE W3ODATMD, only : IAPROC IMPLICIT NONE !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 INTEGER, intent(in) :: ISEA, IP REAL, intent(out) :: DMM(0:NK2) REAL, intent(in) :: DTG REAL, intent(out) :: CAS(NSPEC) REAL :: DB(NK2), DSDD(0:NK+1) REAL :: eDCXDX, eDCXDY, eDCYDX, eDCYDY, eCX, eCY, eDDDX, EDDDY REAL :: DCXX, DCXYYX, DCYY, FKD, FACK REAL :: VELNOFILT, VELFAC, DEPTH REAL :: CFLK(NK2,NTH), FKC(NTH), FKD0 INTEGER :: IK, ITH, ISP, IY, IX !/S CALL STRACE (IENT, 'PROP_FREQ_SHIFT') ! IF (LPDLIB) THEN eDCXDX = DCXDX(1,IP) eDCXDY = DCXDY(1,IP) eDCYDX = DCYDX(1,IP) eDCYDY = DCYDY(1,IP) eDDDX = DDDX(1,IP) eDDDY = DDDY(1,IP) ELSE IX=MAPSF(ISEA,1) IY=MAPSF(ISEA,2) eDCXDX=DCXDX(IY,IX) eDCXDY=DCXDY(IY,IX) eDCYDX=DCYDX(IY,IX) eDCYDY=DCYDY(IY,IX) eDDDX=DDDX(IY,IX) eDDDY=DDDY(IY,IX) ENDIF eCX=CX(ISEA) eCY=CY(ISEA) DCXX = - eDCXDX DCXYYX = - ( eDCXDY + eDCYDX ) DCYY = - eDCYDY FKD = ( eCX*eDDDX + eCY*eDDDY ) FACK = DTG !/DEBUG WRITE(740+IAPROC,*) 'DCXX=', DCXX, ' DCXYYX=', DCXYYX !/DEBUG WRITE(740+IAPROC,*) 'DCYY=', DCYY, ' FKD=', FKD !/DEBUG WRITE(740+IAPROC,*) 'DTG=', DTG !/DEBUG FLUSH(740+IAPROC) DO ITH=1, NTH FKC(ITH) = EC2(ITH)*DCXX + ESC(ITH)*DCXYYX + ES2(ITH)*DCYY END DO DO IK=0, NK DB(IK+1) = DSIP(IK) / CG(IK,ISEA) DMM(IK+1) = DBLE(WN(IK+1,ISEA) - WN(IK,ISEA)) END DO DB(NK+2) = DSIP(NK+1) / CG(NK+1,ISEA) DMM(NK+2) = ZERO DMM(0)=DMM(1) ! DEPTH = MAX ( DMIN , DW(ISEA) ) DO IK=0, NK+1 IF ( DEPTH*WN(IK,ISEA) .LT. 5. ) THEN DSDD(IK) = MAX ( 0. , CG(IK,ISEA)*WN(IK,ISEA)-0.5*SIG(IK) ) / DEPTH ELSE DSDD(IK) = 0. END IF END DO !/DEBUG WRITE(740+IAPROC,*) 'DSDD(min/max)=', minval(DSDD), maxval(DSDD) !/DEBUG FLUSH(740+IAPROC) DO IK=0, NK+1 FKD0 = FKD / CG(IK,ISEA) * DSDD(IK) VELFAC = FACK/DB(IK+1) DO ITH=1, NTH VELNOFILT = ( FKD0 + WN(IK,ISEA)*FKC(ITH) ) * VELFAC CFLK(IK+1,ITH) = VELNOFILT/VELFAC END DO END DO !/DEBUG WRITE(740+IAPROC,*) 'sum(CFLK)=', sum(CFLK) !/DEBUG FLUSH(740+IAPROC) DO IK=1,NK DO ITH=1,NTH ISP=ITH + (IK-1)*NTH CAS(ISP)=DBLE(CFLK(IK,ITH)) END DO END DO !/DEBUG WRITE(740+IAPROC,*) 'sum(abs(CAS))=', sum(abs(CAS)) !/DEBUG FLUSH(740+IAPROC) !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / SUBROUTINE PROP_FREQ_SHIFT_M2(IP, ISEA, CWNB_M2, DWNI_M2, DTG) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Compute freq. shift alternative approach ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / ! !/S USE W3SERVMD, ONLY: STRACE USE CONSTANTS, ONLY : LPDLIB USE W3GDATMD, ONLY: NK, NK2, NTH, NSPEC, SIG, DSIP, ECOS, ESIN, & EC2, ESC, ES2, FACHFA, MAPWN, FLCTH, FLCK, & CTMAX, DMIN, DTH, MAPSF USE W3ADATMD, ONLY: CG, WN, DCXDX, DCXDY, DCYDX, DCYDY, CX, CY, DDDX, DDDY, DW USE W3ODATMD, only : IAPROC IMPLICIT NONE !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 INTEGER, intent(in) :: ISEA, IP REAL, intent(out) :: CWNB_M2(1-NTH:NSPEC) REAL, intent(out) :: DWNI_M2(NK) REAL, intent(in) :: DTG ! REAL :: eDCXDX, eDCXDY, eDCYDX, eDCYDY, eCX, eCY, eDDDX, EDDDY REAL :: DCXX, DCXYYX, DCYY, FKD, FACK REAL :: DEPTH REAL :: FKC(NTH), FKD0 REAL :: VCWN(1-NTH:NSPEC+NTH) REAL :: DSDD(0:NK+1) REAL :: sumDiff, sumDiff1, sumDiff2, sumDiff3 REAL :: sumDiff0, sumDiff4, sumDiff5 INTEGER :: IK, ITH, ISP, IY, IX !/ ------------------------------------------------------------------- / !/S CALL STRACE (IENT, 'PROP_FREQ_SHIFT_M2') !/DEBUGDCXDX WRITE(740+IAPROC,*) 'Now we use DCXDX array in PROP_FREQ_SHIFT_M2' IF (LPDLIB) THEN eDCXDX = DCXDX(1,IP) eDCXDY = DCXDY(1,IP) eDCYDX = DCYDX(1,IP) eDCYDY = DCYDY(1,IP) eDDDX = DDDX(1,IP) eDDDY = DDDY(1,IP) ELSE IX=MAPSF(ISEA,1) IY=MAPSF(ISEA,2) eDCXDX=DCXDX(IY,IX) eDCXDY=DCXDY(IY,IX) eDCYDX=DCYDX(IY,IX) eDCYDY=DCYDY(IY,IX) eDDDX=DDDX(IY,IX) eDDDY=DDDY(IY,IX) ENDIF eCX = CX(ISEA) eCY = CY(ISEA) FACK = DTG DCXX = - FACK * eDCXDX DCXYYX = - FACK * ( eDCXDY + eDCYDX ) DCYY = - FACK * eDCYDY FKD = FACK * ( eCX*eDDDX + eCY*eDDDY ) !/DEBUGDCXDX sumDiff0=0 !/DEBUGDCXDX sumDiff1=0 !/DEBUGDCXDX sumDiff2=0 !/DEBUGDCXDX sumDiff3=0 !/DEBUGDCXDX sumDiff4=0 !/DEBUGDCXDX sumDiff5=0 DO ITH=1, NTH FKC(ITH) = EC2(ITH)*DCXX + ESC(ITH)*DCXYYX + ES2(ITH)*DCYY !/DEBUGDCXDX sumDiff0 = sumDiff0 + MIN(EC2(ITH), ZERO) !/DEBUGDCXDX sumDiff1 = sumDiff1 + MIN(DCXX, ZERO) !/DEBUGDCXDX sumDiff2 = sumDiff2 + MIN(ESC(ITH), ZERO) !/DEBUGDCXDX sumDiff3 = sumDiff3 + MIN(DCXYYX, ZERO) !/DEBUGDCXDX sumDiff4 = sumDiff4 + MIN(ES2(ITH), ZERO) !/DEBUGDCXDX sumDiff5 = sumDiff5 + MIN(DCYY, ZERO) END DO !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff0=', sumDiff0 !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff1=', sumDiff1 !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff2=', sumDiff2 !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff3=', sumDiff3 !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff4=', sumDiff4 !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff5=', sumDiff5 ! DEPTH = MAX ( DMIN , DW(ISEA) ) DO IK=0, NK+1 IF ( DEPTH*WN(IK,ISEA) .LT. 5. ) THEN DSDD(IK) = MAX ( 0. , CG(IK,ISEA)*WN(IK,ISEA)-0.5*SIG(IK) ) / DEPTH ELSE DSDD(IK) = 0. END IF END DO ISP = -NTH !/DEBUGDCXDX sumDiff=0 !/DEBUGDCXDX sumDiff1=0 !/DEBUGDCXDX sumDiff2=0 !/DEBUGDCXDX sumDiff3=0 DO IK=0, NK+1 FKD0 = FKD / CG(IK,ISEA) * DSDD(IK) DO ITH=1, NTH ISP = ISP + 1 VCWN(ISP) = FKD0 + WN(IK,ISEA)*FKC(ITH) !/DEBUGDCXDX sumDiff = sumDiff + MAX(VCWN(ISP),ZERO) !/DEBUGDCXDX sumDiff1 = sumDiff1 + MAX(FKD0,ZERO) !/DEBUGDCXDX sumDiff2 = sumDiff2 + MAX(WN(IK,ISEA),ZERO) !/DEBUGDCXDX sumDiff3 = sumDiff3 + MAX(FKC(ITH),ZERO) END DO END DO !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff=', sumDiff !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff1=', sumDiff1 !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff2=', sumDiff2 !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff3=', sumDiff3 sumDiff=0 DO ISP=1-NTH,NSPEC CWNB_M2(ISP) = DBLE(0.5 * ( VCWN(ISP) + VCWN(ISP+NTH) )) sumDiff = sumDiff + MAX(CWNB_M2(ISP), ZERO) END DO !/DEBUGDCXDX WRITE(740+IAPROC,*) 'sumDiff=', sumDiff DO IK=1,NK DWNI_M2(IK) = DBLE( CG(IK,ISEA) / DSIP(IK) ) END DO !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / SUBROUTINE SYNCHRONIZE_IPGL_ETC_ARRAY(IMOD, IsMulti) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Sync global local arrays ! 2. Method : ! All the process need to have IPGL_tot and IPGL_TO_PROC ! This is especially the case for the output process. ! So we need some painful exportation business ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE !/PDLIB USE yowDatapool, only: istatus !/PDLIB USE yowNodepool, only: np_global !/PDLIB USE W3ODATMD, ONLY: NTPROC, NAPROC, IAPROC !/PDLIB USE W3GDATMD, ONLY: MAPSF, NSEA !/PDLIB USE W3ADATMD, ONLY: MPI_COMM_WAVE, MPI_COMM_WCMP !/PDLIB USE yowRankModule, only : IPGL_TO_PROC, IPGL_tot !/PDLIB USE WMMDATMD, ONLY: MDATAS IMPLICIT NONE !/PDLIB INCLUDE "mpif.h" INTEGER, intent(in) :: IMOD logical, intent(in) :: IsMulti !/S INTEGER, SAVE :: IENT = 0 !/PDLIB INTEGER :: Iarr(1) !/PDLIB INTEGER :: ISEA, IP_glob !/PDLIB INTEGER :: IPROC, IERR_MPI, istat !/S CALL STRACE (IENT, 'SYNCHRONIZE_IPGL_ETC_ARRAY') !/PDLIB IF (IAPROC .le. NAPROC) THEN !/PDLIB IF (IAPROC .eq. 1) THEN !/PDLIB Iarr(1)=np_global !/PDLIB DO IPROC=NAPROC+1,NTPROC !/PDLIB CALL MPI_SEND(Iarr,1,MPI_INT, IPROC-1, 37, MPI_COMM_WAVE, IERR_MPI) !/PDLIB END DO !/PDLIB DO IPROC=NAPROC+1,NTPROC !/PDLIB CALL MPI_SEND(ipgl_tot,np_global,MPI_INT, IPROC-1, 43, MPI_COMM_WAVE, IERR_MPI) !/PDLIB CALL MPI_SEND(ipgl_to_proc,np_global,MPI_INT, IPROC-1, 91, MPI_COMM_WAVE, IERR_MPI) !/PDLIB END DO !/PDLIB END IF !/PDLIB ELSE !/PDLIB CALL MPI_RECV(Iarr,1,MPI_INT, 0, 37, MPI_COMM_WAVE, istatus, IERR_MPI) !/PDLIB np_global=Iarr(1) !/PDLIB allocate(IPGL_tot(np_global), IPGL_TO_PROC(np_global), stat=istat) !/PDLIB CALL MPI_RECV(ipgl_tot,np_global,MPI_INT, 0, 43, MPI_COMM_WAVE, istatus, IERR_MPI) !/PDLIB CALL MPI_RECV(ipgl_to_proc,np_global,MPI_INT, 0, 91, MPI_COMM_WAVE, istatus, IERR_MPI) !/PDLIB END IF !/PDLIB IF (IsMulti) THEN !/PDLIB WRITE(*,*) ' Before allocation of MDATAS % SEA_IPGL, SEA_IPGL_TO_PROC : IMOD=', IMOD, ' NSEA=', NSEA !/PDLIB ALLOCATE(MDATAS(IMOD)%SEA_IPGL(NSEA), MDATAS(IMOD)%SEA_IPGL_TO_PROC(NSEA), STAT=ISTAT) !/PDLIB !CHECK_ALLOC_STATUS ( ISTAT ) !/PDLIB DO ISEA=1,NSEA !/PDLIB IP_glob = MAPSF(ISEA, 1) !/PDLIB MDATAS(IMOD)%SEA_IPGL(ISEA) = IPGL_tot(IP_glob) !/PDLIB MDATAS(IMOD)%SEA_IPGL_TO_PROC(ISEA) = IPGL_TO_PROC(IP_glob) !/PDLIB END DO !/PDLIB END IF !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ....................----------------------------------------------- / SUBROUTINE SET_UP_NSEAL_NSEALM(NSEALout, NSEALMout) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Setup nseal, nsealm in contect of pdlib ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE !/ !/ ------------------------------------------------------------------- / !/PDLIB use yowDatapool, only: istatus !/PDLIB use yowNodepool, only: npa !/PDLIB use yowRankModule, only : rank !/PDLIB USE W3GDATMD, ONLY: GTYPE, UNGTYPE !/MPI USE W3ADATMD, ONLY: MPI_COMM_WAVE, MPI_COMM_WCMP USE CONSTANTS, ONLY : LPDLIB USE W3GDATMD, ONLY: NSEA USE W3ODATMD, ONLY: NTPROC, NAPROC, IAPROC IMPLICIT NONE INTEGER, intent(out) :: NSEALout, NSEALMout !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 !/ !/ ------------------------------------------------------------------- / !/ !/S CALL STRACE (IENT, 'SET_UP_NSEAL_NSEALM') !!/PDLIB WRITE(*,*) 'LPDLIB=', LPDLIB !!/PDLIB WRITE(*,*) 'GTYPE=', GTYPE, ' UNGTYPE=', UNGTYPE !/DEBUG WRITE(740+IAPROC,*) 'SET_UP, PDLIB=', LPDLIB !/DEBUG FLUSH(740+IAPROC) !/SHRD NSEALout = NSEA !/SHRD NSEALMout = NSEA ! !/DIST IF (.NOT. LPDLIB ) THEN !/DIST IF ( IAPROC .LE. NAPROC ) THEN !/DIST NSEALout = 1 + (NSEA-IAPROC)/NAPROC !/DIST ELSE !/DIST NSEALout = 0 !/DIST END IF !/DIST NSEALMout = 1 + (NSEA-1)/NAPROC !/DIST ELSE !/PDLIB IF (GTYPE .eq. UNGTYPE) THEN !/PDLIB NSEALout = PDLIB_NSEAL !/PDLIB NSEALMout = PDLIB_NSEALM !/PDLIB WRITE(*,*) 'PDLIB_NSEAL=', PDLIB_NSEAL, ' PDLIB_NSEALM=', PDLIB_NSEALM !/PDLIB ELSE !/PDLIB IF ( IAPROC .LE. NAPROC ) THEN !/PDLIB NSEALout = 1 + (NSEA-IAPROC)/NAPROC !/PDLIB ELSE !/PDLIB NSEALout = 0 !/PDLIB END IF !/PDLIB NSEALMout = 1 + (NSEA-1)/NAPROC !/PDLIB ENDIF !/DIST ENDIF !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / SUBROUTINE INIT_GET_JSEA_ISPROC(ISEA, JSEA, ISPROC) !/ ------------------------------------------------------------------- / !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Set Jsea for all schemes ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE !/ USE W3ODATMD, ONLY: OUTPTS, IAPROC, NAPROC USE W3GDATMD, ONLY: GTYPE, UNGTYPE, MAPSF USE CONSTANTS, ONLY : LPDLIB !/PDLIB USE yowRankModule, only : IPGL_TO_PROC, IPGL_tot !/PDLIB use yowNodepool, only: ipgl, iplg IMPLICIT NONE !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 !/ !/ ------------------------------------------------------------------- / INTEGER, intent(in) :: ISEA INTEGER, intent(out) :: JSEA, ISPROC INTEGER IP_glob !/S CALL STRACE (IENT, 'INIT_GET_JSEA_ISPROC') !!/DEBUG WRITE(740+IAPROC,*) 'PDLIB=', PDLIB !!/DEBUG WRITE(740+IAPROC,*) 'GTYPE=', GTYPE, ' UNGTYPE=', UNGTYPE !!/DEBUG FLUSH(740+IAPROC) IF (.NOT. LPDLIB) THEN JSEA = 1 + (ISEA-1)/NAPROC ISPROC = ISEA - (JSEA-1)*NAPROC ELSE !/PDLIB IF (GTYPE .ne. UNGTYPE) THEN !/PDLIB JSEA = 1 + (ISEA-1)/NAPROC !/PDLIB ISPROC = ISEA - (JSEA-1)*NAPROC !/PDLIB ELSE !/PDLIB IP_glob = MAPSF(ISEA,1) !/PDLIB IF (IAPROC .le. NAPROC) THEN !/PDLIB JSEA = ISEA_TO_JSEA(ISEA) !/PDLIB ELSE !/PDLIB JSEA = -1 !/PDLIB END IF !/PDLIB ISPROC = IPGL_TO_PROC(IP_glob) !/PDLIB ENDIF ENDIF !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / SUBROUTINE GET_JSEA_IBELONG(ISEA, JSEA, IBELONG) !/ ------------------------------------------------------------------- / !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Set belongings of jsea in context of pdlib ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE !/ USE W3ODATMD, ONLY: OUTPTS, IAPROC, NAPROC USE W3GDATMD, ONLY: GTYPE, UNGTYPE, MAPSF USE CONSTANTS, ONLY : LPDLIB !/PDLIB USE yowRankModule, only : IPGL_TO_PROC, IPGL_tot, IPGL_npa !/PDLIB use yowNodepool, only: ipgl, iplg IMPLICIT NONE !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 !/ !/ ------------------------------------------------------------------- / !/ INTEGER, intent(in) :: ISEA INTEGER, intent(out) :: JSEA, IBELONG INTEGER ISPROC, IX, JX !/S CALL STRACE (IENT, 'GET_JSEA_IBELONG') IF (.NOT. LPDLIB) THEN JSEA = 1 + (ISEA-1)/NAPROC ISPROC = ISEA - (JSEA-1)*NAPROC IF (ISPROC .eq. IAPROC) THEN IBELONG=1 ELSE IBELONG=0 END IF ELSE !/PDLIB IF (GTYPE .ne. UNGTYPE) THEN !/PDLIB JSEA = 1 + (ISEA-1)/NAPROC !/PDLIB ISPROC = ISEA - (JSEA-1)*NAPROC !/PDLIB IF (ISPROC .eq. IAPROC) THEN !/PDLIB IBELONG=1 !/PDLIB ELSE !/PDLIB IBELONG=0 !/PDLIB END IF !/PDLIB ELSE !/PDLIB IF (IAPROC .le. NAPROC) THEN !/PDLIB IX = MAPSF(ISEA,1) !/PDLIB JX = IPGL_npa(IX) !/PDLIB IF (JX .eq. 0) THEN !/PDLIB IBELONG=0 !/PDLIB JSEA=-1 !/PDLIB ELSE !/PDLIB IBELONG=1 !/PDLIB JSEA=JX_TO_JSEA(JX) !/PDLIB END IF !/PDLIB ELSE !/PDLIB IBELONG=0 !/PDLIB JSEA=-1 !/PDLIB END IF !/PDLIB ENDIF ENDIF !/ !/ End of INIT_GET_ISEA ---------------------------------------------- / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / SUBROUTINE INIT_GET_ISEA(ISEA, JSEA) !/ ------------------------------------------------------------------- / !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Set Isea for all schemes ! 2. Method : ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE !/ USE W3ODATMD, ONLY: OUTPTS, IAPROC, NAPROC USE W3GDATMD, ONLY: GTYPE, UNGTYPE USE CONSTANTS, ONLY : LPDLIB !/PDLIB USE YOWNODEPOOL, ONLY: iplg !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/ ------------------------------------------------------------------- / !/ !/ !/ ------------------------------------------------------------------- / !/ !/ ------------------------------------------------------------------- / ! USE W3ODATMD, ONLY: OUTPTS, IAPROC, NAPROC USE W3GDATMD, ONLY: GTYPE, UNGTYPE USE CONSTANTS, ONLY : LPDLIB !/PDLIB USE YOWNODEPOOL, ONLY: iplg IMPLICIT NONE !/S INTEGER, SAVE :: IENT = 0 INTEGER, intent(in) :: JSEA INTEGER, intent(out) :: ISEA !/S CALL STRACE (IENT, 'INIT_GET_ISEA') !/SHRD ISEA = JSEA !/DIST IF (.NOT. LPDLIB) THEN !/DIST ISEA = IAPROC + (JSEA-1)*NAPROC !/DIST ELSE !/PDLIB IF (GTYPE .eq. UNGTYPE) THEN !/PDLIB ISEA = iplg(JSEA) !/PDLIB ELSE !/PDLIB ISEA = IAPROC + (JSEA-1)*NAPROC !/PDLIB ENDIF !/DIST ENDIF !/ !/ End of INIT_GET_ISEA ------------------------------------------------ / !/ END SUBROUTINE !********************************************************************** !* An array of size (NSEA) is send but only the (1:NSEAL) values * !* are correct. The program synchonizes everything on all nodes. * !********************************************************************** SUBROUTINE SYNCHRONIZE_GLOBAL_ARRAY(TheVar) !/ ------------------------------------------------------------------- / !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | | !/ | Aron Roland (BGS IT&E GmbH) | !/ | Mathieu Dutour-Sikiric (IRB) | !/ | | !/ | FORTRAN 90 | !/ | Last update : 01-June-2018 | !/ +-----------------------------------+ !/ !/ 01-June-2018 : Origination. ( version 6.04 ) !/ ! 1. Purpose : Sync global array in context of pdlib ! 2. Method : ! An array of size (NSEA) is send but only the (1:NSEAL) values ! are correct. The program synchonizes everything on all nodes. ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! 7. Remarks ! 8. Structure : ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/S USE W3SERVMD, ONLY: STRACE ! USE W3GDATMD, ONLY: NSEAL, NSEA, NX !/PDLIB USE W3ODATMD, only : IAPROC, NAPROC, NTPROC !/PDLIB USE W3ADATMD, ONLY: MPI_COMM_WCMP !/PDLIB use yowDatapool, only: rtype, istatus !/PDLIB USE yowNodepool, only: npa !/PDLIB use yowNodepool, only: iplg IMPLICIT NONE !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 !/ !/ ------------------------------------------------------------------- / !/ !/MPI INCLUDE "mpif.h" INTEGER ISEA, JSEA, Status(NX), rStatus(NX) INTEGER IPROC, I, ierr, IP, IX, IP_glob REAL*8, intent(inout) :: TheVar(NX) REAL*8 rVect(NX) Status=0 !/S CALL STRACE (IENT, 'SYNCHRONIZE_GLOBAL_ARRAY') !/PDLIB DO IP=1,npa !/PDLIB IP_glob=iplg(IP) !/PDLIB Status(IP_glob)=1 !/PDLIB END DO !/PDLIB IF (IAPROC .eq. 1) THEN !/PDLIB DO iProc=2,NAPROC !/PDLIB CALL MPI_RECV(rVect,NX,rtype, iProc-1, 19, MPI_COMM_WCMP, istatus, ierr) !/PDLIB CALL MPI_RECV(rStatus,NX,MPI_INTEGER, iProc-1, 23, MPI_COMM_WCMP, istatus, ierr) !/PDLIB DO I=1,NX !/PDLIB IF (rStatus(I) .eq. 1) THEN !/PDLIB TheVar(I)=rVect(I) !/PDLIB Status(I)=1 !/PDLIB END IF !/PDLIB END DO !/PDLIB END DO !/PDLIB DO IPROC=2,NAPROC !/PDLIB CALL MPI_SEND(TheVar,NX,rtype, iProc-1, 29, MPI_COMM_WCMP, ierr) !/PDLIB END DO !/PDLIB ELSE !/PDLIB CALL MPI_SEND(TheVar,NX,rtype, 0, 19, MPI_COMM_WCMP, ierr) !/PDLIB CALL MPI_SEND(Status,NX,MPI_INTEGER, 0, 23, MPI_COMM_WCMP, ierr) !/PDLIB CALL MPI_RECV(TheVar,NX,rtype, 0, 29, MPI_COMM_WCMP, istatus, ierr) !/PDLIB END IF !/ !/ End of JACOBI_INIT ------------------------------------------------ / !/ END SUBROUTINE !/ ------------------------------------------------------------------- / END MODULE W3PARALL !/ ------------------------------------------------------------------- /