#include "w3macros.h" !/ ------------------------------------------------------------------- / MODULE W3WAVEMD !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | FORTRAN 90 | !/ | Last update : 27-Aug-2015 | !/ +-----------------------------------+ !/ !/ 04-Feb-2000 : Origination. ( version 2.00 ) !/ For upgrades see subroutines. !/ 14-Feb-2000 : Exact-NL added. ( version 2.01 ) !/ 05-Jan-2001 : Bug fix to allow model to run ( version 2.05 ) !/ without output. !/ 24-Jan-2001 : Flat grid version. ( version 2.06 ) !/ 09-Feb-2001 : Third propagation scheme added. ( version 2.08 ) !/ 23-Feb-2001 : Check for barrier after source !/ terms added ( W3NMIN ). ( delayed version 2.07 ) !/ 16-Mar-2001 : Fourth propagation scheme added. ( version 2.09 ) !/ 30-Mar-2001 : Sub-grid obstacles added. ( version 2.10 ) !/ 23-May-2001 : Clean up and bug fixes. ( version 2.11 ) !/ 10-Dec-2001 : Sub-grid obstacles for UQ schemes. ( version 2.14 ) !/ 11-Jan-2002 : Sub-grid ice. ( version 2.15 ) !/ 24-Jan-2002 : Zero time step dor data ass. ( version 2.17 ) !/ 18-Feb-2002 : Point output diagnostics added. ( version 2.18 ) !/ 30-Apr-2002 : Add field output types 17-18. ( version 2.20 ) !/ 09-May-2002 : Switch clean up. ( version 2.21 ) !/ 13-Nov-2002 : Add stress vector. ( version 3.00 ) !/ 26-Dec-2002 : Moving grid version. ( version 3.02 ) !/ 01-Aug-2003 : Moving grid GSE correction. ( version 3.03 ) !/ 20-Aug-2003 : Output server options added. ( version 3.04 ) !/ 07-Oct-2003 : Output options for NN training. ( version 3.05 ) !/ 29-Dec-2004 : Multiple grid version. ( version 3.06 ) !/ W3INIT, W3MPII-O and WWVER moved to w3initmd.ftn !/ 04-Feb-2005 : Add STAMP to par list of W3WAVE. ( version 3.07 ) !/ 04-May-2005 : Change to MPI_COMM_WAVE. ( version 3.07 ) !/ 28-Jun-2005 : Adding map recalc for W3ULEV call. ( version 3.07 ) !/ 07-Sep-2005 : Updated boundary conditions. ( version 3.08 ) !/ Fix NRQSG1/2 = 0 array bound issue. !/ 13-Jun-2006 : Split STORE in G/SSTORE ( version 3.09 ) !/ 26-Jun-2006 : Add output type 6. ( version 3.09 ) !/ 04-Jul-2006 : Consolidate stress arrays. ( version 3.09 ) !/ 18-Oct-2006 : Partitioned spectral data output. ( version 3.10 ) !/ 02-Feb-2007 : Add FLAGST test. ( version 3.10 ) !/ 02-Apr-2007 : Add partitioned field data. ( version 3.11 ) !/ 07-May-2007 : Bug fix SKIP_O treatment. ( version 3.11 ) !/ 17-May-2007 : Adding NTPROC/NAPROC separation. ( version 3.11 ) !/ 08-Oct-2007 : Adding AS CX-Y to W3SRCE par. list. ( version 3.13 ) !/ 22-Feb-2008 : Initialize VGX-Y properly. ( version 3.13 ) !/ 10-Apr-2008 : Bug fix writing log file (MPI). ( version 3.13 ) !/ 29-May-2009 : Preparing distribution version. ( version 3.14 ) !/ 30-Oct-2009 : Implement run-time grid selection. ( version 3.14 ) !/ (W. E. Rogers & T. J. Campbell, NRL) !/ 30-Oct-2009 : Implement curvilinear grid type. ( version 3.14 ) !/ (W. E. Rogers & T. J. Campbell, NRL) !/ 29-Mar-2010 : Adding coupling, ice in W3SRCE. ( version 3.14_SHOM ) !/ 16-May-2010 : Adding transparencies in W3SCRE ( version 3.14_SHOM ) !/ 23-Jun-2011 : Movable bed bottom friction BT4 ( version 4.04 ) !/ 03-Nov-2011 : Shoreline reflection on unst. grids ( version 4.04 ) !/ 02-Jul-2011 : Update for PALM coupling ( version 4.07 ) !/ 06-Mar-2012 : Initializing ITEST as needed. ( version 4.07 ) !/ 02-Jul-2012 : Update for PALM coupling ( version 4.07 ) !/ 02-Sep-2012 : Clean up of open BC for UG grids ( version 4.08 ) !/ 03-Sep-2012 : Fix format 902. ( version 4.10 ) !/ 07-Dec-2012 : Wrap W3SRCE with TMPn to limit WARN ( version 4.OF ) !/ 10-Dec-2012 : Modify field output MPI for new ( version 4.OF ) !/ structure and smaller memory footprint. !/ 12-Dec-2012 : Adding SMC grid. JG_Li ( version 4.08 ) !/ 26-Dec-2012 : Move FIELD init. to W3GATH. ( version 4.OF ) !/ 16-Sep-2013 : Add Arctic part for SMC grid. ( version 4.11 ) !/ 11-Nov-2013 : SMC and rotated grid incorporated in the main !/ trunk ( version 4.13 ) !/ 14-Nov-2013 : Remove orphaned work arrays. ( version 4.13 ) !/ 27-Nov-2013 : Fixes for OpenMP versions. ( version 4.15 ) !/ 23-May-2014 : Adding ice fluxes to W3SRCE ( version 5.01 ) !/ 27-May-2014 : Move to OMPG/X switch. ( version 5.02 ) !/ 24-Apr-2015 : Adding OASIS coupling calls ( version 5.07 ) !/ (M. Accensi & F. Ardhuin, IFREMER) !/ 27-Aug-2015 : Update for ICEH, ICEF ( version 5.08 ) !/ 14-Sep-2018 : Remove PALM implementation ( version 6.06 ) !/ !/ Copyright 2009-2014 National Weather Service (NWS), !/ National Oceanic and Atmospheric Administration. All rights !/ reserved. WAVEWATCH III is a trademark of the NWS. !/ No unauthorized use without permission. !/ ! 1. Purpose : ! ! 2. Variables and types : ! ! 3. Subroutines and functions : ! ! Name Type Scope Description ! ---------------------------------------------------------------- ! W3WAVE Subr. Public Actual wave model. ! W3GATH Subr. Public Data transpose before propagation. ! W3SCAT Subr. Public Data transpose after propagation. ! W3NMIN Subr. Public Calculate minimum number of sea ! points per processor. ! ---------------------------------------------------------------- ! ! 4. Subroutines and functions used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3SETx Subr. W3xDATMD Point to data structure. ! ! W3UCUR Subr. W3UPDTMD Interpolate current fields in time. ! W3UWND Subr. W3UPDTMD Interpolate wind fields in time. ! W3UINI Subr. W3UPDTMD Update initial conditions if init. ! with initial wind conditions. ! W3UBPT Subr. W3UPDTMD Update boundary points. ! W3UICE Subr. W3UPDTMD Update ice coverage. ! W3ULEV Subr. W3UPDTMD Transform the wavenumber grid. ! W3DDXY Subr. W3UPDTMD Calculate dirivatives of the depth. ! W3DCXY Subr. W3UPDTMD Calculate dirivatives of the current. ! ! W3MAPn Subr. W3PROnMD Preparation for ropagation schemes. ! W3XYPn Subr. W3PROnMD Longitude-latitude ("XY") propagation. ! W3KTPn Subr. W3PROnMD Intra-spectral ("k-theta") propagation. ! ! W3SRCE Subr. W3SRCEMD Source term integration and calculation. ! ! W3IOGR Subr. W3IOGRMD Reading/writing model definition file. ! W3OUTG Subr. W3IOGOMD Generate gridded output fields. ! W3IOGO Subr. W3IOGOMD Read/write gridded output. ! W3IOPE Subr. W3IOPOMD Extract point output. ! W3IOPO Subr. W3IOPOMD Read/write point output. ! W3IOTR Subr. W3IOTRMD Process spectral output along tracks. ! W3IORS Subr. W3IORSMD Read/write restart files. ! W3IOBC Subr. W3IOBCMD Read/write boundary conditions. ! W3CPRT Subr. W3IOSFMD Partition spectra. ! W3IOSF Subr. Id. Write partitioned spectral data. ! ! STRACE Subr. W3SERVMD Subroutine tracing. ! WWTIME Subr. Id. System time in readable format. ! EXTCDE Subr. Id. Program abort. ! ! TICK21 Subr. W3TIMEMD Advance the clock. ! DSEC21 Func. Id. Difference between times. ! STME21 Subr. Id. Time in readable format. ! ! MPI_BARRIER, MPI_STARTALL, MPI_WAITALL ! Subr. Basic MPI routines. ! ---------------------------------------------------------------- ! ! 5. Remarks : ! ! 6. Switches : ! ! !/SHRD Switch for shared / distributed memory architecture. ! !/DIST Id. ! !/MPI Id. ! !/OMPG Id. ! !/OMPX Id. ! ! !/PR1 First order propagation schemes. ! !/PR2 ULTIMATE QUICKEST scheme. ! !/PR3 Averaged ULTIMATE QUICKEST scheme. ! !/PRX User-defined scheme. ! !/SMC UNO2 scheme on SMC grid. ! ! !/S Enable subroutine tracing. ! !/T Test output. ! !/MPIT Test output for MPI specific code. ! ! 7. Source code : ! !/ ------------------------------------------------------------------- / USE W3ADATMD, ONLY: MPIBUF ! PUBLIC !/ CONTAINS !/ ------------------------------------------------------------------- / SUBROUTINE W3WAVE ( IMOD, TEND, STAMP, NO_OUT & ,ID_LCOMM & ) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | FORTRAN 90 | !/ | Last update : 06-Jun-2018 | !/ +-----------------------------------+ !/ !/ 17-Mar-1999 : Distributed FORTRAN 77 version. ( version 1.18 ) !/ 04-Feb-2000 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ Major changes to logistics. !/ 05-Jan-2001 : Bug fix to allow model to run ( version 2.05 ) !/ without output. !/ 24-Jan-2001 : Flat grid version. ( version 2.06 ) !/ 09-Feb-2001 : Third propagation scheme added. ( version 2.08 ) !/ 23-Feb-2001 : Check for barrier after source !/ terms added ( W3NMIN ). ( delayed version 2.07 ) !/ 16-Mar-2001 : Fourth propagation scheme added. ( version 2.09 ) !/ 30-Mar-2001 : Sub-grid obstacles added. ( version 2.10 ) !/ 23-May-2001 : Barrier added for dry run, changed ( version 2.10 ) !/ declaration of FLIWND. !/ 10-Dec-2001 : Sub-grid obstacles for UQ schemes. ( version 2.14 ) !/ 11-Jan-2002 : Sub-grid ice. ( version 2.15 ) !/ 24-Jan-2002 : Zero time step dor data ass. ( version 2.17 ) !/ 09-May-2002 : Switch clean up. ( version 2.21 ) !/ 13-Nov-2002 : Add stress vector. ( version 3.00 ) !/ 26-Dec-2002 : Moving grid version. ( version 3.02 ) !/ 01-Aug-2003 : Moving grid GSE correction. ( version 3.03 ) !/ 07-Oct-2003 : Output options for NN training. ( version 3.05 ) !/ 29-Dec-2004 : Multiple grid version. ( version 3.06 ) !/ 04-Feb-2005 : Add STAMP to par list. ( version 3.07 ) !/ 04-May-2005 : Change to MPI_COMM_WAVE. ( version 3.07 ) !/ 28-Jun-2005 : Adding map recalc for W3ULEV call. ( version 3.07 ) !/ 07-Sep-2005 : Updated boundary conditions. ( version 3.08 ) !/ 26-Jun-2006 : Add output type 6. ( version 3.09 ) !/ 04-Jul-2006 : Consolidate stress arrays. ( version 3.09 ) !/ 18-Oct-2006 : Partitioned spectral data output. ( version 3.10 ) !/ 02-Feb-2007 : Add FLAGST test. ( version 3.10 ) !/ 02-Apr-2007 : Add partitioned field data. ( version 3.11 ) !/ Improve MPI_WAITALL call tests/allocations. !/ 07-May-2007 : Bug fix SKIP_O treatment. ( version 3.11 ) !/ 17-May-2007 : Adding NTPROC/NAPROC separation. ( version 3.11 ) !/ 08-Oct-2007 : Adding AS CX-Y to W3SRCE par. list. ( version 3.13 ) !/ 22-Feb-2008 : Initialize VGX-Y properly. ( version 3.13 ) !/ 10-Apr-2008 : Bug fix writing log file (MPI). ( version 3.13 ) !/ 30-Oct-2009 : Implement run-time grid selection. ( version 3.14 ) !/ (W. E. Rogers & T. J. Campbell, NRL) !/ 30-Oct-2009 : Implement curvilinear grid type. ( version 3.14 ) !/ (W. E. Rogers & T. J. Campbell, NRL) !/ 31-Mar-2010 : Add reflections ( version 3.14.4 ) !/ 29-Oct-2010 : Implement unstructured grids ( version 3.14.4 ) !/ (A. Roland and F. Ardhuin) !/ 06-Mar-2011 : Output of max. CFL (F.Ardhuin) ( version 3.14.4 ) !/ 05-Apr-2011 : Implement iteration for DTMAX <1s ( version 3.14.4 ) !/ 02-Jul-2012 : Update for PALM coupling ( version 4.07 ) !/ 02-Sep-2012 : Clean up of open BC for UG grids ( version 4.08 ) !/ 03-Sep-2012 : Fix format 902. ( version 4.10 ) !/ 10-Dec-2012 : Modify field output MPI for new ( version 4.OF ) !/ structure and smaller memory footprint. !/ 16-Nov-2013 : Allows reflection on curvi. grids ( version 4.13 ) !/ 27-Nov-2013 : Fixes for OpenMP versions. ( version 4.15 ) !/ 23-May-2014 : Adding ice fluxes to W3SRCE ( version 5.01 ) !/ 27-May-2014 : Move to OMPG/X switch. ( version 5.02 ) !/ 24-Apr-2015 : Adding OASIS coupling calls ( version 5.07 ) !/ (M. Accensi & F. Ardhuin, IFREMER) !/ 27-Aug-2015 : Update for ICEH, ICEF ( version 5.10 ) !/ 31-Mar-2016 : Current option for smc grid. ( version 5.18 ) !/ 06-Jun-2018 : Add PDLIB/MEMCHECK/SETUP/NETCDF_QAD/TIMING !/ OASIS/DEBUGINIT/DEBUGSRC/DEBUGRUN/DEBUGCOH !/ DEBUGIOBP/DEBUGIOBC ( version 6.04 ) !/ 14-Sep-2018 : Remove PALM implementation ( version 6.06 ) !/ ! 1. Purpose : ! ! Run WAVEWATCH III for a given time interval. ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! IMOD Int. I Model number. ! TEND I.A. I Ending time of integration. ! STAMP Log. I WRITE(*,*)time stamp (optional, defaults to T). ! NO_OUT Log. I Skip output (optional, defaults to F). ! Skip at ending time only! ! ---------------------------------------------------------------- ! ! Local parameters : Flags ! ---------------------------------------------------------------- ! FLOUTG Log. Flag for running W3OUTG. ! FLPART Log. Flag for running W3CPRT. ! FLZERO Log. Flag for zero time interval. ! FLAG0 Log. Flag for processors without tasks. ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! See module documentation. ! ! 5. Called by : ! ! Any program shell or integrated model which uses WAVEWATCH III. ! ! 6. Error messages : ! ! 7. Remarks : ! ! - Currents are updated before winds as currents are used in wind ! and USTAR processing. ! - Ice and water levels can be updated only once per call. ! - If ice or water level time are undefined, the update ! takes place asap, otherwise around the "half-way point" ! betweem the old and new times. ! - To increase accuracy, the calculation of the intra-spectral ! propagation is performed in two parts around the spatial propagation. ! ! 8. Structure : ! ! ----------------------------------------------------------- ! 0. Initializations ! a Point to data structures ! b Subroutine tracing ! c Local parameter initialization ! d Test output ! 1. Check the consistency of the input. ! a Ending time versus initial time. ! b Water level time. ! c Current time interval. ! d Wind time interval. ! e Ice time. ! 2. Determine next time from ending and output ! time and get corresponding time step. ! 3. Loop over time steps (see below). ! 4. Perform output to file if requested. ! a Check if time is output time. ! b Processing and MPP preparations. ( W3CPRT, W3OUTG ) ! c Reset next output time. ! -------------- loop over output types ------------------ ! d Perform output. ( W3IOxx ) ! e Update next output time. ! -------------------- end loop -------------------------- ! 5. Update log file. ! 6. If time is not ending time, branch back to 2. ! ----------------------------------------------------------- ! ! Section 3. ! ---------------------------------------------------------- ! 3.1 Interpolate winds and currents. ( W3UCUR, W3DCXY ) ! ( W3UWND ) ! ( W3UINI ) ! 3.2 Update boundary conditions. ( W3IOBC, W3UBPT ) ! 3.3 Update ice coverage (if new ice map). ( W3UICE ) ! 3.4 Transform grid (if new water level). ( W3ULEV ) ! 3.5 Update maps and dirivatives. ( W3MAPn, W3DDXY ) ! ( W3NMIN, W3UTRN ) ! Update grid advection vector. ! 3.6 Perform propagation ! a Preparations. ! b Intra spectral part 1. ( W3KTPn ) ! c Longitude-latitude ( W3GATH, W3XYPn W3SCAT ) ! b Intra spectral part 2. ( W3KTPn ) ! 3.7 Calculate and integrate source terms. ( W3SRCE ) ! 3.8 Update global time step. ! ---------------------------------------------------------- ! ! 9. Switches : ! ! See module documentation. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE CONSTANTS !/ USE W3GDATMD USE W3WDATMD USE W3ADATMD USE W3IDATMD USE W3ODATMD !/ USE W3UPDTMD USE W3SRCEMD USE W3PRO3MD ! USE W3PROFSMD !/ USE W3TRIAMD USE W3IOGRMD USE W3IOGOMD USE W3IOPOMD USE W3IOTRMD USE W3IORSMD USE W3IOBCMD USE W3IOSFMD !/ USE W3SERVMD USE W3TIMEMD USE W3PARALL, ONLY : INIT_GET_ISEA USE W3OACPMD, ONLY: ID_OASIS_TIME USE W3OGCMMD, ONLY: SND_FIELDS_TO_OCEAN USE W3AGCMMD, ONLY: SND_FIELDS_TO_ATMOS ! IMPLICIT NONE ! INCLUDE "mpif.h" !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ INTEGER, INTENT(IN) :: IMOD, TEND(2) LOGICAL, INTENT(IN), OPTIONAL :: STAMP, NO_OUT INTEGER, INTENT(IN), OPTIONAL :: ID_LCOMM !/ !/ ------------------------------------------------------------------- / !/ Local parameters : !/ INTEGER :: IP INTEGER :: TCALC(2), IT, IT0, NT, ITEST, & ITLOC, ITLOCH, NTLOC, ISEA, JSEA, & IX, IY, ISPEC, J, TOUT(2), TLST(2), & REFLED(6), IK, ITH, IS, NKCFL INTEGER :: ISP, IP_glob REAL :: ICEDAVE INTEGER :: OASISED INTEGER :: IERR_MPI, NRQMAX INTEGER, ALLOCATABLE :: STATCO(:,:), STATIO(:,:) INTEGER :: IXrel REAL :: DTTST, DTTST1, DTTST2, DTTST3, & DTL0, DTI0, DTI10, DTI50, & DTGA, DTG, DTGpre, DTRES, & FAC, VGX, VGY, FACK, FACTH, & FACX, XXX, REFLEC(4), & DELX, DELY, DELA, DEPTH, D50, PSIC REAL :: VSioDummy(NSPEC), VDioDummy(NSPEC), VAoldDummy(NSPEC) LOGICAL :: SHAVETOTioDummy ! REAL, ALLOCATABLE :: FIELD(:) REAL :: TMP1(4), TMP2(3), TMP3(2), TMP4(2) ! ! Orphaned arrays from old data structure ! REAL, ALLOCATABLE :: TAUWX(:), TAUWY(:) ! LOGICAL :: FLACT, FLZERO, FLFRST, FLMAP, TSTAMP,& SKIP_O, FLAG_O, FLDDIR, READBC, & FLAG0 = .FALSE., FLOUTG, FLPFLD, & FLPART, LOCAL, FLOUTG2, & FLOMP = .FALSE. ! !Li Logical variable to control regular gird lines in conflict with SMC option. !AR SMC option is in conflict with lofical variables for regular grid ... chicken ... egg ... stuff LOGICAL :: RGLGRD = .TRUE., ARCTIC = .FALSE. !!Li LOGICAL :: FLGMPI(0:6) LOGICAL :: UGDTUPDATE ! true if time step should be updated for UG schemes CHARACTER(LEN=8) :: STTIME CHARACTER(LEN=17) :: IDACT CHARACTER(LEN=13) :: OUTID CHARACTER(LEN=23) :: IDTIME INTEGER eIOBP INTEGER ITH_F ! !/ !/ ------------------------------------------------------------------- / ! 0. Initializations ! ! 0.a Set pointers to data structure ! SCREEN = 333 ! IF ( IOUTP .NE. IMOD ) CALL W3SETO ( IMOD, NDSE, NDST ) IF ( IGRID .NE. IMOD ) CALL W3SETG ( IMOD, NDSE, NDST ) IF ( IWDATA .NE. IMOD ) CALL W3SETW ( IMOD, NDSE, NDST ) IF ( IADATA .NE. IMOD ) CALL W3SETA ( IMOD, NDSE, NDST ) IF ( IIDATA .NE. IMOD ) CALL W3SETI ( IMOD, NDSE, NDST ) ! ALLOCATE(TAUWX(NSEAL), TAUWY(NSEAL)) ! IF ( PRESENT(STAMP) ) THEN TSTAMP = STAMP ELSE TSTAMP = .TRUE. END IF ! IF ( PRESENT(NO_OUT) ) THEN SKIP_O = NO_OUT ELSE SKIP_O = .FALSE. END IF ! ! 0.b Subroutine tracing ! ! 0.c Local parameter initialization ! IPASS = IPASS + 1 IDACT = ' ' OUTID = ' ' FLACT = ITIME .EQ. 0 FLMAP = ITIME .EQ. 0 FLDDIR = ITIME .EQ. 0 .AND. ( FLCTH .OR. FSREFRACTION & .OR. FLCK .OR. FSFREQSHIFT ) ! FLPFLD = .FALSE. DO J=1,NOGE(4) FLPFLD = FLPFLD .OR. FLOGRD(4,J) .OR. FLOGR2(4,J) END DO ! IF ( IAPROC .EQ. NAPLOG ) BACKSPACE ( NDSO ) ! IF ( FLCOLD ) THEN DTDYN = 0. FCUT = SIG(NK) * TPIINV END IF ! !!Li ALLOCATE ( FIELD(1-NY:NY*(NX+2)) ) IF( RGLGRD .AND. .NOT. FLOMP ) ALLOCATE ( FIELD(1-NY:NY*(NX+2)) ) ! ! FIELD(:) = 0. ! LOCAL = IAPROC .LE. NAPROC UGDTUPDATE = .FALSE. IF (FLAGLL) THEN FACX = 1./(DERA * RADIUS) ELSE FACX = 1. END IF ! TAUWX = 0. TAUWY = 0. ! ! 0.d Test output ! ! 1. Check the consistency of the input ----------------------------- / ! 1.a Ending time versus initial time ! DTTST = DSEC21 ( TIME , TEND ) FLZERO = DTTST .EQ. 0. IF ( DTTST .LT. 0. ) THEN IF ( IAPROC .EQ. NAPERR ) WRITE (NDSE,1000) CALL EXTCDE ( 1 ) END IF ! ! 1.b Water level time ! IF ( FLLEV ) THEN IF ( TLEV(1) .GE. 0. ) THEN DTL0 = DSEC21 ( TLEV , TLN ) ELSE DTL0 = 1. END IF IF ( DTL0 .LT. 0. ) THEN IF ( IAPROC .EQ. NAPERR ) WRITE (NDSE,1001) CALL EXTCDE ( 2 ) END IF ELSE DTL0 = 0. END IF ! ! 1.c Current interval ! IF ( FLCUR ) THEN DTTST1 = DSEC21 ( TC0 , TCN ) DTTST2 = DSEC21 ( TC0 , TIME ) DTTST3 = DSEC21 ( TEND , TCN ) IF ( DTTST1.LT.0. .OR. DTTST2.LT.0. .OR. DTTST3.LT.0. ) THEN IF ( IAPROC .EQ. NAPERR ) WRITE (NDSE,1002) CALL EXTCDE ( 3 ) END IF IF ( DTTST2.EQ.0..AND. ITIME.EQ.0 ) THEN IDACT(7:7) = 'F' TOFRST = TIME END IF END IF ! ! 1.d Wind interval ! IF ( FLWIND ) THEN DTTST1 = DSEC21 ( TW0 , TWN ) DTTST2 = DSEC21 ( TW0 , TIME ) DTTST3 = DSEC21 ( TEND , TWN ) IF ( DTTST1.LT.0. .OR. DTTST2.LT.0. .OR. DTTST3.LT.0. ) THEN IF ( IAPROC .EQ. NAPERR ) WRITE (NDSE,1003) CALL EXTCDE ( 4 ) END IF IF ( DTTST2.EQ.0..AND. ITIME.EQ.0 ) THEN IDACT(3:3) = 'F' TOFRST = TIME END IF END IF ! ! 1.e Ice concentration interval ! IF ( FLICE ) THEN IF ( TICE(1) .GE. 0 ) THEN DTI0 = DSEC21 ( TICE , TIN ) ELSE DTI0 = 1. END IF IF ( DTI0 .LT. 0. ) THEN IF ( IAPROC .EQ. NAPERR ) WRITE (NDSE,1004) CALL EXTCDE ( 5 ) END IF ELSE DTI0 = 0. END IF ! ! 1.e Ice thickness interval ! IF ( FLIC1 ) THEN IF ( TIC1(1) .GE. 0 ) THEN DTI10 = DSEC21 ( TIC1 , TI1 ) ELSE DTI10 = 1. END IF IF ( DTI10 .LT. 0. ) THEN IF ( IAPROC .EQ. NAPERR ) WRITE (NDSE,1005) CALL EXTCDE ( 5 ) END IF ELSE DTI10 = 0. END IF ! ! 1.e Ice floe interval ! ! 2. Determine next time from ending and output --------------------- / ! time and get corresponding time step. ! FLFRST = .TRUE. DO ! DO JSEA = 1, NSEAL ! DO IS = 1, NSPEC ! IF (VA(IS, JSEA) .LT. 0.) THEN ! WRITE(740+IAPROC,*) 'TEST W3WAVE 2', VA(IS,JSEA) ! CALL FLUSH(740+IAPROC) ! ENDIF ! ENDDO ! ENDDO ! IF (SUM(VA) .NE. SUM(VA)) THEN ! WRITE(740+IAPROC,*) 'NAN in ACTION 2', IX, IY, SUM(VA) ! CALL FLUSH(740+IAPROC) ! STOP ! ENDIF ! ! 2.a Pre-calculate table for IC3 ------------------------------------ / ! 2.b Update group velocity and wavenumber from ice parameters ------- / ! from W3SIC3MD module. ------------------------------------------ / ! Note: "IF FLFRST" can be added for efficiency, but testing req'd JSEA=1 ! no switch (intentional) ISEA = IAPROC + (JSEA-1)*NAPROC ! 2.b.1 Using Cheng method: requires stationary/uniform rheology. ! However, ice thickness may be input by either method ! 2.b.2 If not using Cheng method: require FLIC1 to FLIC4 (not strictly ! necesssary, but makes code simpler) ! IF ( TOFRST(1) .GT. 0 ) THEN DTTST = DSEC21 ( TEND , TOFRST ) ELSE DTTST = 0. ENDIF ! IF ( DTTST.GE.0. ) THEN TCALC = TEND ELSE TCALC = TOFRST END IF ! DTTST = DSEC21 ( TIME , TCALC ) NT = 1 + INT ( DTTST / DTMAX - 0.001 ) DTGA = DTTST / REAL(NT) IF ( DTTST .EQ. 0. ) THEN IT0 = 0 IF ( .NOT.FLZERO ) ITIME = ITIME - 1 NT = 0 ELSE IT0 = 1 END IF ! ! ==================================================================== / ! ! 3. Loop over time steps ! DTRES = 0. ! DO IT=IT0, NT ! copy old values ! ITIME = ITIME + 1 DTG = REAL(NINT(DTGA+DTRES+0.0001)) DTRES = DTRES + DTGA - DTG IF ( ABS(DTRES) .LT. 0.001 ) DTRES = 0. CALL TICK21 ( TIME , DTG ) ! IF ( TSTAMP .AND. SCREEN.NE.NDSO .AND. IAPROC.EQ.NAPOUT ) THEN CALL WWTIME ( STTIME ) CALL STME21 ( TIME , IDTIME ) WRITE (SCREEN,950) IDTIME, STTIME END IF ! VGX = 0. VGY = 0. IF(INFLAGS1(8)) THEN DTTST1 = DSEC21 ( TIME, TGN ) DTTST2 = DSEC21 ( TG0, TGN ) FAC = DTTST1 / MAX ( 1. , DTTST2 ) VGX = (FAC*GA0+(1.-FAC)*GAN) * & COS(FAC*GD0+(1.-FAC)*GDN) VGY = (FAC*GA0+(1.-FAC)*GAN) * & SIN(FAC*GD0+(1.-FAC)*GDN) END IF ! ! 3.1 Interpolate winds and currents. ! (Initialize wave fields with winds) ! IF ( FLCUR ) THEN CALL W3UCUR ( FLFRST ) IF (GTYPE .NE. UNGTYPE) THEN IF( RGLGRD ) THEN CALL W3DZXY(CX(1:UBOUND(CX,1)),'m/s',DCXDX, DCXDY) !CX GRADIENT CALL W3DZXY(CY(1:UBOUND(CY,1)),'m/s',DCYDX, DCYDY) !CY GRADIENT ENDIF ELSE CALL UG_GRADIENTS(CX, DCXDX, DCXDY) CALL UG_GRADIENTS(CY, DCYDX, DCYDY) CALL GET_INTERFACE UGDTUPDATE=.TRUE. CFLXYMAX = 0. END IF ! ELSE IF ( FLFRST ) THEN UGDTUPDATE=.TRUE. CFLXYMAX = 0. CX = 0. CY = 0. END IF ! FLCUR ! IF ( FLWIND ) THEN IF ( FLFRST ) ASF = 1. CALL W3UWND ( FLFRST, VGX, VGY ) ELSE IF ( FLFRST ) THEN U10 = 0.01 U10D = 0. UST = 0.05 USTDIR = 0.05 END IF ! DO JSEA = 1, NSEAL ! DO IS = 1, NSPEC ! IF (VA(IS, JSEA) .LT. 0.) THEN ! WRITE(740+IAPROC,*) 'TEST W3WAVE 5', VA(IS,JSEA) ! CALL FLUSH(740+IAPROC) ! ENDIF ! ENDDO ! ENDDO ! IF (SUM(VA) .NE. SUM(VA)) THEN ! WRITE(740+IAPROC,*) 'NAN in ACTION 5', IX, IY, SUM(VA) ! CALL FLUSH(740+IAPROC) ! STOP ! ENDIF ! IF ( FLIWND .AND. LOCAL ) CALL W3UINI ( VA ) ! ! 3.2 Update boundary conditions if boundary flag is true (FLBPI) ! IF ( FLBPI .AND. LOCAL ) THEN ! DO IF ( TBPIN(1) .EQ. -1 ) THEN READBC = .TRUE. IDACT(1:1) = 'F' ELSE READBC = DSEC21(TIME,TBPIN).LT.0. IF (READBC.AND.IDACT(1:1).EQ.' ') IDACT(1:1) = 'X' END IF FLACT = READBC .OR. FLACT IF ( READBC ) THEN CALL W3IOBC ( 'READ', NDS(9), TBPI0, TBPIN, & ITEST, IMOD ) IF ( ITEST .NE. 1 ) CALL W3UBPT ELSE ITEST = 0 END IF IF ( ITEST .LT. 0 ) IDACT(1:1) = 'L' IF ( ITEST .GT. 0 ) IDACT(1:1) = ' ' IF ( .NOT. (READBC.AND.FLBPI) ) EXIT END DO END IF ! ! 3.3.1 Update ice coverage (if new ice map). ! Need to be run on output nodes too, to update MAPSTx ! IF ( FLICE .AND. DTI0.NE.0. ) THEN ! IF ( TICE(1).GE.0 ) THEN IF ( DTI0 .LT. 0. ) THEN IDACT(9:9) = 'B' ELSE DTTST = DSEC21 ( TIME, TIN ) IF ( DTTST .LE. 0.5*DTI0 ) IDACT(9:9) = 'U' END IF ELSE IDACT(9:9) = 'I' END IF ! IF ( IDACT(9:9).NE.' ' ) THEN CALL W3UICE ( VA, VA ) DTI0 = 0. FLACT = .TRUE. FLMAP = .TRUE. END IF END IF ! ! 3.3.2 Update ice thickness ! IF ( FLIC1 .AND. DTI10.NE.0. ) THEN ! IF ( TIC1(1).GE.0 ) THEN IF ( DTI10 .LT. 0. ) THEN IDACT(11:11) = 'B' ELSE DTTST = DSEC21 ( TIME, TI1 ) IF ( DTTST .LE. 0.5*DTI10 ) IDACT(11:11) = 'U' END IF ELSE IDACT(11:11) = 'I' END IF ! IF ( IDACT(11:11).NE.' ' ) THEN CALL W3UIC1 ( FLFRST ) DTI10 = 0. FLACT = .TRUE. FLMAP = .TRUE. END IF ! END IF ! ! 3.3.3 Update ice floe diameter ! ! ! 3.4 Transform grid (if new water level). ! ! write(740+IAPROC,*) 'TEST ARON', FLLEV, DTL0, TLEV(1), IDACT(5:5), DSEC21 ( TIME, TLN ), TIME, TLN IF ( FLLEV .AND. DTL0 .NE.0. ) THEN ! IF ( TLEV(1) .GE. 0 ) THEN IF ( DTL0 .LT. 0. ) THEN IDACT(5:5) = 'B' ELSE DTTST = DSEC21 ( TIME, TLN ) IF ( DTTST .LE. 0.5*DTL0 ) IDACT(5:5) = 'U' END IF ELSE IDACT(5:5) = 'I' END IF ! IF ( IDACT(5:5).NE.' ' ) THEN CALL W3ULEV ( VA, VA ) UGDTUPDATE=.TRUE. CFLXYMAX = 0. DTL0 = 0. FLACT = .TRUE. FLMAP = .TRUE. FLDDIR = FLDDIR .OR. FLCTH .OR. FSREFRACTION & .OR. FLCK .OR. FSFREQSHIFT END IF END IF ! ! 3.5 Update maps and derivatives. ! IF ( FLMAP .AND. RGLGRD ) THEN CALL W3MAP3 CALL W3UTRN ( TRNX, TRNY ) CALL W3MAPT CALL W3NMIN ( MAPSTA, FLAG0 ) IF ( FLAG0 .AND. IAPROC.EQ.NAPERR ) WRITE (NDSE,1030) IMOD FLMAP = .FALSE. END IF ! IF ( FLDDIR ) THEN IF (GTYPE .NE. UNGTYPE) THEN !!Li CALL W3DZXY(DW(1:NSEA),'m',DDDX,DDDY) !DEPTH (DW) GRADIENT IF( RGLGRD ) CALL W3DZXY(DW(1:UBOUND(DW,1)),'m',DDDX,DDDY) ! ELSE CALL UG_GRADIENTS(DW, DDDX, DDDY) END IF FLDDIR = .FALSE. END IF ! ! Calculate PHASE SPEED GRADIENT. DCDX = 0. DCDY = 0. ! FLIWND = .FALSE. FLFRST = .FALSE. ! ! IF ( FLSOU .and. LPDLIB) THEN ! D50=0.0002 REFLEC(:)=0. REFLED(:)=0 PSIC=0. DO JSEA=1, NSEAL CALL INIT_GET_ISEA(ISEA, JSEA) IX = MAPSF(ISEA,1) IY = MAPSF(ISEA,2) DELA=1. DELX=1. DELY=1. IF (GTYPE.EQ.RLGTYPE) THEN DELX=SX*CLATS(ISEA)/FACX DELY=SY/FACX DELA=DELX*DELY END IF IF (GTYPE.EQ.CLGTYPE) THEN DELX=HPFAC(IY,IX)/ FACX DELY=HQFAC(IY,IX)/ FACX DELA=DELX*DELY END IF ! REFLEC=REFLC(:,ISEA) REFLEC(4)=BERG(ISEA)*REFLEC(4) REFLED=REFLD(:,ISEA) D50=SED_D50(ISEA) PSIC=SED_PSIC(ISEA) REFLEC=REFLC(:,ISEA) REFLEC(4)=BERG(ISEA)*REFLEC(4) REFLED=REFLD(:,ISEA) D50=SED_D50(ISEA) PSIC=SED_PSIC(ISEA) ! IF ( MAPSTA(IY,IX) .EQ. 1 .AND. FLAGST(ISEA)) THEN IF (FSSOURCE) THEN ENDIF ELSE UST (ISEA) = UNDEF USTDIR(ISEA) = UNDEF DTDYN (JSEA) = UNDEF FCUT (JSEA) = UNDEF END IF END DO ! JSEA END IF ! PDLIB IF ( FLZERO ) THEN GOTO 400 END IF IF ( IT.EQ.0 ) THEN DTG = 1. ! DTG = 60. GOTO 370 END IF IF ( FLDRY .OR. IAPROC.GT.NAPROC ) THEN GOTO 380 END IF ! ! Estimation of the local maximum CFL for XY propagation ! IF ( FLOGRD(9,3).AND. UGDTUPDATE ) THEN IF (FSTOTALIMP .eqv. .FALSE.) THEN NKCFL=NK ! DO JSEA=1, NSEAL CALL INIT_GET_ISEA(ISEA, JSEA) IF (GTYPE .EQ. UNGTYPE) THEN IF ( FLOGRD(9,3) ) THEN CALL W3CFLUG ( ISEA, NKCFL, FACX, FACX, DTG, & MAPFS, CFLXYMAX(JSEA), VGX, VGY ) END IF ELSE CALL W3CFLXY ( ISEA, DTG, MAPSTA, MAPFS, & CFLXYMAX(JSEA), VGX, VGY ) END IF END DO ! END IF END IF ! ! ! ! 3.6 Perform Propagation = = = = = = = = = = = = = = = = = = = = = = = ! 3.6.1 Preparations ! NTLOC = 1 + INT( DTG/DTCFLI - 0.001 ) ! FACTH = DTG / (DTH*REAL(NTLOC)) FACK = DTG / REAL(NTLOC) ITLOCH = ( NTLOC + 1 - MOD(ITIME,2) ) / 2 ! ! 3.6.2 Intra-spectral part 1 ! IF ( FLCTH .OR. FLCK ) THEN DO ITLOC=1, ITLOCH ! DO JSEA=1, NSEAL CALL INIT_GET_ISEA(ISEA, JSEA) IX = MAPSF(ISEA,1) IY = MAPSF(ISEA,2) IF ( GTYPE .EQ. UNGTYPE ) THEN IF (IOBP(ISEA) .NE. 1) CYCLE ENDIF IF ( MAPSTA(IY,IX) .EQ. 1 ) THEN DEPTH = MAX ( DMIN , DW(ISEA) ) IF (LPDLIB) THEN IXrel = JSEA ELSE IXrel = IX END IF CALL W3KTP3 ( ISEA, FACTH, FACK, CTHG0S(ISEA), & CG(:,ISEA), WN(:,ISEA), DEPTH, & DDDX(IY,IXrel), DDDY(IY,IXrel), CX(ISEA), & CY(ISEA), DCXDX(IY,IXrel), DCXDY(IY,IXrel), & DCYDX(IY,IXrel), DCYDY(IY,IXrel), & DCDX(:,IY,IXrel), DCDY(:,IY,IXrel), VA(:,JSEA), & CFLTHMAX(JSEA), CFLKMAX(JSEA) ) ! END IF END DO ! END DO END IF ! ! 3.6.3 Longitude-latitude ! (time step correction in routine) ! IF (GTYPE .EQ. UNGTYPE) THEN IF (FLAGLL) THEN FACX = 1./(DERA * RADIUS) ELSE FACX = 1. END IF END IF IF ((GTYPE .EQ. UNGTYPE) .and. LPDLIB) THEN ! ELSE IF (FLCX .or. FLCY) THEN ! IF ( NRQSG1 .GT. 0 ) THEN CALL MPI_STARTALL (NRQSG1, IRQSG1(1,1), IERR_MPI) CALL MPI_STARTALL (NRQSG1, IRQSG1(1,2), IERR_MPI) END IF ! IF ( FLOMP ) ALLOCATE ( FIELD(1-NY:NY*(NX+2)) ) ! DO ISPEC=1, NSPEC IF ( IAPPRO(ISPEC) .EQ. IAPROC ) THEN IF (.NOT.LPDLIB .AND. RGLGRD) CALL W3GATH ( ISPEC, FIELD ) ! IF (GTYPE .NE. UNGTYPE) THEN CALL W3XYP3 ( ISPEC, DTG, MAPSTA, MAPFS, FIELD, VGX, VGY ) ! ELSE IF (GTYPE .EQ. UNGTYPE) THEN IF (.NOT. LPDLIB) THEN CALL W3XYPUG ( ISPEC, FACX, FACX, DTG, & FIELD, VGX, VGY, UGDTUPDATE ) END IF END IF IF (.NOT.LPDLIB .AND. RGLGRD) CALL W3SCAT ( ISPEC, MAPSTA, FIELD ) ! END IF END DO ! IF ( NRQSG1 .GT. 0 ) THEN ALLOCATE ( STATCO(MPI_STATUS_SIZE,NRQSG1) ) CALL MPI_WAITALL (NRQSG1, IRQSG1(1,1), STATCO, & IERR_MPI) CALL MPI_WAITALL (NRQSG1, IRQSG1(1,2), STATCO, & IERR_MPI) DEALLOCATE ( STATCO ) END IF IF ( FLOMP ) DEALLOCATE ( FIELD ) !Li Initialise IK IX IY in case ARC option is not used to avoid warnings. IK=1 IX=1 IY=1 ! IF( ARCTIC ) THEN ISPEC = MOD( IY-1, NAPROC ) JSEA = 1 + (IY - ISPEC - 1)/NAPROC ENDIF ! !!Li Broadcast local SPCBAC(:,IK) to all other PEs. IF( ARCTIC ) THEN CALL MPI_BCAST(SPCBAC(1,IK),NSPEC,MPI_REAL,ISPEC,MPI_COMM_WAVE,IERR_MPI) CALL MPI_BARRIER (MPI_COMM_WAVE,IERR_MPI) ENDIF ! IF( ARCTIC ) THEN ISPEC = MOD( IX-1, NAPROC ) JSEA = 1 + (IX - ISPEC - 1)/NAPROC ENDIF ! ! End of test FLCX.OR.FLCY ! END IF END IF ! ! 3.6.4 Intra-spectral part 2 ! IF ( FLCTH .OR. FLCK ) THEN DO ITLOC=ITLOCH+1, NTLOC ! DO JSEA = 1, NSEAL CALL INIT_GET_ISEA(ISEA, JSEA) IX = MAPSF(ISEA,1) IY = MAPSF(ISEA,2) DEPTH = MAX ( DMIN , DW(ISEA) ) IF ( GTYPE .EQ. UNGTYPE ) THEN IF (IOBP(ISEA) .NE. 1) CYCLE ENDIF IF ( MAPSTA(IY,IX) .EQ. 1 ) THEN IF (LPDLIB) THEN IXrel = JSEA ELSE IXrel = IX END IF CALL W3KTP3 ( ISEA, FACTH, FACK, CTHG0S(ISEA), & CG(:,ISEA), WN(:,ISEA), DEPTH, & DDDX(IY,IXrel), DDDY(IY,IXrel), CX(ISEA), & CY(ISEA), DCXDX(IY,IXrel), DCXDY(IY,IXrel), & DCYDX(IY,IXrel), DCYDY(IY,IXrel), & DCDX(:,IY,IXrel), DCDY(:,IY,IXrel), VA(:,JSEA), & CFLTHMAX(JSEA), CFLKMAX(JSEA) ) ! END IF END DO ! END DO END IF ! UGDTUPDATE = .FALSE. ! ! 3.6 End propapgation = = = = = = = = = = = = = = = = = = = = = = = = ! 3.7 Calculate and integrate source terms. ! 370 CONTINUE IF ( FLSOU ) THEN ! D50=0.0002 REFLEC(:)=0. REFLED(:)=0 PSIC=0. ! DO JSEA=1, NSEAL CALL INIT_GET_ISEA(ISEA, JSEA) IX = MAPSF(ISEA,1) IY = MAPSF(ISEA,2) DELA=1. DELX=1. DELY=1. IF (GTYPE.EQ.RLGTYPE) THEN DELX=SX*CLATS(ISEA)/FACX DELY=SY/FACX DELA=DELX*DELY END IF IF (GTYPE.EQ.CLGTYPE) THEN DELX=HPFAC(IY,IX)/ FACX DELY=HQFAC(IY,IX)/ FACX DELA=DELX*DELY END IF ! REFLEC=REFLC(:,ISEA) REFLEC(4)=BERG(ISEA)*REFLEC(4) REFLED=REFLD(:,ISEA) D50=SED_D50(ISEA) PSIC=SED_PSIC(ISEA) IF ( MAPSTA(IY,IX) .EQ. 1 .AND. FLAGST(ISEA)) THEN TMP1 = WHITECAP(JSEA,1:4) TMP2 = BEDFORMS(JSEA,1:3) TMP3 = TAUBBL(JSEA,1:2) TMP4 = TAUICE(JSEA,1:2) CALL W3SRCE(srce_direct, IT, JSEA, IX, IY, IMOD, & VAoldDummy, VA(:,JSEA), & VSioDummy, VDioDummy, SHAVETOTioDummy, & ALPHA(1:NK,JSEA), WN(1:NK,ISEA), & CG(1:NK,ISEA), DW(ISEA), U10(ISEA), & U10D(ISEA), AS(ISEA), UST(ISEA), & USTDIR(ISEA), CX(ISEA), CY(ISEA), & ICE(ISEA), ICEH(ISEA), ICEF(ISEA), & ICEDMAX(ISEA), & REFLEC, REFLED, DELX, DELY, DELA, & TRNX(IY,IX), TRNY(IY,IX), BERG(ISEA), & FPIS(ISEA), DTDYN(JSEA), & FCUT(JSEA), DTG, TAUWX(JSEA), TAUWY(JSEA), & TAUOX(JSEA), TAUOY(JSEA), TAUWIX(JSEA), & TAUWIY(JSEA), TAUWNX(JSEA), & TAUWNY(JSEA), PHIAW(JSEA), CHARN(JSEA), & TWS(JSEA), PHIOC(JSEA), TMP1, D50, PSIC,TMP2,& PHIBBL(JSEA), TMP3, TMP4 , PHICE(JSEA), & ASF(ISEA)) WHITECAP(JSEA,1:4) = TMP1 BEDFORMS(JSEA,1:3) = TMP2 TAUBBL(JSEA,1:2) = TMP3 TAUICE(JSEA,1:2) = TMP4 ELSE UST (ISEA) = UNDEF USTDIR(ISEA) = UNDEF DTDYN (JSEA) = UNDEF FCUT (JSEA) = UNDEF ! VA(:,JSEA) = 0. END IF END DO ! ! ! This barrier is from older code versions. It has been removed in 3.11 ! to optimize IO2/3 settings. May be needed on some systems still ! !!/MPI IF (FLAG0) CALL MPI_BARRIER (MPI_COMM_WCMP,IERR_MPI) !!/MPI ELSE !!/MPI CALL MPI_BARRIER (MPI_COMM_WCMP,IERR_MPI) ! END IF ! ! End of interations for DTMAX < 1s ! ! 3.8 Update global time step. ! (Branch point FLDRY, IT=0) ! 380 CONTINUE ! IF (IT.NE.NT) THEN DTTST = DSEC21 ( TIME , TCALC ) DTG = DTTST / REAL(NT-IT) END IF ! IF ( FLACT .AND. IT.NE.NT .AND. IAPROC.EQ.NAPLOG ) THEN CALL STME21 ( TIME , IDTIME ) IF ( IDLAST .NE. TIME(1) ) THEN WRITE (NDSO,900) ITIME, IPASS, IDTIME(01:19), & IDACT, OUTID IDLAST = TIME(1) ELSE WRITE (NDSO,901) ITIME, IPASS, IDTIME(12:19), & IDACT, OUTID END IF FLACT = .FALSE. IDACT = ' ' END IF ! END DO ! ! End of loop over time steps ! ==================================================================== / ! 400 CONTINUE ! ! 4. Perform output to file if requested ---------------------------- / ! 4.a Check if time is output time ! Delay if data assimilation time. ! IF ( TOFRST(1) .EQ. -1 ) THEN DTTST = 1. ELSE DTTST = DSEC21 ( TIME, TOFRST ) END IF ! IF ( TDN(1) .EQ. -1 ) THEN DTTST1 = 1. ELSE DTTST1 = DSEC21 ( TIME, TDN ) END IF ! DTTST2 = DSEC21 ( TIME, TEND ) FLAG_O = .NOT.SKIP_O .OR. ( SKIP_O .AND. DTTST2.NE.0. ) ! IF ( DTTST.LE.0. .AND. DTTST1.NE.0. .AND. FLAG_O ) THEN ! ! 4.b Processing and MPP preparations ! IF ( FLOUT(1) ) THEN FLOUTG = DSEC21(TIME,TONEXT(:,1)).EQ.0. ELSE FLOUTG = .FALSE. END IF ! IF ( FLOUT(7) ) THEN FLOUTG2 = DSEC21(TIME,TONEXT(:,7)).EQ.0. ELSE FLOUTG2 = .FALSE. END IF ! FLPART = .FALSE. IF ( FLOUT(1) .AND. FLPFLD ) & FLPART = FLPART .OR. DSEC21(TIME,TONEXT(:,1)).EQ.0. IF ( FLOUT(6) ) & FLPART = FLPART .OR. DSEC21(TIME,TONEXT(:,6)).EQ.0. ! IF ( LOCAL .AND. FLPART ) CALL W3CPRT ( IMOD ) IF ( LOCAL .AND. (FLOUTG .OR. FLOUTG2) ) & CALL W3OUTG ( VA, FLPFLD, FLOUTG, FLOUTG2 ) ! FLGMPI = .FALSE. NRQMAX = 0 ! IF ( ( ( DSEC21(TIME,TONEXT(:,1)).EQ.0. ) .OR. & ( DSEC21(TIME,TONEXT(:,7)).EQ.0. ) ).and. & (FLOUT(1) .OR. FLOUT(7)) ) THEN IF (.NOT. LPDLIB .or. (GTYPE.ne.UNGTYPE)) THEN IF (NRQGO.NE.0 ) THEN CALL MPI_STARTALL ( NRQGO, IRQGO , IERR_MPI ) FLGMPI(0) = .TRUE. NRQMAX = MAX ( NRQMAX , NRQGO ) END IF ! IF (NRQGO2.NE.0 ) THEN CALL MPI_STARTALL ( NRQGO2, IRQGO2, IERR_MPI ) FLGMPI(1) = .TRUE. NRQMAX = MAX ( NRQMAX , NRQGO2 ) END IF ELSE END IF END IF ! IF ( FLOUT(2) .AND. NRQPO.NE.0 ) THEN IF ( DSEC21(TIME,TONEXT(:,2)).EQ.0. ) THEN CALL MPI_STARTALL ( NRQPO, IRQPO1, IERR_MPI ) FLGMPI(2) = .TRUE. NRQMAX = MAX ( NRQMAX , NRQPO ) END IF END IF ! IF ( FLOUT(4) .AND. NRQRS.NE.0 ) THEN IF ( DSEC21(TIME,TONEXT(:,4)).EQ.0. ) THEN CALL MPI_STARTALL ( NRQRS, IRQRS , IERR_MPI ) FLGMPI(4) = .TRUE. NRQMAX = MAX ( NRQMAX , NRQRS ) END IF END IF ! IF ( FLOUT(5) .AND. NRQBP.NE.0 ) THEN IF ( DSEC21(TIME,TONEXT(:,5)).EQ.0. ) THEN CALL MPI_STARTALL ( NRQBP , IRQBP1, IERR_MPI ) FLGMPI(5) = .TRUE. NRQMAX = MAX ( NRQMAX , NRQBP ) END IF END IF ! IF ( FLOUT(5) .AND. NRQBP2.NE.0 .AND. & IAPROC.EQ.NAPBPT) THEN IF ( DSEC21(TIME,TONEXT(:,5)).EQ.0. ) THEN CALL MPI_STARTALL (NRQBP2,IRQBP2,IERR_MPI) NRQMAX = MAX ( NRQMAX , NRQBP2 ) END IF END IF ! IF ( NRQMAX .NE. 0 ) ALLOCATE & ( STATIO(MPI_STATUS_SIZE,NRQMAX) ) ! ! 4.c Reset next output time ! TOFRST(1) = -1 TOFRST(2) = 0 ! DO J=1, NOTYPE IF ( FLOUT(J) ) THEN ! ! 4.d Perform output ! TOUT(:) = TONEXT(:,J) DTTST = DSEC21 ( TIME, TOUT ) ! IF ( DTTST .EQ. 0. ) THEN IF ( ( J .EQ. 1 ) .OR. ( J .EQ. 7 ) ) THEN IF ( IAPROC .EQ. NAPFLD ) THEN IF ( FLGMPI(1) ) CALL MPI_WAITALL & ( NRQGO2, IRQGO2, STATIO, IERR_MPI ) FLGMPI(1) = .FALSE. IF ( J .EQ. 1 ) CALL W3IOGO & ( 'WRITE', NDS(7), ITEST, IMOD ) END IF IF ( J .EQ. 7 ) THEN ! ! Send variables to atmospheric or ocean circulation or ice model ! IF (DTOUT(7).NE.0) THEN IF ( (MOD(ID_OASIS_TIME/DTOUT(7),1.0) .LT. 1.E-7 ) .AND. & (DSEC21 (TIME00, TIME) .GT. 0.0) ) THEN CALL SND_FIELDS_TO_ATMOS() CALL SND_FIELDS_TO_OCEAN() ID_OASIS_TIME = DSEC21 ( TIME00 , TIME ) ENDIF END IF END IF ELSE IF ( J .EQ. 2 ) THEN ! ! Point output ! IF ( IAPROC .EQ. NAPPNT ) THEN ! ! Gets the necessary spectral data ! CALL W3IOPE ( VA ) CALL W3IOPO ( 'WRITE', NDS(8), ITEST, IMOD ) END IF ! ELSE IF ( J .EQ. 3 ) THEN ! ! Track output ! CALL W3IOTR ( NDS(11), NDS(12), VA, IMOD ) ELSE IF ( J .EQ. 4 ) THEN CALL W3IORS ('HOT', NDS(6), XXX, ITEST, IMOD ) ELSE IF ( J .EQ. 5 ) THEN IF ( IAPROC .EQ. NAPBPT ) THEN IF (NRQBP2.NE.0) CALL MPI_WAITALL & ( NRQBP2, IRQBP2,STATIO, IERR_MPI ) CALL W3IOBC ( 'WRITE', NDS(10), & TIME, TIME, ITEST, IMOD ) END IF ELSE IF ( J .EQ. 6 ) THEN CALL W3IOSF ( NDS(13), IMOD ) END IF ! CALL TICK21 ( TOUT, DTOUT(J) ) TONEXT(:,J) = TOUT TLST = TOLAST(:,J) DTTST = DSEC21 ( TOUT , TLST ) FLOUT(J) = DTTST.GE.0. IF ( FLOUT(J) ) THEN OUTID(2*J-1:2*J-1) = 'X' IF ( (DTOUT(7).NE.0) .AND. & (DSEC21(TIME,TIME00).EQ.0 .OR. & DSEC21(TIME,TIMEEND).EQ.0) ) OUTID(13:13) = ' ' ELSE OUTID(2*J-1:2*J-1) = 'L' END IF END IF ! ! 4.e Update next output time ! IF ( FLOUT(J) ) THEN IF ( TOFRST(1).EQ.-1 ) THEN TOFRST = TOUT ELSE DTTST = DSEC21 ( TOUT , TOFRST ) IF ( DTTST.GT.0.) THEN TOFRST = TOUT END IF END IF END IF ! END IF ! END DO ! IF ( FLGMPI(0) ) CALL MPI_WAITALL & ( NRQGO, IRQGO , STATIO, IERR_MPI ) IF ( FLGMPI(2) ) CALL MPI_WAITALL & ( NRQPO, IRQPO1, STATIO, IERR_MPI ) IF ( FLGMPI(4) ) CALL MPI_WAITALL & ( NRQRS, IRQRS , STATIO, IERR_MPI ) IF ( FLGMPI(5) ) CALL MPI_WAITALL & ( NRQBP, IRQBP1, STATIO, IERR_MPI ) IF ( NRQMAX .NE. 0 ) DEALLOCATE ( STATIO ) ! ! This barrier is from older code versions. It has been removed in 3.11 ! to optimize IO2/3 settings. May be needed on some systems still ! !!/MPI IF (FLDRY) CALL MPI_BARRIER (MPI_COMM_WAVE,IERR_MPI) ! END IF ! ! 5. Update log file ------------------------------------------------ / ! IF (MINVAL(VA) .LT. 0.) THEN ! WRITE(740+IAPROC,*) 'TEST W3WAVE 13', SUM(VA), MINVAL(VA), MAXVAL(VA) ! CALL FLUSH(740+IAPROC) ! STOP ! ENDIF ! IF ( IAPROC.EQ.NAPLOG ) THEN ! CALL STME21 ( TIME , IDTIME ) IF ( FLCUR ) THEN DTTST = DSEC21 ( TIME , TCN ) IF ( DTTST .EQ. 0. ) IDACT(7:7) = 'X' END IF IF ( FLWIND ) THEN DTTST = DSEC21 ( TIME , TWN ) IF ( DTTST .EQ. 0. ) IDACT(3:3) = 'X' END IF IF ( TDN(1) .GT. 0 ) THEN DTTST = DSEC21 ( TIME , TDN ) IF ( DTTST .EQ. 0. ) IDACT(17:17) = 'X' END IF ! IF ( IDLAST.NE.TIME(1) ) THEN WRITE (NDSO,900) ITIME, IPASS, IDTIME(1:19), & IDACT, OUTID IDLAST = TIME(1) ELSE WRITE (NDSO,901) ITIME, IPASS, IDTIME(12:19), & IDACT, OUTID END IF ! END IF ! IDACT = ' ' OUTID = ' ' FLACT = .FALSE. ! ! 6. If time is not ending time, branch back to 2 ------------------- / ! DTTST = DSEC21 ( TIME, TEND ) IF ( DTTST .EQ. 0. ) EXIT END DO ! IF ( TSTAMP .AND. SCREEN.NE.NDSO .AND. IAPROC.EQ.NAPOUT ) THEN CALL WWTIME ( STTIME ) WRITE (SCREEN,951) STTIME END IF IF ( IAPROC .EQ. NAPLOG ) WRITE (NDSO,902) ! IF ( .NOT. FLOMP ) DEALLOCATE ( FIELD ) ! DEALLOCATE(TAUWX, TAUWY) ! RETURN ! ! Formats ! 900 FORMAT (4X,I6,'|',I6,'| ', A19 ,' | ',A,' | ',A,' |') 901 FORMAT (4X,I6,'|',I6,'| ',11X,A8,' | ',A,' | ',A,' |') 902 FORMAT (2X,'--------+------+---------------------+' & ,'-------------------+---------------+') ! 950 FORMAT (' WAVEWATCH III calculating for ',A,' at ',A) 951 FORMAT (' WAVEWATCH III reached the end of a computation', & ' loop at ',A) 1000 FORMAT (/' *** WAVEWATCH III ERROR IN W3WAVE :'/ & ' ENDING TIME BEFORE STARTING TIME '/) 1001 FORMAT (/' *** WAVEWATCH III ERROR IN W3WAVE :'/ & ' NEW WATER LEVEL BEFORE OLD WATER LEVEL '/) 1002 FORMAT (/' *** WAVEWATCH III ERROR IN W3WAVE :'/ & ' ILLEGAL CURRENT INTERVAL '/) 1003 FORMAT (/' *** WAVEWATCH III ERROR IN W3WAVE :'/ & ' ILLEGAL WIND INTERVAL '/) 1004 FORMAT (/' *** WAVEWATCH III ERROR IN W3WAVE :'/ & ' NEW ICE FIELD BEFORE OLD ICE FIELD '/) 1005 FORMAT (/' *** WAVEWATCH III ERROR IN W3WAVE :'/ & ' NEW IC1 FIELD BEFORE OLD IC1 FIELD '/) 1030 FORMAT (/' *** WAVEWATCH III WARING IN W3WAVE :'/ & ' AT LEAST ONE PROCESSOR HAS 0 ACTIVE POINTS', & ' IN GRID',I3) ! !/ !/ End of W3WAVE ----------------------------------------------------- / !/ END SUBROUTINE W3WAVE !/ ------------------------------------------------------------------- / SUBROUTINE W3GATH ( ISPEC, FIELD ) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | FORTRAN 90 | !/ | Last update : 26-Dec-2012 | !/ +-----------------------------------+ !/ !/ 04-Jan-1999 : Distributed FORTRAN 77 version. ( version 1.18 ) !/ 13-Jan-2000 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ Major changes to logistics. !/ 29-Dec-2004 : Multiple grid version. ( version 3.06 ) !/ 13-Jun-2006 : Split STORE in G/SSTORE ( version 3.09 ) !/ 26-Dec-2012 : Move FIELD init. to W3GATH. ( version 4.OF ) !/ ! 1. Purpose : ! ! Gather spectral bin information into a propagation field array. ! ! 2. Method : ! ! Direct copy or communication calls (MPP version). ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ISPEC Int. I Spectral bin considered. ! FIELD R.A. O Full field to be propagated. ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ! MPI_STARTALL, MPI_WAITALL ! Subr. mpif.h MPI persistent comm. routines (!/MPI). ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3WAVE Subr. W3WAVEMD Actual wave model routine. ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! - The field is extracted but not converted. ! - MPI version requires posing of send and receive calls in ! W3WAVE to match local calls. ! - MPI version does not require an MPI_TESTALL call for the ! posted gather operation as MPI_WAITALL is mandatory to ! reset persistent communication for next time step. ! - MPI version allows only two new pre-fetch postings per ! call to minimize chances to be slowed down by gathers that ! are not yet needed, while maximizing the pre-loading ! during the early (low-frequency) calls to the routine ! where the amount of calculation needed for proagation is ! the largest. ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/SHRD Switch for message passing method. ! !/MPI Id. ! ! !/S Enable subroutine tracing. ! !/MPIT MPI test output. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/ USE W3GDATMD, ONLY: NSPEC, NX, NY, NSEA, NSEAL, MAPSF, DMIN USE W3PARALL, ONLY: INIT_GET_ISEA USE W3WDATMD, ONLY: A => VA USE W3ADATMD, ONLY: MPIBUF, BSTAT, IBFLOC, ISPLOC, BISPL, & NSPLOC, NRQSG2, IRQSG2, GSTORE USE W3ODATMD, ONLY: NDST, IAPROC, NAPROC, NOTYPE !/ IMPLICIT NONE ! INCLUDE "mpif.h" !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ INTEGER, INTENT(IN) :: ISPEC REAL, INTENT(OUT) :: FIELD(1-NY:NY*(NX+2)) !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ INTEGER :: STATUS(MPI_STATUS_SIZE,NSPEC), & IOFF, IERR_MPI, JSEA, ISEA, & IXY, IS0, IB0, NPST, J !/ !/ ------------------------------------------------------------------- / !/ ! FIELD = 0. ! ! 1. Shared memory version ------------------------------------------ / ! ! 2. Distributed memory version ( MPI ) ----------------------------- / ! 2.a Update counters ! ISPLOC = ISPLOC + 1 IBFLOC = IBFLOC + 1 IF ( IBFLOC .GT. MPIBUF ) IBFLOC = 1 ! ! 2.b Check status of present buffer ! 2.b.1 Scatter (send) still in progress, wait to end ! IF ( BSTAT(IBFLOC) .EQ. 2 ) THEN IOFF = 1 + (BISPL(IBFLOC)-1)*NRQSG2 IF ( NRQSG2 .GT. 0 ) CALL & MPI_WAITALL ( NRQSG2, IRQSG2(IOFF,2), & STATUS, IERR_MPI ) BSTAT(IBFLOC) = 0 END IF ! ! 2.b.2 Gather (recv) not yet posted, post now ! IF ( BSTAT(IBFLOC) .EQ. 0 ) THEN BSTAT(IBFLOC) = 1 BISPL(IBFLOC) = ISPLOC IOFF = 1 + (ISPLOC-1)*NRQSG2 IF ( NRQSG2 .GT. 0 ) CALL & MPI_STARTALL ( NRQSG2, IRQSG2(IOFF,1), IERR_MPI ) END IF ! ! 2.c Put local spectral densities in store ! DO JSEA=1, NSEAL CALL INIT_GET_ISEA(ISEA, JSEA) GSTORE(ISEA,IBFLOC) = A(ISPEC,JSEA) END DO ! ! 2.d Wait for remote spectral densities ! IOFF = 1 + (BISPL(IBFLOC)-1)*NRQSG2 IF ( NRQSG2 .GT. 0 ) CALL & MPI_WAITALL ( NRQSG2, IRQSG2(IOFF,1), STATUS, IERR_MPI ) ! ! 2.e Convert storage array to field. ! DO ISEA=1, NSEA IXY = MAPSF(ISEA,3) FIELD(IXY) = GSTORE(ISEA,IBFLOC) END DO ! ! 2.f Pre-fetch data in available buffers ! IS0 = ISPLOC IB0 = IBFLOC NPST = 0 ! DO J=1, MPIBUF-1 IS0 = IS0 + 1 IF ( IS0 .GT. NSPLOC ) EXIT IB0 = 1 + MOD(IB0,MPIBUF) IF ( BSTAT(IB0) .EQ. 0 ) THEN BSTAT(IB0) = 1 BISPL(IB0) = IS0 IOFF = 1 + (IS0-1)*NRQSG2 IF ( NRQSG2 .GT. 0 ) CALL & MPI_STARTALL ( NRQSG2, IRQSG2(IOFF,1), IERR_MPI ) NPST = NPST + 1 END IF IF ( NPST .GE. 2 ) EXIT END DO ! ! 2.g Test output ! RETURN ! ! Formats ! !/ !/ End of W3GATH ----------------------------------------------------- / !/ END SUBROUTINE W3GATH !/ ------------------------------------------------------------------- / SUBROUTINE W3SCAT ( ISPEC, MAPSTA, FIELD ) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | FORTRAN 90 | !/ | Last update : 13-Jun-2006 | !/ +-----------------------------------+ !/ !/ 04-Jan-1999 : Distributed FORTRAN 77 version. ( version 1.18 ) !/ 13-Jan-2000 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ Major changes to logistics. !/ 28-Dec-2004 : Multiple grid version. ( version 3.06 ) !/ 07-Sep-2005 : Updated boundary conditions. ( version 3.08 ) !/ 13-Jun-2006 : Split STORE in G/SSTORE ( version 3.09 ) !/ ! 1. Purpose : ! ! 'Scatter' data back to spectral storage after propagation. ! ! 2. Method : ! ! Direct copy or communication calls (MPP version). ! See also W3GATH. ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ISPEC Int. I Spectral bin considered. ! MAPSTA I.A. I Status map for spatial grid. ! FIELD R.A. I Full field to be propagated. ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ! MPI_STARTALL, MPI_WAITALL, MPI_TESTALL ! Subr. mpif.h MPI persistent comm. routines (!/MPI). ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3WAVE Subr. W3WAVEMD Actual wave model routine. ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! - The field is put back but not converted ! ! - MPI persistent communication calls initialize in W3MPII. ! - See W3GATH and W3MPII for additional comments on data ! buffering. ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/SHRD Switch for message passing method. ! !/MPI Id. ! ! !/S Enable subroutine tracing. ! !/MPIT MPI test output. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/ USE W3GDATMD, ONLY: NSPEC, NX, NY, NSEA, NSEAL, MAPSF USE W3WDATMD, ONLY: A => VA USE W3ADATMD, ONLY: MPIBUF, BSTAT, IBFLOC, ISPLOC, BISPL, & NSPLOC, NRQSG2, IRQSG2, SSTORE USE W3ODATMD, ONLY: NDST USE W3ODATMD, ONLY: IAPROC, NAPROC USE CONSTANTS, ONLY : LPDLIB USE W3PARALL, only: INIT_GET_ISEA !/ IMPLICIT NONE ! INCLUDE "mpif.h" !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ INTEGER, INTENT(IN) :: ISPEC, MAPSTA(NY*NX) REAL, INTENT(IN) :: FIELD(1-NY:NY*(NX+2)) !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ INTEGER :: ISEA, IXY, IOFF, IERR_MPI, J, & STATUS(MPI_STATUS_SIZE,NSPEC), & JSEA, IB0 LOGICAL :: DONE !/ !/ ------------------------------------------------------------------- / !/ ! ! 1. Shared memory version ------------------------------------------ * ! ! 2. Distributed memory version ( MPI ) ----------------------------- * ! 2.a Initializations ! ! 2.b Convert full grid to sea grid, active points only ! DO ISEA=1, NSEA IXY = MAPSF(ISEA,3) IF ( MAPSTA(IXY) .GE. 1 ) SSTORE(ISEA,IBFLOC) = FIELD(IXY) END DO ! ! 2.c Send spectral densities to appropriate remote ! IOFF = 1 + (ISPLOC-1)*NRQSG2 IF ( NRQSG2 .GT. 0 ) CALL & MPI_STARTALL ( NRQSG2, IRQSG2(IOFF,2), IERR_MPI ) BSTAT(IBFLOC) = 2 ! ! 2.d Save locally stored results ! DO JSEA=1, NSEAL CALL INIT_GET_ISEA(ISEA, JSEA) IXY = MAPSF(ISEA,3) IF (MAPSTA(IXY) .GE. 1) A(ISPEC,JSEA) = SSTORE(ISEA,IBFLOC) END DO ! ! 2.e Check if any sends have finished ! IB0 = IBFLOC ! DO J=1, MPIBUF IB0 = 1 + MOD(IB0,MPIBUF) IF ( BSTAT(IB0) .EQ. 2 ) THEN IOFF = 1 + (BISPL(IB0)-1)*NRQSG2 IF ( NRQSG2 .GT. 0 ) THEN CALL MPI_TESTALL ( NRQSG2, IRQSG2(IOFF,2), DONE, & STATUS, IERR_MPI ) ELSE DONE = .TRUE. END IF IF ( DONE .AND. NRQSG2.GT.0 ) CALL & MPI_WAITALL ( NRQSG2, IRQSG2(IOFF,2), & STATUS, IERR_MPI ) IF ( DONE ) THEN BSTAT(IB0) = 0 END IF END IF END DO ! ! 2.f Last component, finish message passing, reset buffer control ! IF ( ISPLOC .EQ. NSPLOC ) THEN ! DO IB0=1, MPIBUF IF ( BSTAT(IB0) .EQ. 2 ) THEN IOFF = 1 + (BISPL(IB0)-1)*NRQSG2 IF ( NRQSG2 .GT. 0 ) CALL & MPI_WAITALL ( NRQSG2, IRQSG2(IOFF,2), & STATUS, IERR_MPI ) BSTAT(IB0) = 0 END IF END DO ! ISPLOC = 0 IBFLOC = 0 ! END IF ! ! 2.g Test output ! RETURN ! ! Formats ! !/ !/ End of W3SCAT ----------------------------------------------------- / !/ END SUBROUTINE W3SCAT !/ ------------------------------------------------------------------- / SUBROUTINE W3NMIN ( MAPSTA, FLAG0 ) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | FORTRAN 90 | !/ | Last update : 28-Dec-2004 | !/ +-----------------------------------+ !/ !/ 23-Feb-2001 : Origination. ( version 2.07 ) !/ 28-Dec-2004 : Multiple grid version. ( version 3.06 ) !/ ! 1. Purpose : ! ! Check minimum number of active sea points at given processor to ! evaluate the need for a MPI_BARRIER call. ! ! 2. Method : ! ! Evaluate mapsta. ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! MAPSTA I.A. I Status map for spatial grid. ! FLAG0 log. O Flag to identify 0 as minimum. ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3WAVE Subr. W3WAVEMD Actual wave model routine. ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/S Enable subroutine tracing. ! !/T Test output. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/ USE W3GDATMD, ONLY: NX, NY, NSEA, MAPSF USE W3ODATMD, ONLY: NDST, NAPROC USE W3PARALL, ONLY: INIT_GET_JSEA_ISPROC !/ IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ INTEGER, INTENT(IN) :: MAPSTA(NY*NX) LOGICAL, INTENT(OUT) :: FLAG0 !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ INTEGER :: NMIN, IPROC, NLOC, ISEA, IXY INTEGER :: JSEA, ISPROC !/ !/ ------------------------------------------------------------------- / !/ ! NMIN = NSEA ! DO IPROC=1, NAPROC NLOC = 0 DO ISEA=1, NSEA CALL INIT_GET_JSEA_ISPROC(ISEA, JSEA, ISPROC) IF (ISPROC .eq. IPROC) THEN IXY = MAPSF(ISEA,3) IF ( MAPSTA(IXY) .EQ. 1 ) NLOC = NLOC + 1 END IF END DO ! NMIN = MIN ( NMIN , NLOC ) END DO ! FLAG0 = NMIN .EQ. 0 ! RETURN ! ! Formats ! !/ !/ End of W3NMIN ----------------------------------------------------- / !/ END SUBROUTINE W3NMIN !/ !/ End of module W3WAVEMD -------------------------------------------- / !/ END MODULE W3WAVEMD