! ! SWAN/SWREAD file 2 of 2 ! ! Contents of this file: ! SPROUT: Reading and processing of the user output commands ! SWREPS: Reading and processing of the commands defining output points ! SWREOQ: Reading and processing of the output requests ! SIRAY : Searching the first point on a ray where the depth is DP ! SWNMPS ! SVARTP ! SWBOUN 40.00 ! BCFILE 40.00 ! BCWAMN 40.00 ! BCWW3N ! SWBCPT ! RETSTP 40.00 ! !*********************************************************************** ! * SUBROUTINE SPROUT (FOUND, BOTLEV, WATLEV) 40.31 ! * !*********************************************************************** ! USE OCPCOMM4 40.41 USE SWCOMM3 40.41 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 30.70: Nico Booij ! 30.72: IJsbrand Haagsma ! 30.81: Annette Kieftenburg ! 30.82: IJsbrand Haagsma ! 30.90: IJsbrand Haagsma (Equivalence version) ! 32.02: Roeland Ris & Cor van der Schelde (1D version) ! 34.01: Jeroen Adema ! 40.02: IJsbrand Haagsma ! 40.03: Nico Booij ! 40.31: Marcel Zijlema ! 40.41: Marcel Zijlema ! ! 1. Updates ! ! 100.04, Nov. 92: Filename of plotfile will be given by user ! 30.70, Nov. 97: Arguments BOTLEV and WATLEV added ! 32.02, Feb. 98: 1D version introduced ! 30.72, Feb. 98: Introduced generic names XCGRID, YCGRID and SPCSIG for SWAN ! 30.82, Apr. 98: Removed reference to commons KAART and KAR ! 30.90, Oct. 98: Introduced EQUIVALENCE POOL-arrays ! 30.81, Nov. 98: Replaced variable STATUS by IERR (because STATUS is a ! reserved word) ! 30.81, Jan. 99: Replaced variable FROM by FROM_ (because FROM is a ! reserved word) ! 34.01, Feb. 99: Introducing STPNOW ! 40.03, Sep. 00: inconsistency with manual corrected ! 40.02, Oct. 00: Initialisation of IERR ! 40.31, Nov. 03: removing POOL construction and HPGL functionality ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. Purpose ! ! Reading and processing of the user output commands ! ! 3. Method ! ! If the first characters of the last read command are equal to a ! given string (KEYWIS ('STRING')), the keywords and varia- ! bles of this command are further read and processed ! ! 4. Argument variables ! ! i BOTLEV: Bottom levels 30.70 ! i WATLEV: Water levels 30.70 ! REAL BOTLEV(*) 30.70 REAL WATLEV(*) 30.70 ! ! 6. Local variables ! ! 8. Subroutines used ! ! MSGERR ! SWREPS ! SWREOQ ! STPNOW ! LOGICAL STPNOW 34.01 ! ! 9. Subroutines calling ! ! SWREAD ! ! 10. Error messages ! ! --- ! ! 11. Remarks ! ! --- ! ! 12. Structure ! ! ---------------------------------------------------------------- ! Most of the source code will be clear with the ! aid of the user manual, the system documentation ! and the additional comments in the source code. ! ---------------------------------------------------------------- ! ! 13. Source text ! LOGICAL FOUND LOGICAL KEYWIS SAVE IENT DATA IENT/0/ 40.31 30.81 CALL STRACE (IENT,'SPROUT') ! FOUND = .FALSE. ! ! definition of output point sets ! CALL SWREPS ( FOUND, BOTLEV, WATLEV ) 40.31 IF (STPNOW()) RETURN 34.01 IF (FOUND) RETURN ! ! output requests ! CALL SWREOQ ( FOUND ) 40.31 30.90 IF (STPNOW()) RETURN 34.01 IF (FOUND) RETURN ! IF (KEYWIS('SIT') .OR. KEYWIS('PLA')) THEN CALL MSGERR(2,'Keyword SITES is no longer maintained') 40.31 GOTO 800 40.31 ENDIF ! IF (KEYWIS ('LIN')) THEN CALL MSGERR(2,'Keyword LINE is no longer maintained') 40.31 GOTO 800 40.31 ENDIF ! ------------------------------------------------------------------- ! ***** command name not found ***** RETURN ! 800 FOUND = .TRUE. RETURN ! * end of subroutine SPROUT * END !*********************************************************************** ! * SUBROUTINE SWREPS ( FOUND, BOTLEV, WATLEV ) 40.31 ! * !*********************************************************************** ! USE OCPCOMM1 40.41 USE OCPCOMM3 40.41 USE OCPCOMM4 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 USE SWCOMM4 40.41 USE OUTP_DATA 40.31 USE M_PARALL 40.31 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 30.70, 40.03, 40.13: Nico Booij ! 30.72: IJsbrand Haagsma ! 30.81: Annette Kieftenburg ! 30.82: IJsbrand Haagsma ! 32.02: Roeland Ris & Cor van der Schelde (1D version) ! 34.01: Jeroen Adema ! 40.30: Marcel Zijlema ! 40.31: Marcel Zijlema ! 40.41: Marcel Zijlema ! ! 1. Updates ! ! 30.72, Sept 97: Changed DO-block with one CONTINUE to DO-block with ! two CONTINUE's ! 30.70, Nov. 97: comm ISO, inquire pointer added to get correct value ! for IADRAY ! 30.70, Nov. 97: comm ISO, offset origin added in message concerning rays ! declaration INT SIRAY added ! 30.70, Nov. 97: arguments BOTLEV and WATLEV added ! 30.72, Feb. 98: Declaration of Argument variables updated ! 32.02, Feb. 98: 1D version introduced ! 30.72, Feb. 98: Introduced generic names XCGRID, YCGRID and SPCSIG for SWAN ! 30.82, Apr. 98: removed reference to commons KAART and KAR ! 30.81, Nov. 98: Replaced variable STATUS by IERR (because STATUS is a ! reserved word) ! 34.01, Feb. 99: Introducing STPNOW ! 40.01, Sep. 99: XASM and YASM replace fixed numbers ! 33.09, Sep. 00: modifications in view of spherical coordinates ! 40.03, Sep. 00: inconsistency with manual corrected ! 40.13, Sep. 01: nesting in curvilinear grid: division by 0 prevented ! 40.30, May 03: introduction distributed-memory approach using MPI ! 40.31, Dec. 03: removing POOL-mechanism ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. PURPOSE ! ! Reading and processing of the commands defining output points ! ! 4. Argument variables (updated 30.72) ! ! BOTLEV: input bottom levels 30.70 ! WATLEV: input water levels 30.70 ! REAL BOTLEV(*), WATLEV(*) 30.70 ! ! FOUND : output parameter indicating whether command ! being processed is found (value True) ! or not (False) ! LOGICAL FOUND ! ! Local variables REAL :: XPFR, YPFR, XLENFR, YLENFR 40.31 TYPE(OPSDAT), POINTER :: OPSTMP, ROPS 40.31 TYPE XYPT 40.31 REAL :: X, Y, XQ, YQ TYPE(XYPT), POINTER :: NEXTXY END TYPE XYPT TYPE(XYPT), TARGET :: FRST 40.31 TYPE(XYPT), POINTER :: CURR, TMP 40.31 ! 8. Subroutines used ! ! command reading routines ! (all Ocean Pack) LOGICAL :: STPNOW 34.01 LOGICAL :: EQREAL ! if True the two (real) arguments are equal 33.09 ! 9. Subroutines calling ! ! SPROUT ! ! 10. Error messages ! ! --- ! ! 13. Source text ! LOGICAL PP 30.72 INTEGER IERR, SIRAY 30.81 30.72 CHARACTER PSNAME*16, STYPE*1, PRNAME*16 40.31 30.21 LOGICAL KEYWIS, BOTDEP 30.70 SAVE IENT DATA IENT/0/ CALL STRACE (IENT,'SWREPS') ! ! -------------------------------------------------------------------- ! FRAME 'sname' [xpfr] [ypfr] [alpfr] [xlenfr] [ylenfr] & ! [mxfr] [myfr] ! -------------------------------------------------------------------- IF (KEYWIS ('FRA')) THEN ! IF (ONED) THEN 32.02 CALL MSGERR (2,' Illegal keyword (FRA) in combination with'// 32.02 & ' 1D-computation') 32.02 GOTO 800 32.02 ELSE 32.02 ! ver 30.20: names of input variables changed, order of data changed ALLOCATE(OPSTMP) 40.31 CALL INCSTR ('SNAME',PSNAME,'REQ',' ') IF (LENCST.GT.8) CALL MSGERR (2, 'SNAME is too long') OPSTMP%PSNAME = PSNAME 40.31 CALL READXY ('XPFR', 'YPFR', XPFR, YPFR, 'REQ', 0., 0.) 40.31 30.20 OPSTMP%OPR(1) = XPFR 40.31 OPSTMP%OPR(2) = YPFR 40.31 CALL INREAL('ALPFR',ALPK,'REQ',0.) 30.20 IF (KSPHER.GT.0 .AND. .NOT.EQREAL(ALPK,0.)) CALL MSGERR (2, & '[alpfr] must be 0 with spherical coordinates') 33.09 CALL INREAL('XLENFR', XLENFR,'REQ',0.) 40.31 30.20 CALL INREAL('YLENFR', YLENFR,'REQ',0.) 40.31 30.20 OPSTMP%OPR(3) = XLENFR 40.31 OPSTMP%OPR(4) = YLENFR 40.31 OPSTMP%OPR(5) = PI2 * (ALPK/360.-NINT(ALPK/360.)) ! ***** the user gives number of meshes along each side ***** ! ***** program uses the number of points ***** CALL ININTG ('MXFR',MXK,'STA',20) 30.20 CALL ININTG ('MYFR',MYK,'STA',20) 30.20 OPSTMP%PSTYPE = 'F' 40.31 OPSTMP%OPI(1) = MXK+1 40.31 OPSTMP%OPI(2) = MYK+1 40.31 ALLOCATE(OPSTMP%XP(0)) 40.31 ALLOCATE(OPSTMP%YP(0)) 40.31 NULLIFY(OPSTMP%NEXTOPS) 40.31 IF ( .NOT.LOPS ) THEN 40.31 FOPS = OPSTMP 40.31 COPS => FOPS 40.31 LOPS = .TRUE. 40.31 ELSE 40.31 COPS%NEXTOPS => OPSTMP 40.31 COPS => OPSTMP 40.31 END IF 40.31 GOTO 800 ENDIF 32.02 ENDIF ! ! ------------------------------------------------------------------ ! GROUP 'sname' SUBGRID [ix1] [ix2] [iy1] [iy2] ! ------------------------------------------------------------------ IF (KEYWIS('GROUP') .OR. KEYWIS ('SUBG')) THEN 970221 ! IF (ONED) THEN 32.02 CALL MSGERR (2,' Illegal keyword (GROUP) in combination'// 32.02 & ' with 1D-computation') 32.02 GOTO 800 32.02 ELSE 32.02 ! mod 970221: GROUP is introduced as a new command instead of ! an option SUBG within the Frame command ALLOCATE(OPSTMP) 40.31 CALL INCSTR ('SNAME',PSNAME,'REQ',' ') IF (LENCST.GT.8) CALL MSGERR (2, 'SNAME is too long') OPSTMP%PSNAME = PSNAME 40.31 CALL INKEYW ('STA', ' ') CALL IGNORE ('SUBG') 970221 CALL ININTG ('IX1', IX1, 'REQ', 0) CALL ININTG ('IX2', IX2, 'REQ', 0) CALL ININTG ('IY1', IY1, 'REQ', 0) CALL ININTG ('IY2', IY2, 'REQ', 0) IF (IX1 .LT. 0 .OR. IX2 .GT. MXCGL-1 .OR. IX1 .GT. IX2 .OR. & IY1 .LT. 0 .OR. IY2 .GT. MYCGL-1 .OR. IY1 .GT. IY2) THEN CALL MSGERR (3, 'Check corners of GROUP (SUBGRID) command') CALL MSGERR (3, ' .........the values should be.........') CALL MSGERR (3, 'ix1 FOPS 40.31 LOPS = .TRUE. 40.31 ELSE 40.31 COPS%NEXTOPS => OPSTMP 40.31 COPS => OPSTMP 40.31 END IF 40.31 GOTO 800 ENDIF 32.01 ENDIF ! ! ------------------------------------------------------------------ ! CURVE 'sname' [xp1] [yp1] < [int] [xp] [yp] > ! ------------------------------------------------------------------ ! command CURVE, MIP is number of point on a curve IF (KEYWIS ('CURV')) THEN ALLOCATE(OPSTMP) 40.31 CALL INCSTR('SNAME',PSNAME,'REQ',' ') IF (LENCST.GT.8) CALL MSGERR (2, 'SNAME is too long') OPSTMP%PSNAME = PSNAME 40.31 OPSTMP%PSTYPE = 'C' 40.31 MIP = 0 OPSTMP%MIP = MIP 40.31 ! ***** first point of a curve ***** 30 CALL NWLINE IF (STPNOW()) RETURN 34.01 CALL READXY ('XP1', 'YP1', XP, YP, 'REQ', 0., 0.) FRST%X = XP 40.31 FRST%Y = YP 40.31 NULLIFY(FRST%NEXTXY) 40.31 CURR => FRST 40.31 MIP = 1 ! ***** interval and next corner point ***** 33 CALL ININTG ('INT',INTV,'REP',-1) IF (INTV .NE. -1) THEN IF (INTV .LE. 0) THEN CALL MSGERR (2,'INT is negative or zero') INTV = 1 ENDIF XP1 = XP YP1 = YP CALL READXY ('XP', 'YP', XP, YP, 'REQ', 0., 0.) 40.03 IF (ITEST .GE. 200 .OR. INTES .GE. 20) THEN 30.21 WRITE(PRINTF, 31) PSNAME 31 FORMAT ('COORDINATES OF OUTPUT POINTS FOR CURVE : ', A) ENDIF DO 36 JJ=1,INTV MIP = MIP+1 ALLOCATE(TMP) 40.31 TMP%X = XP1+REAL(JJ)*(XP-XP1)/REAL(INTV) 40.31 TMP%Y = YP1+REAL(JJ)*(YP-YP1)/REAL(INTV) 40.31 NULLIFY(TMP%NEXTXY) 40.31 CURR%NEXTXY => TMP 40.31 CURR => TMP 40.31 36 CONTINUE GOTO 33 ENDIF ALLOCATE(OPSTMP%XP(MIP)) 40.31 ALLOCATE(OPSTMP%YP(MIP)) 40.31 CURR => FRST 40.31 DO JJ = 1, MIP 40.31 OPSTMP%XP(JJ) = CURR%X 40.31 OPSTMP%YP(JJ) = CURR%Y 40.31 IF (ITEST .GE. 200 .OR. INTES .GE. 50) THEN 40.31 WRITE(PRINTF,32) JJ, CURR%X, CURR%Y 40.31 32 FORMAT(' POINT(',I4,')',' (IX,IY) -> ',2F10.2) 40.31 ENDIF 40.31 CURR => CURR%NEXTXY 40.31 END DO 40.31 DEALLOCATE(TMP) 40.31 ! ***** store number of points of the curve ***** OPSTMP%MIP = MIP 40.31 IF (MIP .EQ. 0) CALL MSGERR(1,'No output points found') NULLIFY(OPSTMP%NEXTOPS) 40.31 IF ( .NOT.LOPS ) THEN 40.31 FOPS = OPSTMP 40.31 COPS => FOPS 40.31 LOPS = .TRUE. 40.31 ELSE 40.31 COPS%NEXTOPS => OPSTMP 40.31 COPS => OPSTMP 40.31 END IF 40.31 GOTO 800 ENDIF ! ! ------------------------------------------------------------------ ! POINTS 'sname' < [xp] [yp] > / FILE 'fname' ! ------------------------------------------------------------------ ! POINTS set of individual output points IF (KEYWIS ('POIN')) THEN ALLOCATE(OPSTMP) 40.31 CALL INCSTR('SNAME',PSNAME,'REQ',' ') IF (LENCST.GT.8) CALL MSGERR (2, 'SNAME is too long') OPSTMP%PSNAME = PSNAME 40.31 OPSTMP%PSTYPE = 'P' 40.31 MIP = 0 OPSTMP%MIP = MIP 40.31 CALL INKEYW ('STA', ' ') IF (KEYWIS('FILE')) THEN IOSTAT = 0 NDS = 0 PP = .TRUE. CALL INCSTR ('FNAME', FILENM, 'REQ', ' ') CALL FOR (NDS, FILENM, 'OF', IOSTAT) 10.31 IF (STPNOW()) RETURN 34.01 ELSE PP = .FALSE. ENDIF FRST%X = 0. 40.31 FRST%Y = 0. 40.31 NULLIFY(FRST%NEXTXY) 40.31 CURR => FRST 40.31 DO IF (PP) THEN IERR = 0 10.18 CALL REFIXY (NDS, XP, YP, IERR) 10.18 IF (IERR.EQ.-1) GOTO 47 IF (IERR.EQ.-2) THEN CALL MSGERR (2, 'Error reading point coord. from file') GOTO 800 ENDIF ELSE CALL READXY ('XP', 'YP', XP, YP, 'REP', -1.E10, -1.E10) IF (XP .LT. -0.9E10) GOTO 47 ENDIF MIP = MIP+1 ALLOCATE(TMP) 40.31 TMP%X = XP 40.31 TMP%Y = YP 40.31 NULLIFY(TMP%NEXTXY) 40.31 CURR%NEXTXY => TMP 40.31 CURR => TMP 40.31 ENDDO 47 ALLOCATE(OPSTMP%XP(MIP)) 40.31 ALLOCATE(OPSTMP%YP(MIP)) 40.31 CURR => FRST%NEXTXY 40.31 DO JJ = 1, MIP 40.31 OPSTMP%XP(JJ) = CURR%X 40.31 OPSTMP%YP(JJ) = CURR%Y 40.31 CURR => CURR%NEXTXY 40.31 END DO 40.31 DEALLOCATE(TMP) 40.31 ! ***** store number of output points ***** OPSTMP%MIP = MIP 40.31 IF (MIP .EQ. 0) CALL MSGERR (2, 'No output points found') 10.32 NULLIFY(OPSTMP%NEXTOPS) 40.31 IF ( .NOT.LOPS ) THEN 40.31 FOPS = OPSTMP 40.31 COPS => FOPS 40.31 LOPS = .TRUE. 40.31 ELSE 40.31 COPS%NEXTOPS => OPSTMP 40.31 COPS => OPSTMP 40.31 END IF 40.31 GOTO 800 ENDIF ! ! ------------------------------------------------------------------- ! RAY 'rname' [xp1] [yp1] [xq1] [yq1] & ! < [int] [xp] [yp] [xq] [yq] > ! ------------------------------------------------------------------- ! RAY set of rays IF (KEYWIS ('RAY')) THEN ! IF (ONED) THEN 32.02 CALL MSGERR (2,' Illegal keyword (RAY) in combination'// 32.02 & ' with 1D-computation') 32.02 GOTO 800 32.02 ELSE 32.02 ALLOCATE(OPSTMP) 40.31 CALL INCSTR('RNAME',PSNAME,'REQ',' ') IF (LENCST.GT.8) CALL MSGERR (2, 'RNAME is too long') OPSTMP%PSNAME = PSNAME 40.31 OPSTMP%PSTYPE = 'R' 40.31 MIP = 1 OPSTMP%MIP = MIP 40.31 ! first ray CALL NWLINE IF (STPNOW()) RETURN 34.01 CALL READXY ('XP1', 'YP1', XP, YP, 'REQ', 0., 0.) CALL READXY ('XQ1', 'YQ1', XQ, YQ, 'REQ', 0., 0.) FRST%X = XP 40.31 FRST%Y = YP 40.31 FRST%XQ = XQ 40.31 FRST%YQ = YQ 40.31 NULLIFY(FRST%NEXTXY) 40.31 CURR => FRST 40.31 ! following rays 110 CALL ININTG ('INT',INTD,'REP',-1) IF (INTD .NE. -1) THEN IF (INTD .LE. 0) THEN CALL MSGERR(2, 'INT negative or zero') INTD = 1 ENDIF XP1 = XP YP1 = YP XQ1 = XQ YQ1 = YQ CALL READXY ('XP', 'YP', XP, YP, 'REQ', 0., 0.) CALL READXY ('XQ', 'YQ', XQ, YQ, 'REQ', 0., 0.) DO 115 JJ=1,INTD MIP = MIP+1 ALLOCATE(TMP) 40.31 TMP%X = XP1 + REAL(JJ)*(XP-XP1)/REAL(INTD) 40.31 TMP%Y = YP1 + REAL(JJ)*(YP-YP1)/REAL(INTD) 40.31 TMP%XQ = XQ1 + REAL(JJ)*(XQ-XQ1)/REAL(INTD) 40.31 TMP%YQ = YQ1 + REAL(JJ)*(YQ-YQ1)/REAL(INTD) 40.31 NULLIFY(TMP%NEXTXY) 40.31 CURR%NEXTXY => TMP 40.31 CURR => TMP 40.31 115 CONTINUE GOTO 110 ENDIF ALLOCATE(OPSTMP%XP(MIP)) 40.31 ALLOCATE(OPSTMP%YP(MIP)) 40.31 ALLOCATE(OPSTMP%XQ(MIP)) 40.31 ALLOCATE(OPSTMP%YQ(MIP)) 40.31 CURR => FRST 40.31 DO JJ = 1, MIP 40.31 OPSTMP%XP(JJ) = CURR%X 40.31 OPSTMP%YP(JJ) = CURR%Y 40.31 OPSTMP%XQ(JJ) = CURR%XQ 40.31 OPSTMP%YQ(JJ) = CURR%YQ 40.31 CURR => CURR%NEXTXY 40.31 END DO 40.31 DEALLOCATE(TMP) 40.31 ! ! ***** termination ***** OPSTMP%MIP = MIP 40.31 IF (MIP .EQ. 1) CALL MSGERR (1,'Only one ray is defined') NULLIFY(OPSTMP%NEXTOPS) 40.31 IF ( .NOT.LOPS ) THEN 40.31 FOPS = OPSTMP 40.31 COPS => FOPS 40.31 LOPS = .TRUE. 40.31 ELSE 40.31 COPS%NEXTOPS => OPSTMP 40.31 COPS => OPSTMP 40.31 END IF 40.31 GOTO 800 ENDIF 32.02 ENDIF ! ! ------------------------------------------------------------------- ! ISOline 'sname' 'rname' DEPTH / BOTTOM [dep] ! ------------------------------------------------------------------- ! ISO depth or bottom level contour line IF (KEYWIS ('ISO')) THEN ! IF (ONED) THEN 32.02 CALL MSGERR (2,' Illegal keyword (ISO) in combination'// 32.02 & ' with 1D-computation') 32.02 GOTO 800 32.02 ELSE 32.02 ALLOCATE(OPSTMP) 40.31 CALL INCSTR ('SNAME',PSNAME,'REQ',' ') IF (LENCST.GT.8) CALL MSGERR (2, 'SNAME is too long') OPSTMP%PSNAME = PSNAME 40.31 CALL INCSTR ('RNAME', PRNAME, 'REQ', ' ') IF (LENCST.GT.8) CALL MSGERR (2, 'RNAME is too long') CALL INKEYW ('STA', 'DEP') IF (KEYWIS ('BOT')) THEN BOTDEP = .TRUE. 30.70 ELSE CALL IGNORE ('DEP') BOTDEP = .FALSE. 30.70 IF (DYNDEP) CALL MSGERR (2,'depths will vary with time') 40.00 ENDIF CALL INREAL ('DEP',DP,'REP',-1.E10) ROPS => FOPS 40.31 DO 40.31 IF (ROPS%PSNAME.EQ.PRNAME) EXIT 40.31 IF (.NOT.ASSOCIATED(ROPS%NEXTOPS)) THEN 40.31 CALL MSGERR(2,'Set of rays not defined') 40.31 GOTO 800 40.31 END IF 40.31 ROPS => ROPS%NEXTOPS 40.31 END DO 40.31 STYPE = ROPS%PSTYPE 40.31 IF (STYPE .NE. 'R') THEN 40.31 CALL MSGERR 40.31 & (2,'Ray name set assigned to set of output locations') 40.31 GOTO 800 40.31 END IF 40.31 MIPR = ROPS%MIP 40.31 OPSTMP%PSTYPE = 'C' 40.31 MIP = 0 OPSTMP%MIP = MIP 40.31 FRST%X = 0. 40.31 FRST%Y = 0. 40.31 NULLIFY(FRST%NEXTXY) 40.31 CURR => FRST 40.31 DO 125 IK=1,MIPR XP = ROPS%XP(IK) 40.31 YP = ROPS%YP(IK) 40.31 XQ = ROPS%XQ(IK) 40.31 YQ = ROPS%YQ(IK) 40.31 II = SIRAY (DP, XP, YP, XQ, YQ, XX, YY, BOTDEP, 30.70 & BOTLEV, WATLEV) 30.70 IF (II.EQ.0) THEN WRITE (PRINTF, 6120) DP, XP+XOFFS, YP+YOFFS, 30.70 & XQ+XOFFS, YQ+YOFFS 30.70 6120 FORMAT(' No point with depth ',F5.2, & ' is found in ray :',4F10.2) ELSE MIP = MIP+1 ALLOCATE(TMP) 40.31 TMP%X = XX 40.31 TMP%Y = YY 40.31 NULLIFY(TMP%NEXTXY) 40.31 CURR%NEXTXY => TMP 40.31 CURR => TMP 40.31 ENDIF 125 CONTINUE ALLOCATE(OPSTMP%XP(MIP)) 40.31 ALLOCATE(OPSTMP%YP(MIP)) 40.31 CURR => FRST%NEXTXY 40.31 DO IK = 1, MIP 40.31 OPSTMP%XP(IK) = CURR%X 40.31 OPSTMP%YP(IK) = CURR%Y 40.31 CURR => CURR%NEXTXY 40.31 END DO 40.31 DEALLOCATE(TMP) 40.31 IF (MIP.EQ.0) CALL MSGERR & (2, 'No points with valid depth found') ! ***** store number of points of the curve ***** OPSTMP%MIP = MIP 40.31 NULLIFY(OPSTMP%NEXTOPS) 40.31 IF ( .NOT.LOPS ) THEN 40.31 FOPS = OPSTMP 40.31 COPS => FOPS 40.31 LOPS = .TRUE. 40.31 ELSE 40.31 COPS%NEXTOPS => OPSTMP 40.31 COPS => OPSTMP 40.31 END IF 40.31 GOTO 800 ENDIF 32.02 ENDIF ! ! ------------------------------------------------------------------- ! NGRID 'sname' [xpn] [ypn] [alpn] [xlenn] [ylenn] [mxn] [myn] ! ------------------------------------------------------------------- ! NGRID computational grid for nesting 20.63 IF (KEYWIS ('NGR')) THEN ! ! =============================================================== ! ! NGRid 'sname' [xpn] [ypn] [alpn] [xlenn] [ylenn] [mxn] [myn] 40.00 ! ! ================================================================ ! IF (ONED) THEN 32.02 CALL MSGERR (2,' Illegal keyword (NGR) in combination'// 32.02 & ' with 1D-computation') 32.02 GOTO 800 32.02 ELSE 32.02 ! ver 30.20: names changed, order changed ALLOCATE(OPSTMP) 40.31 CALL INCSTR('SNAME',PSNAME,'REQ',' ') IF (LENCST.GT.8) CALL MSGERR (2, 'SNAME is too long') OPSTMP%PSNAME = PSNAME 40.31 OPSTMP%PSTYPE = 'N' 40.31 CALL READXY ('XPN', 'YPN', XPCN, YPCN , 'REQ', 0., 0.) 30.20 CALL INREAL('ALPN',ALPCN,'REQ',0.) 30.20 CALL INREAL('XLENN',XNLEN,'REQ',0.) 30.20 CALL INREAL('YLENN',YNLEN,'REQ',0.) 30.20 ALTNP = ALPCN / 360. ANG = PI2 * (ALTNP - NINT(ALTNP)) ! estimate step size for output 40.00 IF (OPTG.EQ.1) THEN 40.13 COSA2 = (COS(ANG-ALPC))**2 40.00 SINA2 = (SIN(ANG-ALPC))**2 40.00 DXN = DX*COSA2 + DY*SINA2 40.00 DYN = DX*SINA2 + DY*COSA2 40.00 ELSE 40.13 ! curvilinear grid, DXN and DYN are average step size 40.13 DXN = (XCLEN+YCLEN)/REAL(MXCGL+MYCGL) 40.31 40.13 DYN = DXN 40.13 ENDIF 40.13 CALL ININTG ('MXN',MXN,'STA',MAX(1,NINT(XNLEN/DXN))) 40.00 CALL ININTG ('MYN',MYN,'STA',MAX(1,NINT(YNLEN/DYN))) 40.00 MIP = 0 ALLOCATE(OPSTMP%XP(2*(MXN+MYN))) 40.31 ALLOCATE(OPSTMP%YP(2*(MXN+MYN))) 40.31 XF=XPCN YF=YPCN ! ***** start to calculate the positions ******** ! ***** of the boundary point in the four sides ******** DO 50 I=1,4 INTE=MXN ALON=XNLEN ANGLE=ANG+PI2*(90.*REAL(I-1))/360. IF (I .EQ. 2 .OR. I .EQ. 4) THEN INTE=MYN ALON=YNLEN ENDIF XI=XF YI=YF XF=XI+ALON*COS(ANGLE) YF=YI+ALON*SIN(ANGLE) DO 55 KK=1,INTE MIP=MIP+1 OPSTMP%XP(MIP)=XI+REAL(KK)*(XF-XI)/REAL(INTE) OPSTMP%YP(MIP)=YI+REAL(KK)*(YF-YI)/REAL(INTE) 55 CONTINUE 50 CONTINUE ! ***** store number of points of the curve ***** OPSTMP%MIP = MIP 40.31 OPSTMP%OPR(1) = XNLEN 40.31 OPSTMP%OPR(2) = YNLEN 40.31 OPSTMP%OPR(3) = XPCN 40.31 OPSTMP%OPR(4) = YPCN 40.31 OPSTMP%OPR(5) = PI2 * (ALPCN/360.-NINT(ALPCN/360.)) 40.31 OPSTMP%OPI(1) = MXN 40.31 OPSTMP%OPI(2) = MYN 40.31 IF (MIP .EQ. 0) CALL MSGERR(1,'No output points found') NULLIFY(OPSTMP%NEXTOPS) 40.31 IF ( .NOT.LOPS ) THEN 40.31 FOPS = OPSTMP 40.31 COPS => FOPS 40.31 LOPS = .TRUE. 40.31 ELSE 40.31 COPS%NEXTOPS => OPSTMP 40.31 COPS => OPSTMP 40.31 END IF 40.31 GOTO 800 ENDIF 32.02 ENDIF ! --------------------------------------------------------- ! command not found: RETURN 800 FOUND = .TRUE. RETURN !* end of subroutine SWREPS ** END !*********************************************************************** ! * SUBROUTINE SWREOQ ( FOUND ) 40.31 ! * !*********************************************************************** ! USE OCPCOMM1 40.41 USE OCPCOMM2 40.41 USE OCPCOMM3 40.41 USE OCPCOMM4 40.41 USE SWCOMM1 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 USE SWCOMM4 40.41 USE OUTP_DATA 40.13 USE M_PARALL 40.31 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 30.72: IJsbrand Haagsma ! 30.81: Annette Kieftenburg ! 30.82: IJsbrand Haagsma ! 32.02: Roeland Ris & Cor van der Schelde (1D version) ! 34.01: Jeroen Adema ! 40.03, 40.13: Nico Booij ! 40.14: Annette Kieftenburg ! 40.30: Marcel Zijlema ! 40.31: Marcel Zijlema ! 40.41: Marcel Zijlema ! ! 1. Updates ! ! 30.50 : option COORD added in command PLOT ! option STAR added in command PLOT ! 32.02, Feb. 98: 1D version introduced ! 30.72, Feb. 98: Introduced generic names XCGRID, YCGRID and SPCSIG for SWAN ! 30.82, Apr. 98: removed reference to commons KAART and KAR ! 30.81, Nov. 98: Replaced variable STATUS by IERR (because STATUS is a ! reserved word) ! 30.81, Jan. 99: Replaced variable TO by TO_ (because TO is a reserved ! word) ! 34.01, Feb. 99: Introducing STPNOW ! 40.03, Nov. 99: in case SPEC2D the value of MXOUTAR is increased by ! 6*MIP ! 40.03, Mar. 00: NQUA increased in case of Isoline plot ! Sep. 00: inconsistency with manual corrected ! 40.13, Mar. 01: option BLOCKed added in plot of problem points ! Aug. 01: array for NESTOUT request extended to 20 (in view of SETUP) ! 40.13, Oct. 01: filenames are stored in array OUTP_FILES ! not any more in array containing output request parameters ! 40.14, Dec. 01: format for setup corrected. ! 40.30, Mar. 03: introduction distributed-memory approach using MPI ! 40.31, Nov. 03: removing POOL construction and HPGL funcationality ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. Purpose ! ! Reading and processing of the output requests ! ! 3. Method ! ! --- ! ! 4. Argument variables ! ! FOUND : output parameter indicating whether command ! being processed is found (value True) ! or not (False) ! LOGICAL FOUND 30.72 ! ! 8. Subroutines used ! ! command reading routines ! (all Ocean Pack) ! LOGICAL STPNOW 34.01 ! ! 9. Subroutines calling ! ! SPROUT ! ! 10. Error messages ! ! --- ! ! 13. Source text ! INTEGER IERR 40.31 CHARACTER PSNAME *16, STYPE *1, RTYPE *4 40.31 LOGICAL KEYWIS 40.31 TYPE(ORQDAT), POINTER :: ORQTMP 40.31 TYPE(ORQDAT), SAVE, POINTER :: CORQ 40.31 LOGICAL, SAVE :: LORQ = .FALSE. 40.31 TYPE AUXT 40.31 INTEGER :: I REAL :: R TYPE(AUXT), POINTER :: NEXTI END TYPE AUXT TYPE(AUXT), TARGET :: FRST 40.31 TYPE(AUXT), POINTER :: CURR, TMP 40.31 SAVE IENT DATA IENT/0/ CALL STRACE (IENT,'SWREOQ') ! ! -------------------------------------------------------------------------- ! BLOCK 'sname' HEADER / NOHEADER 'fname' (LAY-OUT [idla]) & ! < DSPR/HSIGN/DIR/PDIR/TDIR/TM01/RTM01/RTP/TM02/FSPR/DEPTH/VEL/ & ! FRCOEFF/WIND/DISSIP/QB/TRANSP/FORCE/UBOT/URMS/WLEN/STEEPNESS/ & ! DHSIGN/DRTM01/LEAK/XP/YP/DIST/SETUP/TMM10/RTMM10/TMBOT/ & ! WATLEV/BOTLEV/TPS/DISBOT/DISSURF/DISWCAP > ([unit]) & ! (OUTPUT [tbegblk] [deltblk] SEC/MIN/HR/DAY) ! -------------------------------------------------------------------------- ! BLO block type output IF (KEYWIS ('BLO')) THEN ! IF (ONED) THEN 32.02 CALL MSGERR (2,' Illegal keyword (BLO) in combination'// 32.02 & ' with 1D-computation') 32.02 GOTO 800 32.02 ELSE 32.02 CALL SWNMPS (PSNAME, STYPE, MIP, IERR) 40.31 IF (IERR.NE.0) GOTO 800 IF (STYPE .NE. 'F' .AND. STYPE .NE. 'H') THEN 30.21 CALL MSGERR(2,'Set of output locations is not type frame') GOTO 800 ENDIF ! ! output frame exists ! ALLOCATE(ORQTMP) 40.31 NREOQ = NREOQ + 1 40.31 IF (NREOQ.GT.MAX_OUTP_REQ) CALL MSGERR (2, 40.31 40.13 & 'too many output requests') 40.13 IDLAO = 1 ! CALL INKEYW ('REQ',' ') IF (KEYWIS('NOHEAD') .OR. KEYWIS ('FIL')) THEN 30.20 CALL INCSTR ('FNAME', FILENM, 'REQ', ' ') DFAC = 1. CALL INKEYW ('STA', ' ') IF (KEYWIS('LONG')) THEN ! option disabled 40.13 CALL MSGERR (2, 'option LONG disabled; use OUTP OPT') 40.13 ENDIF 40.13 RTYPE = 'BLKD' ELSE CALL IGNORE ('HEAD') 30.20 CALL IGNORE ('PAP') DFAC = -1. CALL INCSTR ('FNAME', FILENM, 'STA', ' ') 30.20 RTYPE = 'BLKP' IF ( INDEX( FILENM, '.MAT' ).NE.0 .OR. 40.41 40.30 & INDEX (FILENM, '.mat' ).NE.0 ) THEN 40.41 40.30 CALL MSGERR(4,'No header allowed for Matlab files') 40.30 RETURN 40.30 END IF 40.30 END IF IF (FILENM .EQ. ' ') THEN 24/FEB NREF = PRINTF ELSE NREF = 0 ! --- append node number to FILENM in case of 40.30 ! parallel computing 40.30 IF ( PARLL ) THEN 40.30 ILPOS = INDEX ( FILENM, ' ' )-1 40.30 WRITE(FILENM(ILPOS+1:ILPOS+4),33) INODE 40.30 33 FORMAT('-',I3.3) 40.30 END IF 40.30 ENDIF CALL INKEYW ('STA', ' ') IF (KEYWIS('LAY')) THEN CALL ININTG ('IDLA', IDLAO, 'REQ', 0) CALL INKEYW ('REQ', ' ') IF (IDLAO.NE.1 .AND. IDLAO.NE.3 .AND. IDLAO.NE.4) & CALL MSGERR (2, 'Illegal value for IDLA') ENDIF ORQTMP%OQR(1) = -1. 40.31 ORQTMP%OQR(2) = -1. 40.31 ORQTMP%RQTYPE = RTYPE 40.31 ORQTMP%PSNAME = PSNAME 40.31 ORQTMP%OQI(1) = NREF 40.31 ORQTMP%OQI(2) = NREOQ 40.31 NVAR = 0 40.31 ORQTMP%OQI(3) = NVAR 40.31 ORQTMP%OQI(4) = IDLAO 40.31 OUTP_FILES(NREOQ) = FILENM 40.13 ! ! read types of output quantities ! FRST%I = 0 40.31 FRST%R = 0. 40.31 NULLIFY(FRST%NEXTI) 40.31 CURR => FRST 40.31 70 CALL SVARTP (IVTYPE) IF (IVTYPE .EQ. 98) GOTO 91 30.00 IF (IVTYPE .NE. 99) THEN CALL INREAL ('UNIT', DFAC, 'STA', -1.) IF (OVSVTY(IVTYPE).EQ.5) THEN CALL MSGERR (2, & 'Type of output not allowed for this quantity') WRITE (PRINTF, *) ' -> ', OVSNAM(IVTYPE) ELSE IF (IVTYPE .GT. 0) THEN NVAR = NVAR+1 ALLOCATE(TMP) 40.31 TMP%I = IVTYPE 40.31 TMP%R = DFAC 40.31 NULLIFY(TMP%NEXTI) 40.31 CURR%NEXTI => TMP 40.31 CURR => TMP 40.31 IF (IVTYPE.EQ.6) IUBOTR = 1 ENDIF GOTO 70 ENDIF 91 IF (NVAR.GT.0) THEN 40.31 ALLOCATE(ORQTMP%IVTYP(NVAR)) 40.31 ALLOCATE(ORQTMP%FAC(NVAR)) 40.31 CURR => FRST%NEXTI 40.31 DO JJ = 1, NVAR 40.31 ORQTMP%IVTYP(JJ) = CURR%I 40.31 ORQTMP%FAC (JJ) = CURR%R 40.31 CURR => CURR%NEXTI 40.31 END DO 40.31 DEALLOCATE(TMP) 40.31 END IF 40.31 ! IF (IVTYPE .EQ. 98) THEN 30.00 IF (NSTATM.EQ.0) CALL MSGERR (3, & 'time information not allowed in stationary mode') NSTATM = 1 CALL INCTIM (ITMOPT, 'TBEG', ORQTMP%OQR(1), 'REQ', 0.) 40.31 30.00 CALL ININTV ('DELT', ORQTMP%OQR(2), 'REQ', 0.) 40.31 30.00 ENDIF ! ORQTMP%OQI(3) = NVAR 40.31 NULLIFY(ORQTMP%NEXTORQ) 40.31 IF ( .NOT.LORQ ) THEN 40.31 FORQ = ORQTMP 40.31 CORQ => FORQ 40.31 LORQ = .TRUE. 40.31 ELSE 40.31 CORQ%NEXTORQ => ORQTMP 40.31 CORQ => ORQTMP 40.31 END IF 40.31 GOTO 800 ENDIF 32.02 ENDIF ! -------------------------------------------------------------------------- ! TABLE 'sname' HEADER / NOHEADER / INDEXED 'fname' & ! < DSPR/HSIGN/DIR/PDIR/TDIR/TM01/RTM01/RTP/TM02/FSPR/DEPTH/VEL/ & ! FRCOEFF/WIND/DISSIP/QB/TRANSP/FORCE/UBOT/URMS/WLEN/STEEPNESS/ & ! DHSIGN/DRTM01/LEAK/XP/YP/DIST/SETUP/TMM10/RTMM10/TMBOT/ & ! WATLEV/BOTLEV/TPS/DISBOT/DISSURF/DISWCAP > ([unit]) & ! (OUTPUT [tbegtbl] [delttbl] SEC/MIN/HR/DAY) ! -------------------------------------------------------------------------- ! TABLE output in the form of a table IF (KEYWIS ('TAB')) THEN CALL SWNMPS (PSNAME, STYPE, MIP, IERR) 40.31 IF (IERR.NE.0) GOTO 800 ! ! output points exist ! ALLOCATE(ORQTMP) 40.31 NREOQ = NREOQ + 1 40.31 IF (NREOQ.GT.MAX_OUTP_REQ) CALL MSGERR (2, 40.31 40.13 & 'too many output requests') 40.13 ! CALL INKEYW ('STA','HEAD') 20.67 IF (KEYWIS('NOHEAD') .OR. KEYWIS ('FIL')) THEN 20.67 RTYPE = 'TABD' ELSE IF (KEYWIS ('IND')) THEN 30.50 RTYPE = 'TABI' ELSE IF (KEYWIS ('SWAN')) THEN 40.00 RTYPE = 'TABS' ELSE IF (KEYWIS ('STAB')) THEN 40.00 RTYPE = 'TABT' ELSE CALL IGNORE ('HEAD') 20.67 CALL IGNORE ('PAP') RTYPE = 'TABP' END IF ORQTMP%OQR(1) = -1. 40.31 ORQTMP%OQR(2) = -1. 40.31 ORQTMP%RQTYPE = RTYPE 40.31 ! unit reference number NREF is 0, will be determined in output module CALL INCSTR ('FNAME', FILENM, 'STA', ' ') IF (FILENM .NE. ' ') THEN NREF = 0 ! --- append node number to FILENM in case of 40.30 ! parallel computing 40.30 IF ( PARLL ) THEN 40.30 ILPOS = INDEX ( FILENM, ' ' )-1 40.30 WRITE(FILENM(ILPOS+1:ILPOS+4),33) INODE 40.30 END IF 40.30 ELSE NREF = PRINTF ENDIF ORQTMP%PSNAME = PSNAME 40.31 ORQTMP%OQI(1) = NREF 40.31 ORQTMP%OQI(2) = NREOQ 40.31 OUTP_FILES(NREOQ) = FILENM 40.31 40.13 ! NVAR = 0 ORQTMP%OQI(3) = NVAR 40.31 ! read types of variables to be printed in the table FRST%I = 0 40.31 NULLIFY(FRST%NEXTI) 40.31 CURR => FRST 40.31 80 CALL SVARTP (IVTYPE) IF (IVTYPE .EQ. 98) GOTO 90 30.00 IF (IVTYPE .NE. 99) THEN IF (OVSVTY(IVTYPE).EQ.5) THEN CALL MSGERR (2, & 'Type of output not allowed for this quantity') WRITE (PRINTF, *) ' -> ', OVSNAM(IVTYPE) ELSE IF (IVTYPE .GT. 0) THEN NVAR = NVAR+1 ALLOCATE(TMP) 40.31 TMP%I = IVTYPE 40.31 NULLIFY(TMP%NEXTI) 40.31 CURR%NEXTI => TMP 40.31 CURR => TMP 40.31 IF (IVTYPE.EQ.18) IUBOTR = 1 CALL INKEYW ('STA', ' ') 40.00 IF (KEYWIS('UNIT')) THEN CALL MSGERR (1, 'UNIT is ignored in this version') 40.00 ENDIF ENDIF GOTO 80 ENDIF 90 IF (NVAR.GT.0) THEN 40.31 ALLOCATE(ORQTMP%IVTYP(NVAR)) 40.31 CURR => FRST%NEXTI 40.31 DO JJ = 1, NVAR 40.31 ORQTMP%IVTYP(JJ) = CURR%I 40.31 CURR => CURR%NEXTI 40.31 END DO 40.31 DEALLOCATE(TMP) 40.31 END IF 40.31 ALLOCATE(ORQTMP%FAC(0)) 40.31 IF (IVTYPE .EQ. 98) THEN 30.00 IF (NSTATM.EQ.0) CALL MSGERR (3, & 'time information not allowed in stationary mode') NSTATM = 1 CALL INCTIM (ITMOPT, 'TBEG', ORQTMP%OQR(1), 'REQ', 0.) 40.31 30.00 CALL ININTV ('DELT', ORQTMP%OQR(2), 'REQ', 0.) 40.31 30.00 ENDIF ORQTMP%OQI(3) = NVAR 40.31 NULLIFY(ORQTMP%NEXTORQ) 40.31 IF ( .NOT.LORQ ) THEN 40.31 FORQ = ORQTMP 40.31 CORQ => FORQ 40.31 LORQ = .TRUE. 40.31 ELSE 40.31 CORQ%NEXTORQ => ORQTMP 40.31 CORQ => ORQTMP 40.31 END IF 40.31 GOTO 800 ENDIF ! ! PLOT plot iso lines and/or vector fields ! ! -------------------------------------------------------------------------- IF (KEYWIS ('PLO')) THEN CALL MSGERR(2,'Keyword PLO... is no longer maintained') 40.31 GOTO 800 ENDIF ! ! -------------------------------------------------------------------------- ! SPECout 'sname' SPEC1D/SPEC2D ABS/REL 'fname' & ! (OUTPUT [tbegspc] [deltspc] SEC/MIN/HR/DAY) ! -------------------------------------------------------------------------- ! SPEC output of spectra IF (KEYWIS ('SPEC')) THEN CALL SWNMPS (PSNAME, STYPE, MIP, IERR) 40.31 IF (IERR.NE.0) GOTO 800 ! ! output points exist ! ALLOCATE(ORQTMP) 40.31 NREOQ = NREOQ + 1 40.31 IF (NREOQ.GT.MAX_OUTP_REQ) CALL MSGERR (2, 40.31 40.13 & 'too many output requests') 40.13 ! CALL INKEYW ('STA', 'SPEC2D') 20.67 IF (KEYWIS ('FS1D') .OR. KEYWIS('SPEC1D')) THEN 20.67 RTYPE = 'SPE1' 20.28 ELSE CALL IGNORE ('SFD') 20.28 CALL IGNORE ('SPEC2D') 20.67 RTYPE = 'SPEC' ENDIF CALL INKEYW ('STA', 'ABS') 40.03 IF (KEYWIS ('REL')) THEN 40.03 RTYPE(3:3) = 'R' 20.28 ELSE CALL IGNORE ('ABS') 40.03 ENDIF ORQTMP%OQR(1) = -1. 40.31 ORQTMP%OQR(2) = -1. 40.31 ORQTMP%RQTYPE = RTYPE 40.31 CALL INCSTR ('FNAME', FILENM, 'REQ', ' ') NREF = 0 ! --- append node number to FILENM in case of 40.30 ! parallel computing 40.30 IF ( PARLL ) THEN 40.30 ILPOS = INDEX ( FILENM, ' ' )-1 40.30 WRITE(FILENM(ILPOS+1:ILPOS+4),33) INODE 40.30 END IF 40.30 ORQTMP%PSNAME = PSNAME 40.31 ORQTMP%OQI(1) = NREF 40.31 ORQTMP%OQI(2) = NREOQ 40.31 OUTP_FILES(NREOQ) = FILENM 40.31 40.13 ! NVAR = 0 40.31 ORQTMP%OQI(3) = NVAR 40.31 ALLOCATE(ORQTMP%IVTYP(0)) 40.31 ALLOCATE(ORQTMP%FAC(0)) 40.31 ! read types of variables to be printed in the table CALL INKEYW ('STA', ' ') 30.00 IF (KEYWIS ('OUT')) THEN 40.03 IF (NSTATM.EQ.0) CALL MSGERR (3, & 'time information not allowed in stationary mode') NSTATM = 1 CALL INCTIM (ITMOPT, 'TBEG', ORQTMP%OQR(1), 'REQ', 0.) 40.31 30.00 CALL ININTV ('DELT', ORQTMP%OQR(2), 'REQ', 0.) 40.31 30.00 IF (NSTATM.EQ.0) CALL MSGERR (2, & 'time input not allowed in stationary mode') ENDIF NULLIFY(ORQTMP%NEXTORQ) 40.31 IF ( .NOT.LORQ ) THEN 40.31 FORQ = ORQTMP 40.31 CORQ => FORQ 40.31 LORQ = .TRUE. 40.31 ELSE 40.31 CORQ%NEXTORQ => ORQTMP 40.31 CORQ => ORQTMP 40.31 END IF 40.31 GOTO 800 END IF ! -------------------------------------------------------------------------- ! NESTout 'sname' 'fname' & ! (OUTPUT [tbegnst] [deltnst] SEC/MIN/HR/DAY) ! -------------------------------------------------------------------------- ! NEST output for nesting of models VER. 20.63 IF (KEYWIS ('NEST')) THEN 40.00 ! ! ====================================================================== ! ! NESTout 'sname' 'fname' & ! ! | -> Sec | ! OUTput [tbegnst] [deltnst] < MIn > ! | HR | ! | DAy | ! ! ======================================================================= ! IF (ONED) THEN 32.02 CALL MSGERR (2,' Illegal keyword (NEST) in'// 32.02 & ' combination with 1D-computation') 32.02 GOTO 800 32.02 ELSE 32.02 CALL SWNMPS (PSNAME, STYPE, MIP, IERR) 40.31 IF (IERR.NE.0) GOTO 800 IF (STYPE .NE. 'N') THEN CALL MSGERR(2,'Set of output locations is not correct type') GOTO 800 ENDIF ! output points exist ALLOCATE(ORQTMP) 40.31 NREOQ = NREOQ + 1 40.31 IF (NREOQ.GT.MAX_OUTP_REQ) CALL MSGERR (2, 40.31 40.13 & 'too many output requests') 40.13 ORQTMP%OQR(1) = -1. 40.31 ORQTMP%OQR(2) = -1. 40.31 RTYPE = 'SPRC' 40.00 ORQTMP%RQTYPE = RTYPE 40.31 CALL INCSTR ('FNAME', FILENM, 'REQ', ' ') NREF = 0 ! --- append node number to FILENM in case of 40.30 ! parallel computing 40.30 IF ( PARLL ) THEN 40.30 ILPOS = INDEX ( FILENM, ' ' )-1 40.30 WRITE(FILENM(ILPOS+1:ILPOS+4),33) INODE 40.30 END IF 40.30 ORQTMP%PSNAME = PSNAME 40.31 ORQTMP%OQI(1) = NREF 40.31 ORQTMP%OQI(2) = NREOQ 40.31 OUTP_FILES(NREOQ) = FILENM 40.31 40.13 NVAR = 0 40.31 ORQTMP%OQI(3) = NVAR 40.31 ALLOCATE(ORQTMP%IVTYP(0)) 40.31 ALLOCATE(ORQTMP%FAC(0)) 40.31 ! CALL INKEYW ('STA', ' ') 30.00 IF (KEYWIS ('OUT')) THEN 40.03 IF (NSTATM.EQ.0) CALL MSGERR (3, & 'time information not allowed in stationary mode') NSTATM = 1 CALL INCTIM (ITMOPT, 'TBEG', ORQTMP%OQR(1), 'REQ', 0.) 40.31 30.00 CALL ININTV ('DELT', ORQTMP%OQR(2), 'REQ', 0.) 40.31 30.00 IF (NSTATM.EQ.0) CALL MSGERR (2, & 'time input not allowed in stationary mode') ENDIF ! NULLIFY(ORQTMP%NEXTORQ) 40.31 IF ( .NOT.LORQ ) THEN 40.31 FORQ = ORQTMP 40.31 CORQ => FORQ 40.31 LORQ = .TRUE. 40.31 ELSE 40.31 CORQ%NEXTORQ => ORQTMP 40.31 CORQ => ORQTMP 40.31 END IF 40.31 GOTO 800 ENDIF 32.02 ENDIF ! ------------------------------------------------------- ! command not found: RETURN 800 FOUND = .TRUE. RETURN !* end of subroutine SWREOQ ** END !*********************************************************************** ! * INTEGER FUNCTION SIRAY (DP, XP1, YP1, XP2, YP2, XX, YY, BOTDEP, 30.70 & BOTLEV, WATLEV) 30.70 ! * !*********************************************************************** ! USE OCPCOMM4 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. AUTHORS ! ! 30.72: IJsbrand Haagsma ! 40.41: Marcel Zijlema ! ! 1. UPDATE ! ! 00.00, Mar. 87: heading added, name of routine changed from ! IRAAI in SIRAY ! 30.72, Oct. 97: logical function EQREAL introduced for floating point ! comparisons ! 30.70, Nov. 97: changed into INTEGER function ! test output added ! arguments BOTDEP, BOTLEV, WATLEV added ! 40.03, Nov. 99: X2= etc. moved out of IF-ENDIF group ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. PURPOSE ! ! Searching the first point on a ray where the depth is DP ! ! 3. METHOD ! ! --- ! ! 4. PARAMETERLIST ! ! DP REAL input depth ! XP1 REAL input X-coordinate start point of ray ! YP1 REAL input Y-coordinate start point of ray ! XP2 REAL input X-coordinate end point of ray ! YP2 REAL input Y-coordinate end point of ray ! XX REAL input X-coordinate point with depth DP ! YY REAL input Y-coordinate point with depth DP ! ! 5. SUBROUTINES CALLING ! ! SWREPS (SWAN/READ) ! ! 6. SUBROUTINES USED ! ! SVALQI (SWAN/READ) ! ! 7. ERROR MESSAGES ! ! --- ! ! 8. REMARKS ! ! --- ! ! 9. STRUCTURE ! ! ---------------------------------------------------------------- ! Give SIRAY initial value 0 ! Compute stepsize, raylength and number of steps along the ray ! Compute bottom coordinates of startpoint as number of meshes ! Call SVALQI to interpolate depth in startpoint of ray ! For every step along the ray do ! Compute coordinates of the intermediate point in problem ! grid and bottom grid ! Call SVALQI to interpolate the depth for this point ! If the required depth is in the interval, then ! Compute coordinates of the point with depth DP ! Set SIRAY 1 ! Else ! Coordinates and depth at start of new interval are values ! at end of old interval ! ---------------------------------------------------------------- ! 10. SOURCE TEXT ! LOGICAL EQREAL, BOTDEP 30.72 REAL BOTLEV(*), WATLEV(*) 30.70 SAVE IENT DATA IENT/0/ CALL STRACE (IENT,'SIRAY') ! SIRAY = 0 DIFDEP = 1e+10 DSTEP = MIN(DXG(1), DYG(1)) RAYLEN = SQRT ((XP2-XP1)*(XP2-XP1) + (YP2-YP1)*(YP2-YP1)) NSTEP = 1 + INT(1.5*RAYLEN/DSTEP + 0.5) 30.70 ! DO 10 JJ = 0, NSTEP 30.70 X3 = XP1 + REAL(JJ)*(XP2-XP1)/REAL(NSTEP) Y3 = YP1 + REAL(JJ)*(YP2-YP1)/REAL(NSTEP) IF (BOTDEP) THEN D3 = SVALQI (X3, Y3, 1, BOTLEV, 1, 0, 0) 30.70 ELSE D3 = SVALQI (X3, Y3, 1, BOTLEV, 1, 0, 0) + WLEV 30.70 IF (LEDS(7).GE.2) 30.70 & D3 = D3 + SVALQI (X3, Y3, 7, WATLEV, 1, 0, 0) 30.70 ENDIF IF (ITEST.GE.160) WRITE (PRTEST, 14) X3+XOFFS, Y3+YOFFS, D3 30.70 14 FORMAT (' SIRAY, scan point', 2(1X,F8.0), 1X, F8.2) 30.70 IF (ABS(D3-DP).LT.DIFDEP) THEN 10.20 DIFDEP=ABS(D3-DP) 10.20 JDMINMAX=JJ 10.20 ENDIF IF (JJ.GT.0) THEN IF ((DP-D2)*(DP-D3).LE.0) THEN 40.03 IF (EQREAL(D2,D3)) THEN 30.72 XX = X2 YY = Y2 ELSE XX = X2+(X3-X2)*(D2-DP)/(D2-D3) YY = Y2+(Y3-Y2)*(D2-DP)/(D2-D3) ENDIF SIRAY = 1 GOTO 20 ENDIF 40.03 ENDIF 40.03 X2 = X3 Y2 = Y3 D2 = D3 10 CONTINUE ! ! exact depth not found, take closest value: ! X3 = XP1 + REAL(JDMINMAX)*(XP2-XP1)/REAL(NSTEP) 10.20 Y3 = YP1 + REAL(JDMINMAX)*(YP2-YP1)/REAL(NSTEP) 10.20 XX = X3 10.20 YY = Y3 10.20 ! 20 IF (ITEST.GE.140) WRITE (PRTEST, 24) XX+XOFFS, YY+YOFFS 30.70 24 FORMAT (' SIRAY, result ', 2(1X,F8.0)) 30.70 RETURN ! * end of function SIRAY * END !************************************************************************ ! * SUBROUTINE SWNMPS (PSNAME, PSTYPE, MIP, IERR) 40.31 ! * !************************************************************************ ! USE OCPCOMM1 40.41 USE OCPCOMM4 40.41 USE SWCOMM1 40.41 USE OUTP_DATA 40.31 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! 40.41: Marcel Zijlema ! ! 1. UPDATE ! ! Oct. 1996, ver. 30.50: new subr. ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. PURPOSE ! ! Read name of set of output points; get type and number of ! points in the set ! ! 3. METHOD ! ! ! 4. PARAMETERLIST ! ! PSNAME char output name ! PSTYPE char output type ! MIP int output number of points ! ! 5. SUBROUTINES CALLING ! ! SPREOQ ! ! 6. SUBROUTINES USED ! ! INCSTR (Ocean Pack) ! ! 7. ERROR MESSAGES ! ! --- ! ! 8. REMARKS ! ! --- ! ! 9. STRUCTURE ! ! ! 10. SOURCE TEXT ! INTEGER MIP 40.31 CHARACTER PSNAME *(*), PSTYPE *1 40.31 TYPE(OPSDAT), POINTER :: CUOPS 40.31 SAVE IENT DATA IENT/0/ CALL STRACE (IENT,'SWNMPS') ! IERR = 0 CALL INCSTR ('SNAME', PSNAME, 'STA', 'BOTTGRID') IF (LENCST.GT.8) CALL MSGERR (2, 'SNAME is too long') CUOPS => FOPS 40.31 DO 40.31 IF (CUOPS%PSNAME.EQ.PSNAME) EXIT 40.31 IF (.NOT.ASSOCIATED(CUOPS%NEXTOPS)) THEN 40.31 CALL MSGERR(2, 'Set of output locations is not known') 40.31 GOTO 900 40.31 END IF 40.31 CUOPS => CUOPS%NEXTOPS 40.31 END DO 40.31 PSTYPE = CUOPS%PSTYPE 40.31 IF (PSTYPE.EQ.'F' .OR. PSTYPE.EQ.'H') THEN 40.31 MIP = CUOPS%OPI(1) * CUOPS%OPI(2) 40.31 ! get direction of frame in case of coordinates plotting ALPQ = CUOPS%OPR(5) 40.31 ELSE MIP = CUOPS%MIP 40.31 ALPQ = 0. 40.31 ENDIF 800 IF (ITEST.GE.100) WRITE (PRTEST, 802) PSNAME, PSTYPE, MIP 802 FORMAT (' exit SWNMPS, name:', A8, ' type:', A1, & ' num of p:', I5) RETURN 900 PSTYPE = ' ' MIP = 0 IERR = 1 RETURN END !************************************************************************ ! * SUBROUTINE SVARTP (IVTYPE) ! * !************************************************************************ ! USE SWCOMM1 40.41 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 32.02: Roeland Ris & Cor van der Schelde ! 40.03: Nico Booij ! 40.04: Annette Kieftenburg ! ! 1. Updates ! ! 10.09, Aug. 94: output quantity RPER added ! 20.61, Sep. 95: quantities TM02 and FWID added ! 20.67, Dec. 95: FWID renamed FSPR (freq. spread) ! 32.02, Feb. 98: 1D-version introduced ! 40.00, Apr. 98: subr simplified using new array OVKEYW ! 40.03, Sep. 00: inconsistency with manual corrected ! ! 2. Purpose ! ! Converting keyword into integer ! ! 3. Method ! ! This subroutine determines an integer value indicating the ! required output variable from the keyword denoting the same ! for storage in array with output requests. ! ! 4. PARAMETERLIST ! ! IVTYPE INT output type number output variable ! ! 5. Subroutines calling ! ! SPROUT ! ! 6. Subroutines used ! ! KEYWIS (Ocean Pack) ! ! 10. Error messages ! ! --- ! ! 11. Remarks ! ! --- ! ! 12. Structure ! ! ----------------------------------------------------------------- ! If the keyword is equal to given string, then ! IVTYPE is given integer value ! ----------------------------------------------------------------- ! ! 13. Source text ! LOGICAL KEYWIS SAVE IENT DATA IENT/0/ CALL STRACE (IENT,'SVARTP') ! IVTYPE = 0 ! CALL INKEYW ('STA', 'ZZZZ') ! check if given keyword corresponds to output quantity DO IVT = NMOVAR, 1, -1 ! loop in reverse order to check more specific names first ! e.g. HSWE before HS IF (KEYWIS (OVKEYW(IVT))) THEN 40.00 IVTYPE = IVT 40.00 GOTO 40 ENDIF ENDDO ! aliases: IF (KEYWIS ('PPER')) IVTYPE = 12 40.00 IF (KEYWIS ('RPER')) IVTYPE = 28 40.00 IF (KEYWIS ( 'DTM')) IVTYPE = 31 40.00 IF (KEYWIS ('FWID')) IVTYPE = 33 40.00 ! keyword OUTPUT means that output times will be entered 40.00 IF (KEYWIS ('OUT')) IVTYPE = 98 40.03 ! keyword ZZZZ means end of list of output quantities 40.00 IF (KEYWIS ('ZZZZ')) IVTYPE = 99 ! IF (IVTYPE .EQ. 0) CALL WRNKEY ! 40 RETURN ! end of subroutine SVARTP * END !************************************************************************ ! * SUBROUTINE SWBOUN ( XCGRID, YCGRID, KGRPNT, XYTST, KGRBND ) 40.31 ! * !************************************************************************ ! USE OCPCOMM2 40.41 USE OCPCOMM4 40.41 USE SWCOMM1 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 USE M_BNDSPEC 40.31 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 30.73: Nico Booij ! 30.81: Annette Kieftenburg ! 30.82: IJsbrand Haagsma ! 30.90: IJsbrand Haagsma (Equivalence version) ! 34.01: Jeroen Adema ! 40.02: IJsbrand Haagsma ! 40.03, 40.13: Nico Booij ! 40.05: Ekaterini E. Kriezi ! 40.31: Marcel Zijlema ! 40.41: Marcel Zijlema ! ! 1. Updates ! ! 30.73, Nov. 97: New subroutine, replacing code in subr. SWREAD (file SWANPRE1) ! 30.90, Oct. 98: Introduced EQUIVALENCE POOL-arrays ! 30.82, Oct. 98: Updated description several arrays ! 30.81, Nov. 98: Adjustment for 1-D case of new boundary conditions ! 34.01, Feb. 99: Introducing STPNOW ! 30.81, Apr. 99: Prevent negative powers for cosine directional spreading (DSPR); ! prevented DSPR > 360 and DSPR < 0 (except for exception value). ! 30.82, July 99: Used EQREAL for real equality comparisons ! 40.05, Aug 00: WW3 boundary nesting command, in Swan nesting option ! adding of a new option (same as WW3 command) ! 40.03, Sep. 00: inconsistency with manual corrected ! 40.02, Oct. 00: WWIII added as keyword (will appear in the manual) ! 40.13, Nov. 01: determination of side corrected (iside=3) ! 40.31, Nov. 03: removing POOL-mechanism, reconsideration of this subroutine ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. Purpose ! ! Reading and processing BOUNDARY command ! ! 3. Method ! ! ! 4. Argument variables ! ! i XCGRID: Coordinates of computational grid in x-direction 30.82 ! i YCGRID: Coordinates of computational grid in y-direction 30.82 ! REAL XCGRID(MXC,MYC), YCGRID(MXC,MYC) 30.82 ! ! KGRBND int inp grid indices of boundary points 40.31 ! KGRPNT int inp indirect addresses of grid points ! XYTST int inp ix, iy of test points ! INTEGER KGRPNT(MXC,MYC) INTEGER XYTST(*), KGRBND(*) ! ! 5. Parameter variables ! ! ! 6. Local variables ! INTEGER IENT,KOUNTR,IX1,IY1,IX2,IY2 INTEGER MM,IX,IY,ISIDM,ISIDE,KC,KC2,KC1,IX3,IY3,MP INTEGER IP,II,NBSPSS,NFSEQ,IKO,IKO2,IBSPC1,IBSPC2 REAL CRDP, CRDM, SOMX, SOMY REAL XP,YP,XC,YC,RR,DIRSI,COSDIR,SINDIR,DIRSID,DIRREF REAL RLEN1,RDIST,RLEN2,XC1,YC1,XC2,YC2,W1 LOGICAL KEYWIS, LOCGRI, CCW, BPARF, BOUNPT,DONALL LOGICAL LFRST1, LFRST2, LFRST3 40.31 INTEGER UPLO, NUMP LOGICAL, SAVE :: LBFILS = .FALSE. 40.31 LOGICAL, SAVE :: LBS = .FALSE. 40.31 LOGICAL, SAVE :: LBGP = .FALSE. 40.31 TYPE(BSPCDAT), POINTER :: BFLTMP 40.31 TYPE(BSPCDAT), SAVE, POINTER :: CUBFL 40.31 TYPE(BSDAT), POINTER :: BSTMP 40.31 TYPE(BSDAT), SAVE, POINTER :: CUBS 40.31 TYPE(BGPDAT), POINTER :: BGPTMP 40.31 TYPE XYPT 40.31 INTEGER :: JX, JY TYPE(XYPT), POINTER :: NEXTXY END TYPE XYPT TYPE(XYPT), TARGET :: FRST 40.31 TYPE(XYPT), POINTER :: CURR, TMP 40.31 ! ! 8. Subroutines used ! ! Ocean Pack command reading routines ! BOUNPT ! LOGICAL STPNOW 34.01 LOGICAL EQREAL ! ! 9. Subroutines calling ! ! SWREAD ! ! 10. Error messages ! ! ! 11. Remarks ! ! data concerning boundary files are stored in array BFILES ! see subr BCFILE for details ! ! 12. Structure ! ! ----------------------------------------------------------------- ! Read keyword ! Case keyword = ! 'SHAPE': Read spectral shape parameters ! 'WAMN': Read filename ! Open WAM nesting file ! Put file characteristics into array BFILES ! 'WW3N': Read filename ! Open WW3 nesting file ! Put file characteristics into array BFILES ! 'NEST': Read filename ! Call BCFILE to obtain file characteristics ! 'SIDE': read side of boundary ! 'SEG': read set of points on boundary of comp. grid ! Read keyword ! If keyword is 'UNIF' ! Then Read keyword ! If keyword is 'PAR' ! Then read integral wave parameters ! Call SSHAPE to generate spectrum ! put spectrum into array BSPECS ! Else {keyword is 'FILE'} ! Read filename ! Call BCFILE to obtain file characteristics ! Else {keyword is 'VAR'} ! Read keyword ! If keyword is 'PAR' ! Then Repeat until list is exhausted ! Read length and integral wave parameters ! Call SSHAPE to generate spectrum ! put spectrum into array BSPECS ! Else {keyword is 'FILE'} ! Repeat until list is exhausted ! Read filename ! Call BCFILE to obtain file characteristics ! ----------------------------------------------------------------- ! ! 13. Source text SAVE IENT DATA IENT/0/ CALL STRACE (IENT,'SWBOUN') ! CALL INKEYW ('REQ',' ') IF (KEYWIS ('SHAP')) THEN ! ! specification of the spectral shape ! ! ========================================================================= ! ! | JONswap [gamma] | ! | | | -> PEAK | ! BOUNdspec SHAPe < PM > < > & ! | | | MEAN | ! | GAUSs [sigfr] | ! | | ! | BIN | ! ! | DEGRees | ! DSPR < > ! | -> POWer | ! ! ========================================================================= ! CALL INKEYW ('STA', 'JON') IF (KEYWIS ('JON')) THEN FSHAPE = 2 CALL INREAL ('GAMMA', PSHAPE(1), 'STA', 3.3) 40.00 ELSE IF (KEYWIS ('BIN')) THEN FSHAPE = 3 ELSE IF (KEYWIS ('PM')) THEN FSHAPE = 1 ELSE IF (KEYWIS ('GAUS')) THEN FSHAPE = 4 CALL INREAL ('SIGFR', SIGMAG, 'STA', 0.01) ! convert from Hz to rad/s: PSHAPE(2) = PI2 * SIGMAG 40.00 ENDIF ! PEAK or MEAN frequency CALL INKEYW ('STA', ' ') IF (KEYWIS('MEAN')) THEN FSHAPE = -FSHAPE ELSE CALL IGNORE ('PEAK') ENDIF ! directional distribution given by DEGR or by POWER CALL IGNORE ('DSPR') CALL INKEYW ('STA', 'POW') IF (KEYWIS('DEGR')) THEN DSHAPE = 1 ELSE CALL IGNORE ('POW') DSHAPE = 2 ENDIF IF (ITEST.GE.30) WRITE (PRINTF,6100) FSHAPE, DSHAPE 6100 FORMAT (' Shape of inc. spectrum, Freq:', I2, ' ; Dir:', I2) ! ELSE IF (KEYWIS ('WAMN')) THEN ! ! SWAN in WAM nesting 30.04 ! ! ======================================================================= ! ! |-> CRAY | ! | UNFormatted < > | ! | | WKstat | | ! | | ! BOUNdnest2 WAMNest 'fname' < > [xgc] [ygc] ! | | ! | FREE | ! ! ======================================================================= ! IF (MXC .LE. 0) THEN CALL MSGERR(3, & ' Computational grid must be given previously') GOTO 900 ENDIF IF (MCGRD .LE. 1) THEN CALL MSGERR(3, & ' Read depth before entering boundary values') GOTO 900 ENDIF ! NBFILS = NBFILS + 1 ALLOCATE(BFLTMP) 40.31 CALL INCSTR ('FNAME',FILENM,'REQ', ' ') CALL BCWAMN (FILENM, 'NEST', BFLTMP, LBGP, 40.31 & XCGRID, YCGRID, KGRPNT, XYTST) 40.31 IF (STPNOW()) RETURN 34.01 NULLIFY(BFLTMP%NEXTBSPC) 40.31 IF ( .NOT.LBFILS ) THEN 40.31 FBNDFIL = BFLTMP 40.31 CUBFL => FBNDFIL 40.31 LBFILS = .TRUE. 40.31 ELSE 40.31 CUBFL%NEXTBSPC => BFLTMP 40.31 CUBFL => BFLTMP 40.31 END IF 40.31 ! ELSE IF (KEYWIS('WW3').OR.KEYWIS('WWIII')) THEN 40.02 ! ! SWAN in WaveWatch nesting 40.05 ! ! ======================================================================= ! | -> CLOS | ! BOUNdnest2 WWIII 'fname' < > [xgc] [ygc] 40.02 ! | OPEN | ! ======================================================================= ! IF (MXC .LE. 0) THEN CALL MSGERR(3, & ' Computational grid must be given previously') GOTO 900 ENDIF IF (MCGRD .LE. 1) THEN CALL MSGERR(3, & ' Read depth before entering boundary values') GOTO 900 ENDIF ! NBFILS = NBFILS + 1 ALLOCATE(BFLTMP) 40.31 CALL INCSTR ('FNAME',FILENM,'STA', 'nest.ww3') ! if keyword is OPEN (boundary is not a closed contour) then ! DONALL is TRUE and the nesting boundary remain open ! else (default case) ! DONALL is FALSE and boundary is close and interpolation ! between the last and the first point will be done CALL INKEYW ('STA', 'CLOS') 40.05 IF (KEYWIS('OPEN')) THEN 40.05 DONALL = .TRUE. 40.05 ELSE IF (KEYWIS('CLOS')) THEN 40.05 DONALL = .FALSE. 40.05 ELSE 40.05 CALL WRNKEY 40.05 ENDIF 40.05 CALL BCWW3N (FILENM, 'NEST', BFLTMP, LBGP, 40.31 40.05 & XCGRID, YCGRID, KGRPNT, XYTST, KGRBND, 40.31 40.05 & DONALL) 40.31 40.05 IF (STPNOW()) RETURN NULLIFY(BFLTMP%NEXTBSPC) 40.31 IF ( .NOT.LBFILS ) THEN 40.31 FBNDFIL = BFLTMP 40.31 CUBFL => FBNDFIL 40.31 LBFILS = .TRUE. 40.31 ELSE 40.31 CUBFL%NEXTBSPC => BFLTMP 40.31 CUBFL => BFLTMP 40.31 END IF 40.31 ! ELSE IF (KEYWIS ('NE')) THEN ! ! Nesting SWAN model in larger SWAN model ! ========================================== ! | -> CLOS | 40.05 ! BOUNdnest1 NEST 'fname' < > ! | OPEN | 40.05 ! ========================================== IF (MXC .LE. 0) THEN CALL MSGERR(3, & ' Computational grid must be given previously') GOTO 900 ENDIF IF (MCGRD .LE. 1) THEN CALL MSGERR(3, & ' Read depth before entering boundary values') GOTO 900 ENDIF NBFILS = NBFILS + 1 ALLOCATE(BFLTMP) 40.31 CALL INCSTR ('FNAME',FILENM,'REQ', ' ') ! if keyword is OPEN then ! DONALL is TRUE and the nesting boundary remain open ! else (default case) ! DONALL is FALSE and boundary is close and interpolation between ! the last and the first point will be done CALL INKEYW ('STA', 'CLOS') 40.05 IF (KEYWIS('OPEN')) THEN 40.05 DONALL = .TRUE. 40.05 ELSE IF (KEYWIS('CLOS')) THEN 40.05 DONALL = .FALSE. 40.05 ELSE 40.05 CALL WRNKEY 40.05 ENDIF 40.05 CALL BCFILE (FILENM, 'NEST', BFLTMP, LBGP, 40.31 & XCGRID, YCGRID, KGRPNT, XYTST, KGRBND, 40.31 & DONALL) 40.05 IF (STPNOW()) RETURN 34.01 NULLIFY(BFLTMP%NEXTBSPC) 40.31 IF ( .NOT.LBFILS ) THEN 40.31 FBNDFIL = BFLTMP 40.31 CUBFL => FBNDFIL 40.31 LBFILS = .TRUE. 40.31 ELSE 40.31 CUBFL%NEXTBSPC => BFLTMP 40.31 CUBFL => BFLTMP 40.31 END IF 40.31 ! ELSE ! ! parametric or file boundary condition ! ! ====================================================================== ! ! | North | ! | NW | ! | West | ! | SW | | -> CCW | ! | -> SIDE < South > < > | ! | | SE | | CLOCKWise | | ! | | East | | ! | | NE | | ! BOUNdary < > & ! | | -> XY < [x] [y] > | | ! | SEGment < > | ! | IJ < [i] [j] > | ! ! ! | PAR [hs] [per] [dir] [dd] | ! | UNIForm < > | ! | | FILE 'fname' [seq] | | ! < > ! | | PAR < [len] [hs] [per] [dir] [dd] > | | ! | VARiable < > | ! | FILE < [len] 'fname' [seq] > | ! ! ====================================================================== ! IF (MXC .LE. 0) THEN CALL MSGERR(3, & ' Computational grid must be given previously') GOTO 900 ENDIF IF (MCGRD .LE. 1) THEN CALL MSGERR(3, & ' Read depth before entering boundary values') GOTO 900 ENDIF ! ! first define side or segment ! ! *** definition of boundary segment *** ! CALL INKEYW ('REQ',' ') IF (KEYWIS ('STAT')) THEN CALL MSGERR (1, 'keyword STAT ignored') CALL INKEYW ('REQ',' ') ENDIF KOUNTR = 0 FRST%JX = 0 40.31 FRST%JY = 0 40.31 NULLIFY(FRST%NEXTXY) 40.31 CURR => FRST 40.31 IF (KEYWIS ('SEG')) THEN IERR = 0 CALL INKEYW ('STA','XY') IF (KEYWIS('XY') .OR. KEYWIS ('LOC')) THEN LOCGRI = .TRUE. ELSE IF (KEYWIS('IJ') .OR. KEYWIS ('GRI')) THEN LOCGRI = .FALSE. ELSE CALL WRNKEY ENDIF ! loop over points describing the segment IX1 = 1 IY1 = 1 LFRST1 = .TRUE. 40.31 DO IF (LOCGRI) THEN CALL READXY ('XP','YP',XP,YP, 'REP', -1.E10, -1.E10) IF (XP.LT.-.9E10) GOTO 42 CALL CVMESH (XP, YP, XC, YC, KGRPNT, XCGRID, YCGRID, & KGRBND) 40.00 IX2 = NINT(XC) + 1 IY2 = NINT(YC) + 1 40.00 IF (.NOT.BOUNPT(IX2,IY2,KGRPNT)) THEN CALL MSGERR (2, 'invalid boundary point') WRITE (PRTEST, 38) XP+XOFFS, YP+YOFFS, XC, YC, IX2, IY2 40.00 38 FORMAT (' segment point ', 2F10.2, & ' grid ', 2F8.2, 2I4) ENDIF ELSE CALL ININTG ('I' , IX2, 'REP', -1) 40.03 IF (IX2 .LT. 0) GOTO 42 40.00 CALL ININTG ('J' , IY2, 'REQ', 0) 40.03 IX2 = IX2 + 1 40.00 IY2 = IY2 + 1 ENDIF IF (ITEST.GE.80) WRITE (PRTEST, 38) XCGRID(IX2,IY2)+XOFFS, 40.00 & YCGRID(IX2,IY2)+YOFFS, XC, YC, IX2-1, IY2-1 40.00 ! ! --- generate intermediate points on the segment IF (IX2 .GT. 0 .AND. IX2 .LE. MXC .AND. & IY2 .GT. 0 .AND. IY2 .LE. MYC) THEN IF (LFRST1) THEN MM = 1 LFRST1 = .FALSE. ELSE MM = MAX (ABS(IX2-IX1), ABS(IY2-IY1)) ENDIF DO IP = 1, MM RR = REAL(IP) / REAL(MM) ! IF (.NOT. ONED) THEN 30.81 IX = IX1 + NINT(RR*REAL(IX2-IX1)) IY = IY1 + NINT(RR*REAL(IY2-IY1)) ELSE 30.81 IX = IX1 + NINT(RR*REAL(IX2-IX1)) 30.81 IY = IY1 30.81 END IF 30.81 ! IF (ITEST.GE.80) WRITE (PRTEST, *) ' b. point ', & RR, IX, IY IF (KGRPNT(IX,IY) .GT. 1) THEN KOUNTR = KOUNTR + 1 ALLOCATE(TMP) 40.31 TMP%JX = IX 40.31 TMP%JY = IY 40.31 NULLIFY(TMP%NEXTXY) 40.31 CURR%NEXTXY => TMP 40.31 CURR => TMP 40.31 ENDIF ENDDO ELSE CALL MSGERR (2, 'Boundary point outside comp. grid') ENDIF IX1 = IX2 IY1 = IY2 END DO 42 IF (KOUNTR.EQ.0) & CALL MSGERR(1,'No points of the boundary found') ELSE ! boundary condition on one side of the computational grid IF (OPTG.NE.1) THEN 40.31 !CJK:03-03-2004:Only a warning must be given, swan must not stop on keyword side !CJK request of Arjen L., so error code is reduced from 3 to 1 !CJK CALL MSGERR(3, 40.31 !CJK & ' keyword SIDE can be used for regular grid only') 40.31 CALL MSGERR(1, 40.31 & ' keyword SIDE can be used for regular grid only') 40.31 !CJK:END:03-03-2004 END IF 40.31 CALL IGNORE ('SIDE') CALL INKEYW ('REQ',' ') ! *** specification of side for which boundary *** ! *** condition is given *** UPLO = 0 IF (KEYWIS ('NW')) THEN DIRSI = 45. ELSE IF (KEYWIS ('SW')) THEN DIRSI = 135. ELSE IF (KEYWIS ('SE')) THEN DIRSI = -135. ELSE IF (KEYWIS ('NE')) THEN DIRSI = -45. ELSE IF (KEYWIS ('N')) THEN DIRSI = 0. ELSE IF (KEYWIS ('W')) THEN DIRSI = 90. ELSE IF (KEYWIS ('S')) THEN DIRSI = 180. ELSE IF (KEYWIS ('E')) THEN DIRSI = -90. ELSE IF (KEYWIS ('LO')) THEN UPLO = -1 ELSE IF (KEYWIS ('UP')) THEN UPLO = 1 ELSE CALL WRNKEY ENDIF IF (UPLO.NE.0) THEN ! older, obsolete option CALL MSGERR (1, 'UPper and LOwer are becoming obsolete') CALL INKEYW ('REQ', ' ') IF (KEYWIS('Y')) THEN IF (UPLO.GT.0) THEN ISIDM = 3 CCW = .FALSE. ELSE ISIDM = 1 CCW = .TRUE. ENDIF ELSE CALL IGNORE ('X') IF (UPLO.GT.0) THEN ISIDM = 2 CCW = .TRUE. ELSE ISIDM = 4 CCW = .FALSE. ENDIF ENDIF CALL INKEYW ('STA', ' ') IF (KEYWIS ('JON')) THEN FSHAPE = 2 CALL INREAL ('GAMMA', PSHAPE(1), 'STA', 3.3) ELSE IF (KEYWIS ('PM')) THEN FSHAPE = 1 ELSE IF (KEYWIS ('GAUS')) THEN FSHAPE = 4 CALL INREAL ('SIGFR', SIGMAG, 'STA', 0.01) ! convert from Hz to rad/s: PSHAPE(2) = PI2 * SIGMAG ENDIF DSHAPE = 2 GOTO 90 ENDIF ! ! select side in the chosen direction ! CRDM = -1.E10 ISIDM = 0 IF (ONED) THEN 40.00 COSDIR = COS(PI*(DNORTH+DIRSI)/180.) SINDIR = SIN(PI*(DNORTH+DIRSI)/180.) DO ISIDE = 1, 4 SOMX = 0. SOMY = 0. NUMP = 0 IF (ISIDE.EQ.2) THEN KC = KGRPNT(MXC,1) IF (KC.GT.1) THEN SOMX = XCGRID(MXC,1) SOMY = YCGRID(MXC,1) NUMP = 1 ENDIF ELSE IF (ISIDE.EQ.4) THEN KC = KGRPNT(1,1) IF (KC.GT.1) THEN SOMX = XCGRID(1,1) SOMY = YCGRID(1,1) NUMP = 1 ENDIF ENDIF IF (NUMP.GT.0) THEN CRDP = COSDIR*SOMX + SINDIR*SOMY ! side with largest CRDP is the one selected IF (CRDP.GT.CRDM) THEN CRDM = CRDP ISIDM = ISIDE ENDIF ENDIF ENDDO ELSE 40.00 DO ISIDE = 1, 4 SOMX = 0. SOMY = 0. NUMP = 0 IF (ISIDE.EQ.1) THEN DO IX = 1, MXC KC2 = KGRPNT(IX,1) IF (IX.GT.1) THEN IF (KC1.GT.1 .AND. KC2.GT.1) THEN 40.00 ! if both grid points at ends of a step are valid, then ! take DX and DY into account when determining direction SOMX = SOMX + XCGRID(IX,1)-XCGRID(IX-1,1) SOMY = SOMY + YCGRID(IX,1)-YCGRID(IX-1,1) NUMP = NUMP + 1 ENDIF ENDIF KC1 = KC2 40.03 ENDDO ELSE IF (ISIDE.EQ.2) THEN DO IY = 1, MYC KC2 = KGRPNT(MXC,IY) IF (IY.GT.1) THEN IF (KC1.GT.1 .AND. KC2.GT.1) THEN 40.00 SOMX = SOMX + XCGRID(MXC,IY)-XCGRID(MXC,IY-1) SOMY = SOMY + YCGRID(MXC,IY)-YCGRID(MXC,IY-1) NUMP = NUMP + 1 ENDIF ENDIF KC1 = KC2 40.03 ENDDO ELSE IF (ISIDE.EQ.3) THEN 40.00 DO IX = 1, MXC KC2 = KGRPNT(IX,MYC) IF (IX.GT.1) THEN IF (KC1.GT.1 .AND. KC2.GT.1) THEN 40.00 SOMX = SOMX + XCGRID(IX-1,MYC)-XCGRID(IX,MYC) 40.13 SOMY = SOMY + YCGRID(IX-1,MYC)-YCGRID(IX,MYC) 40.13 NUMP = NUMP + 1 ENDIF ENDIF KC1 = KC2 40.03 ENDDO ELSE IF (ISIDE.EQ.4) THEN DO IY = 1, MYC KC2 = KGRPNT(1,IY) IF (IY.GT.1) THEN IF (KC1.GT.1 .AND. KC2.GT.1) THEN 40.00 SOMX = SOMX + XCGRID(1,IY-1)-XCGRID(1,IY) SOMY = SOMY + YCGRID(1,IY-1)-YCGRID(1,IY) NUMP = NUMP + 1 ENDIF ENDIF KC1 = KC2 40.03 ENDDO ENDIF IF (NUMP.GT.0) THEN DIRSID = ATAN2(SOMY,SOMX) DIRREF = PI*(DNORTH+DIRSI)/180. IF (CVLEFT) THEN CRDP = COS(DIRSID - 0.5*PI - DIRREF) ELSE CRDP = COS(DIRSID + 0.5*PI - DIRREF) ENDIF ! side with largest CRDP is the one selected IF (CRDP.GT.CRDM) THEN CRDM = CRDP ISIDM = ISIDE ENDIF ENDIF IF (ITEST.GE.60) WRITE (PRTEST, 151) ISIDE, NUMP, & SOMX, SOMY, DIRSID*180/PI, DIRREF*180/PI, CRDP, CVLEFT 40.13 151 FORMAT (' side ', 2I4, 2(1X,E11.4), 2(1X,F5.0), 2X, & F6.3, 2X, L1) 40.13 ENDDO ENDIF 40.00 IF (ISIDM.EQ.0) THEN CALL MSGERR (2, 'No open boundary found') ENDIF ! ! --- go along boundary clockwise or counterclockwise (default) CALL INKEYW ('STA', 'CCW') IF (KEYWIS('CLOCKW')) THEN CCW = .FALSE. ELSE CALL IGNORE ('CCW') CCW = .TRUE. ENDIF ! 90 IF (ISIDM.EQ.1) THEN IX1 = 1 IY1 = 1 IX2 = MXC IY2 = 1 ELSE IF (ISIDM.EQ.2) THEN IX1 = MXC IY1 = 1 IX2 = MXC IY2 = MYC ELSE IF (ISIDM.EQ.3) THEN IX1 = MXC IY1 = MYC IX2 = 1 IY2 = MYC ELSE IF (ISIDM.EQ.4) THEN IX1 = 1 IY1 = MYC IX2 = 1 IY2 = 1 ENDIF IF (.NOT.CCW .EQV. CVLEFT) THEN ! swap end points IX3 = IX1 IY3 = IY1 IX1 = IX2 IY1 = IY2 IX2 = IX3 IY2 = IY3 ENDIF IF (ITEST.GE.50) WRITE (PRINTF, 112) ISIDM, & IX1-1, IY1-1, XCGRID(IX1,IY1)+XOFFS, YCGRID(IX1,IY1)+YOFFS, 40.00 & IX2-1, IY2-1, XCGRID(IX2,IY2)+XOFFS, YCGRID(IX2,IY2)+YOFFS 40.00 112 FORMAT (' Selected side:', I2, ' from ', 2I4, 2F9.0, 40.00 & ' to ', 2I4, 2F9.0) 40.00 MP = MAX(ABS(IX2-IX1),ABS(IY2-IY1)) DO IP = 0, MP IF (MP.EQ.0) THEN 40.00 RR = 0. ELSE RR = REAL(IP) / REAL(MP) ENDIF IX = IX1 + NINT(RR*REAL(IX2-IX1)) IY = IY1 + NINT(RR*REAL(IY2-IY1)) IF (KGRPNT(IX,IY) .GT. 1) THEN KOUNTR = KOUNTR + 1 ALLOCATE(TMP) 40.31 TMP%JX = IX 40.31 TMP%JY = IY 40.31 NULLIFY(TMP%NEXTXY) 40.31 CURR%NEXTXY => TMP 40.31 CURR => TMP 40.31 ENDIF ENDDO ENDIF ! ! *** boundary condition from file, 1-d or 2-d spectrum ! CURR => FRST%NEXTXY 40.31 CALL INKEYW ('REQ',' ') IF (KEYWIS('UNIF') .OR. KEYWIS('CON') .OR. KEYWIS('PAR')) THEN CALL INKEYW('STA', 'PAR') IF (KEYWIS('PAR')) THEN CALL INREAL ('HS', SPPARM(1), 'REQ', 0.) CALL INKEYW ('STA', ' ') CALL INREAL ('PER', SPPARM(2), 'REQ', 0.) CALL INREAL ('DIR', SPPARM(3), 'REQ', 0.) IF (DSHAPE.EQ.1) THEN CALL INREAL ('DD', SPPARM(4), 'STA', 30.) IF ((SPPARM(4).GT.360. .OR. SPPARM(4).LT. 0.).AND. 30.81 & .NOT.(EQREAL(SPPARM(4),OVEXCV(16)))) THEN 30.82 CALL MSGERR (2,'Directional spreading is less than '// 30.81 & '0 or larger than 360 degrees, and no '// 30.81 & 'exception value') 30.81 END IF 30.81 ELSE CALL INREAL ('DD', SPPARM(4), 'STA', 2.) IF (SPPARM(4).LE. 0.) THEN 30.81 CALL MSGERR (2, 30.81 & 'Power of cosine is less or equal to zero') 30.81 END IF 30.81 IF (SPPARM(4)*DDIR**2/2. .GT. 1.) THEN 40.03 CALL MSGERR (2, 40.03 & 'distribution too narrow to be represented properly') 40.03 WRITE (PRINTF, 142) SQRT(2./SPPARM(4))*180./PI 40.03 142 FORMAT (' Advise: choose Dtheta < ', F8.3, ' degr') 40.03 END IF 40.03 ENDIF NBSPEC = NBSPEC + 1 IF (ITEST.GE.80) WRITE (PRTEST,*) ' bound. spectr.', & NBSPEC, (SPPARM(II), II=1,4) NBSPSS = NBSPEC ALLOCATE(BSTMP) 40.31 BSTMP%NBS = NBSPEC 40.31 BSTMP%FSHAPE = FSHAPE 40.31 BSTMP%DSHAPE = DSHAPE 40.31 BSTMP%SPPARM(1:4) = SPPARM(1:4) 40.31 NULLIFY(BSTMP%NEXTBS) 40.31 IF ( .NOT.LBS ) THEN 40.31 FBS = BSTMP 40.31 CUBS => FBS 40.31 LBS = .TRUE. 40.31 ELSE 40.31 CUBS%NEXTBS => BSTMP 40.31 CUBS => BSTMP 40.31 END IF 40.31 ELSE IF (KEYWIS('FILE') .OR. KEYWIS('SPEC')) THEN CALL INCSTR ('FNAME',FILENM,'REQ', ' ') ! generate new set of file data NBFILS = NBFILS + 1 NBSPSS = NBSPEC ALLOCATE(BFLTMP) 40.31 CALL BCFILE (FILENM, 'PNTS', BFLTMP, LBGP, 40.31 & XCGRID, YCGRID, KGRPNT, XYTST, KGRBND, 40.31 & DONALL) 40.05 IF (STPNOW()) RETURN 34.01 NULLIFY(BFLTMP%NEXTBSPC) 40.31 IF ( .NOT.LBFILS ) THEN 40.31 FBNDFIL = BFLTMP 40.31 CUBFL => FBNDFIL 40.31 LBFILS = .TRUE. 40.31 ELSE 40.31 CUBFL%NEXTBSPC => BFLTMP 40.31 CUBFL => BFLTMP 40.31 END IF 40.31 CALL ININTG('SEQ', NFSEQ, 'STA', 1) NBSPSS = NBSPSS + NFSEQ ENDIF DO IKO = 1, KOUNTR IX = CURR%JX 40.31 IY = CURR%JY 40.31 CURR => CURR%NEXTXY 40.31 ALLOCATE(BGPTMP) 40.31 BGPTMP%BGP(1) = KGRPNT(IX,IY) 40.31 BGPTMP%BGP(2) = 1 40.31 BGPTMP%BGP(3) = 1000 40.31 BGPTMP%BGP(4) = NBSPSS 40.31 BGPTMP%BGP(5) = 0 40.31 BGPTMP%BGP(6) = 1 40.31 NULLIFY(BGPTMP%NEXTBGP) 40.31 IF ( .NOT.LBGP ) THEN 40.31 FBGP = BGPTMP 40.31 CUBGP => FBGP 40.31 LBGP = .TRUE. 40.31 ELSE 40.31 CUBGP%NEXTBGP => BGPTMP 40.31 CUBGP => BGPTMP 40.31 END IF 40.31 ENDDO NBGRPT = NBGRPT + KOUNTR ELSE IF (KEYWIS('VAR')) THEN CALL INKEYW('STA', 'PAR') IF (KEYWIS('PAR')) THEN BPARF = .TRUE. ELSE IF (KEYWIS('FILE')) THEN BPARF = .FALSE. ENDIF RLEN1 = -1.E20 IKO = 1 RDIST = 0. IBSPC1 = 1 LFRST1 = .TRUE. LFRST2 = .TRUE. DO IF (LFRST1) THEN CALL INREAL('LEN', RLEN2, 'REQ', 0.) LFRST1 = .FALSE. ELSE CALL INREAL('LEN', RLEN2, 'STA', 1.E20) ENDIF IF (RLEN2.LT.0.9E20) THEN IF (IKO.GT.KOUNTR) THEN CALL MSGERR(1, & 'Length of segment short, boundary values ignored') 40.00 WRITE (PRINTF, 332) RDIST, RLEN2 332 FORMAT (' segment length=', F9.2, '; [len]=', F9.2) ENDIF IF (BPARF) THEN CALL INREAL ('HS', SPPARM(1), 'REQ', 0.) CALL INKEYW ('STA', ' ') CALL INREAL ('PER', SPPARM(2), 'REQ', 0.) CALL INREAL ('DIR', SPPARM(3), 'REQ', 0.) IF (DSHAPE.EQ.1) THEN CALL INREAL ('DD', SPPARM(4), 'STA', 30.) IF ((SPPARM(4).GT.360. .OR. SPPARM(4).LT. 0.).AND. 30.81 & .NOT.(EQREAL(SPPARM(4),OVEXCV(16)))) THEN 30.82 CALL MSGERR (2,'Directional spreading is less ' // 30.81 & 'than 0 or larger than 360 '// 30.81 & 'degrees and no exception value') 30.81 END IF 30.81 ELSE CALL INREAL ('DD', SPPARM(4), 'STA', 2.) IF (SPPARM(4).LE. 0.) THEN 30.81 CALL MSGERR (2,'Power of cosine is less or equal '// 30.81 & 'to zero') 30.81 END IF 30.81 ENDIF NBSPEC = NBSPEC + 1 IBSPC2 = NBSPEC ALLOCATE(BSTMP) 40.31 BSTMP%NBS = NBSPEC 40.31 BSTMP%FSHAPE = FSHAPE 40.31 BSTMP%DSHAPE = DSHAPE 40.31 BSTMP%SPPARM(1:4) = SPPARM(1:4) 40.31 NULLIFY(BSTMP%NEXTBS) 40.31 IF ( .NOT.LBS ) THEN 40.31 FBS = BSTMP 40.31 CUBS => FBS 40.31 LBS = .TRUE. 40.31 ELSE 40.31 CUBS%NEXTBS => BSTMP 40.31 CUBS => BSTMP 40.31 END IF 40.31 ELSE IF (LFRST2) THEN CALL INCSTR ('FNAME', FILENM, 'REQ', ' ') LFRST2 = .FALSE. ELSE CALL INCSTR ('FNAME', FILENM, 'STA', ' ') ENDIF IF (FILENM.NE.' ') THEN ! generate new set of file data NBFILS = NBFILS + 1 NBSPSS = NBSPEC ALLOCATE(BFLTMP) 40.31 CALL BCFILE (FILENM, 'PNTS', BFLTMP, LBGP, 40.31 & XCGRID, YCGRID, KGRPNT, XYTST, KGRBND, 40.31 & DONALL) 40.05 IF (STPNOW()) RETURN 34.01 NULLIFY(BFLTMP%NEXTBSPC) 40.31 IF ( .NOT.LBFILS ) THEN 40.31 FBNDFIL = BFLTMP 40.31 CUBFL => FBNDFIL 40.31 LBFILS = .TRUE. 40.31 ELSE 40.31 CUBFL%NEXTBSPC => BFLTMP 40.31 CUBFL => BFLTMP 40.31 END IF 40.31 ENDIF CALL ININTG ('SEQ', NFSEQ, 'STA', 1) IBSPC2 = NBSPSS + NFSEQ IF (IBSPC2.GT.NBSPEC) & CALL MSGERR (1,'too large value for SEQ') ENDIF ELSE IF (IKO.GT.KOUNTR) GOTO 360 40.00 ENDIF LFRST3 = .TRUE. DO IX = CURR%JX 40.31 IY = CURR%JY 40.31 XC2 = XCGRID(IX,IY) YC2 = YCGRID(IX,IY) IF (.NOT.LFRST3) THEN RDIST = RDIST + SQRT ((XC2-XC1)**2 + (YC2-YC1)**2) ENDIF LFRST3 = .FALSE. XC1 = XC2 YC1 = YC2 IF (RDIST.GT.RLEN2) GOTO 340 ALLOCATE(BGPTMP) 40.31 BGPTMP%BGP(1) = KGRPNT(IX,IY) 40.31 BGPTMP%BGP(2) = 1 40.31 W1 = (RLEN2-RDIST)/(RLEN2-RLEN1) BGPTMP%BGP(3) = NINT(1000.*W1) 40.31 BGPTMP%BGP(4) = IBSPC1 40.31 BGPTMP%BGP(5) = NINT(1000.*(1.-W1)) 40.31 BGPTMP%BGP(6) = IBSPC2 40.31 NULLIFY(BGPTMP%NEXTBGP) 40.31 IF ( .NOT.LBGP ) THEN 40.31 FBGP = BGPTMP 40.31 CUBGP => FBGP 40.31 LBGP = .TRUE. 40.31 ELSE 40.31 CUBGP%NEXTBGP => BGPTMP 40.31 CUBGP => BGPTMP 40.31 END IF 40.31 IKO = IKO + 1 IF (IKO.GT.KOUNTR) GOTO 340 40.00 IF (.NOT.ASSOCIATED(CURR%NEXTXY)) EXIT 40.31 CURR => CURR%NEXTXY 40.31 ENDDO ! boundary values have been assigned, read new parameters 340 IF (RLEN2.GT.0.9E20) GOTO 360 40.00 RLEN1 = RLEN2 IBSPC1 = IBSPC2 ENDDO ! update NBGRPT = number of boundary grid points 360 NBGRPT = NBGRPT + KOUNTR 40.00 ELSE CALL WRNKEY ENDIF IF (ASSOCIATED(TMP)) DEALLOCATE(TMP) 40.31 ENDIF 900 RETURN END !********************************************************************* ! * SUBROUTINE BCFILE (FBCNAM, BCTYPE, BSPFIL, LBGP, 40.31 & XCGRID, YCGRID, KGRPNT, 40.31 & XYTST, KGRBND, DONALL) 40.31 40.05 ! * !********************************************************************* ! USE OCPCOMM1 40.41 USE OCPCOMM2 40.41 USE OCPCOMM4 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 USE SWCOMM4 40.41 USE M_BNDSPEC 40.31 ! IMPLICIT NONE ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 30.73: Nico Booij ! 30.90: IJsbrand Haagsma (Equivalence version) ! 34.01: Jeroen Adema ! 40.03, 40.13: Nico Booij ! 40.05: Ekaterini E. Kriezi ! 40.31: Marcel Zijlema ! 40.41: Marcel Zijlema ! ! 1. Updates ! ! 30.73, Dec. 97: New subroutine ! 30.90, Oct. 98: Introduced EQUIVALENCE POOL-arrays ! 34.01, Feb. 99: Introducing STPNOW ! 40.03, June 00: function EQCSTR used to compare strings ! July 00: option LONLAT for location coordinates introduced ! 40.05, Aug. 00: replace the source text related with the grid points ! interpolation coef. with a new subroutine BC_POINTS ! 40.13, Jan. 01: ! is now allowed as comment sign in a boundary file ! checking coordinates only for nesting situation ! remove declarations of unused variables ! Nov. 01: initial size of BSPAUX array enlarged ! 40.31, Nov. 03: removing POOL-mechanism, reconsideration of this ! subroutine ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! 40.41, Nov. 04: small corrections ! ! 2. Purpose ! ! Reads file data for boundary condition ! ! 3. Method ! ! ! 4. Argument variables ! REAL XCGRID(MXC,MYC), YCGRID(MXC,MYC) ! INTEGER KGRPNT(MXC,MYC) 40.31 INTEGER XYTST(*), KGRBND(*) 40.31 ! ! FBCNAM char inp filename of boundary data file ! BCTYPE char inp if value is "NEST": nesting b.c. ! XCGRID real inp x-coordinate of computational grid points ! YCGRID real inp y-coordinate of computational grid points ! KGRPNT int inp indirect addresses of grid points ! XYTST int inp ix, iy of test points ! ! DONALL: logic arguments declare if the nesting boundary is open or close ! it is defined by the users ! LOGICAL, INTENT(INOUT) :: DONALL 40.05 ! CHARACTER FBCNAM *(*), BCTYPE *(*) TYPE(BSPCDAT) :: BSPFIL 40.31 LOGICAL :: LBGP 40.31 ! ! 5. Parameter variables ! ! ! 6. Local variables ! INTEGER :: ISTATF, NDSL, NDSD, IOSTAT, IERR, NBOUNC, NANG, NFRE INTEGER :: IBOUNC, DORDER INTEGER :: IENT,IOPTT INTEGER :: NHEDF, NHEDT, NHEDS, IFRE , IANG INTEGER :: NQUANT, IQUANT, IBC, II, NBGRPT_PREV,IIPT2 REAL :: XP, YP, XP2, YP2 REAL :: FREQHZ, DIRDEG, DIRRD1,DIRRAD, EXCV CHARACTER BTYPE *4, HEDLIN *80 ! ! NBGRPT_PREV is the prevous number of NBGRPT ! IIPT2 counter use for the chekinf if there are grid points on nested boundary ! ! 8. Subroutines Used ! ! Ocean Pack command reading routines ! SWBCPT : boundary points interpolation ! LOGICAL STPNOW, EQCSTR 40.03 ! ! 9. Subroutines calling ! ! SWREAD ! ! 10. Error messages ! ! ! 11. Remarks ! ! ! This subroutine reads the heading of the file to determine locations ! of boundary spectra, spectral frequencies and directions etc. ! Reading and processing of spectral energy densities is done during ! computation by subroutine RESPEC (file Swanmain.for) ! ! data concerning boundary files are stored in array BFILED ! there is a subarray for each file; it contains: ! 1. status; 0: stationary, 1: nonstat, -1: exhausted ! 2. time of boundary values read one before last ! 3. time of boundary values read last ! 4. NDSL: unit ref. num. of file containing filenames ! 5. NDSD: unit ref. num. of file containing data ! 6. time coding option for reading time from data file ! 8. number of locations for which spectra are in the file ! 9. order of reading directional information ! 10. number of spectral directions of spectra on file ! 12. number of spectral frequencies ! 14. number of heading lines per file ! 15. number of heading lines per time step ! 16. number of heading lines per spectrum ! 17. =1: energy dens., =2: variance density ! 18. =1: Cartesian direction, =2: Nautical dir. 40.00 ! 19. =1: direction spread in degr, =2: Power of Cos. 40.00 ! ! ! 12. Structure ! ! ----------------------------------------------------------------- ! Open boundary condition data file ! Read type of file from first line of file ! Case filetype is: ! TPAR: make filetype TPAR ! SWAN: make filetype SWAN ! If b.c. type is NEST ! Then calculate data on grid points ! put into array BGRIDP ! ----------------------------------------------------------- ! Read spectral directions from file into array BSPDIR ! Read spectral frequencies from file into array BSPFRQ ! ----------------------------------------------------------------- ! Put file characteristics into array BFILED ! ----------------------------------------------------------------- ! ! 13. Source text ! LOGICAL CCOORD 40.03 SAVE IENT DATA IENT /0/ CALL STRACE (IENT, 'BCFILE') ! NDSL = 0 IIPT2 = 0 40.05 ! open data file NDSD = 0 IOSTAT = 0 CALL FOR (NDSD, FILENM, 'OF', IOSTAT) IF (STPNOW()) RETURN 34.01 ! ! --- initialize array BFILED of BSPFIL 40.31 BSPFIL%BFILED = 0 40.31 ! ! start reading from the data file READ (NDSD, '(A)') HEDLIN IF (EQCSTR(HEDLIN,'TPAR')) THEN 40.03 BTYPE = 'TPAR' ISTATF = 1 IOPTT = 1 NBOUNC = 1 NANG = 0 NFRE = 0 NHEDF = 0 NHEDT = 0 NHEDS = 0 DORDER = 0 ALLOCATE(BSPFIL%BSPFRQ(NFRE)) 40.41 ALLOCATE(BSPFIL%BSPDIR(NANG)) 40.41 IF (NSTATM.EQ.0) CALL MSGERR (3, & 'time information not allowed in stationary mode') NSTATM = 1 ELSE IF (EQCSTR(HEDLIN,'SWAN')) THEN 40.03 NHEDF = 0 10 READ (NDSD, '(A)') HEDLIN IF (ITEST.GE.60) WRITE (PRTEST,11) HEDLIN 40.00 ! skip heading lines starting with comment sign ($ as in command file) IF (HEDLIN(1:1).EQ.COMID .OR. HEDLIN(1:1).EQ.'!') GOTO 10 40.13 IF (EQCSTR(HEDLIN,'TIME')) THEN 40.03 IF (NSTATM.EQ.0) CALL MSGERR (3, & 'nonstationary boundary condition not allowed '// & 'in stationary mode') NSTATM = 1 ISTATF = 1 BTYPE = 'SWNT' READ (NDSD, *) IOPTT READ (NDSD, '(A)') HEDLIN IF (ITEST.GE.60) WRITE (PRTEST,11) HEDLIN 11 FORMAT (' heading line: ', A) NHEDF = 2 NHEDT = 1 ELSE ISTATF = 0 BTYPE = 'SWNS' NHEDT = 0 ENDIF ! ! read geographical locations ! ! read number of boundary points CCOORD = .TRUE. 40.03 IF (EQCSTR(HEDLIN,'LOC')) THEN 40.03 IF (BCTYPE.EQ.'NEST' .AND. KSPHER.EQ.1) CALL MSGERR (3, 40.13 & 'Boundary locations are Cartesian, while comp. is spherical') 40.03 ELSE IF (EQCSTR(HEDLIN,'LONLAT')) THEN 40.03 IF (BCTYPE.EQ.'NEST' .AND. KSPHER.EQ.0) CALL MSGERR (3, 40.13 & 'Boundary locations are spherical, while comp. is Cartesian') 40.03 ELSE ! set CCOORD to False to indicate that no locations are defined 40.03 CCOORD = .FALSE. 40.03 ENDIF IF (CCOORD) THEN READ (NDSD, *) NBOUNC DO IBOUNC = 1, NBOUNC IERR = 0 CALL REFIXY (NDSD, XP, YP, IERR) IF (ITEST.GE.80) THEN WRITE (PRTEST, *) ' B. spectrum ', IBOUNC, XP+XOFFS, & YP+YOFFS, IERR 40.03 ENDIF ! in case of nesting coordinates on file are used to ! determine interpolation coefficients ! in other cases coordinates are ignored IF (BCTYPE .EQ. 'NEST') THEN XP2 = XP YP2 = YP ! ! --- interpolate the boundaries points to the grid points of ! the SWAN computational grid ! NBGRPT_PREV = NBGRPT 40.05 CALL SWBCPT ( LBGP, XCGRID, YCGRID, 40.41 40.31 40.05 & KGRPNT, XYTST, KGRBND,XP2,YP2,IBOUNC, 40.05 & NBOUNC,DONALL ) 40.05 ! check if the grid points are on nested boundary. ! if not, stop the calculation and give an error message IF (NBGRPT.NE.NBGRPT_PREV) THEN 40.05 IIPT2 = IIPT2+1 40.05 ENDIF 40.05 ENDIF ENDDO ! IF ((BCTYPE .EQ. 'NEST').AND. (IIPT2.EQ.0)) CALL MSGERR (2, & 'no grid points on nested boundary') 40.05 ! NHEDF = NHEDF + 2 + NBOUNC IF (ITEST.GE.60) WRITE (PRTEST,16) NBOUNC 16 FORMAT (I6, ' boundary locations') READ (NDSD, '(A)') HEDLIN IF (ITEST.GE.60) WRITE (PRTEST,11) HEDLIN ELSE IF (BCTYPE .EQ. 'NEST') THEN CALL MSGERR (3, 'this file is not a true nesting file') ENDIF NBOUNC = 1 ENDIF ! ! read spectral resolution information ! ! number of spectral frequencies IF (EQCSTR(HEDLIN(2:5),'FREQ')) THEN 40.03 READ (NDSD, *) NFRE ALLOCATE(BSPFIL%BSPFRQ(NFRE)) 40.31 DO IFRE = 1, NFRE ! read frequency in Hz and convert to radians/sec READ (NDSD, *) FREQHZ BSPFIL%BSPFRQ(IFRE) = PI2 * FREQHZ 40.31 ENDDO READ (NDSD, '(A)') HEDLIN IF (ITEST.GE.60) WRITE (PRTEST,11) HEDLIN 40.00 NHEDF = NHEDF + 2 + NFRE ELSE NFRE = 0 IF (BCTYPE.EQ.'NEST') THEN CALL MSGERR (3, 'file is not a true nesting file') ENDIF ENDIF IF (ITEST.GE.60) WRITE (PRTEST,19) NFRE 40.00 19 FORMAT (I6, ' boundary frequencies') ! number of spectral directions IF (EQCSTR(HEDLIN(2:4),'DIR')) THEN 40.03 READ (NDSD, *) NANG ALLOCATE(BSPFIL%BSPDIR(NANG)) 40.31 DO IANG = 1, NANG ! read direction in degr and convert to radians READ (NDSD, *) DIRDEG IF (EQCSTR(HEDLIN,'N')) THEN 40.03 DIRDEG = 180. + DNORTH - DIRDEG ENDIF DIRRAD = DIRDEG * PI / 180. ! reverse order if second direction is smaller than first IF (IANG.EQ.1) THEN DIRRD1 = DIRRAD ELSE IF (IANG.EQ.2) THEN IF (DIRRAD.LT.DIRRD1) THEN DORDER = -1 BSPFIL%BSPDIR(NANG) = DIRRD1 40.31 ELSE DORDER = 1 ENDIF DIRRD1 = DIRRAD ELSE IF (DORDER.LT.0.) THEN IF (DIRRAD.GT.DIRRD1) CALL MSGERR (3, & 'spectral directions in file not in right order') ELSE IF (DIRRAD.LT.DIRRD1) CALL MSGERR (3, & 'spectral directions in file not in right order') ENDIF DIRRD1 = DIRRAD ENDIF IF (DORDER.LT.0) THEN BSPFIL%BSPDIR(NANG+1-IANG) = DIRRAD 40.31 ELSE BSPFIL%BSPDIR(IANG) = DIRRAD 40.31 ENDIF ENDDO READ (NDSD, '(A)') HEDLIN IF (ITEST.GE.60) WRITE (PRTEST,11) HEDLIN NHEDF = NHEDF + 2 + NANG NHEDS = 1 ELSE NANG = 0 ALLOCATE(BSPFIL%BSPDIR(NANG)) 40.41 NHEDS = 1 DORDER = 0 ENDIF IF (ITEST.GE.60) WRITE (PRTEST,23) NANG 40.00 23 FORMAT (I6, ' boundary directions') ! ! read quantities (name, unit, exc. value) ! IF (EQCSTR(HEDLIN,'QUANT')) THEN READ (NDSD, *) NQUANT IF (.NOT.((NQUANT.EQ.1 .AND. NANG.GT.0) .OR. & (NQUANT.EQ.3 .AND. NANG.EQ.0))) THEN CALL MSGERR (2, 'incompatible data on b.c. file') WRITE (PRINTF, 31) NQUANT, NANG 31 FORMAT (I3, ' quantities; ', I5, ' directions') ENDIF DO IQUANT = 1, NQUANT READ (NDSD, '(A)') HEDLIN ! if first quantity is 'EnDens' divide by Rho*Grav IF (IQUANT.EQ.1) THEN IF ( EQCSTR(HEDLIN,'ENDENS')) THEN 40.03 ! quantity on file is energy density BSPFIL%BFILED(17) = 1 40.31 ELSE IF ( EQCSTR(HEDLIN,'VADENS')) THEN 40.03 ! quantity on file is variance density BSPFIL%BFILED(17) = 2 40.31 ELSE CALL MSGERR (2, & 'instead of EN/VADENS b.c.file has: ' // HEDLIN(1:10)) 40.03 BSPFIL%BFILED(17) = 2 40.31 ENDIF ELSE IF (IQUANT.EQ.2) THEN 40.00 ! if second quantity is 'NDIR' transform from Nautical to Cartesian dir. IF ( EQCSTR(HEDLIN,'NDIR')) THEN 40.03 ! quantity on file is Nautical direction BSPFIL%BFILED(18) = 2 40.31 ELSE IF (EQCSTR(HEDLIN,'CDIR')) THEN 40.03 ! quantity on file is Cartesian direction BSPFIL%BFILED(18) = 1 40.31 ELSE CALL MSGERR (2, & 'instead of NDIR,CDIR b.c.file has: ' // HEDLIN(1:10)) 40.03 BSPFIL%BFILED(18) = 1 40.31 ENDIF ELSE IF (IQUANT.EQ.3) THEN 40.00 ! if third quantity is 'DSPRP' or 'POWER' power is given, ! otherwise calculate power from dir. spread in degrees IF (EQCSTR(HEDLIN,'DSPRP') .OR. & EQCSTR(HEDLIN,'POWER')) THEN 40.03 ! quantity on file is power of cos BSPFIL%BFILED(19) = 2 40.31 ELSE IF (EQCSTR(HEDLIN,'DSPR') .OR. & EQCSTR(HEDLIN,'DEGR')) THEN 40.03 ! quantity on file is Directional spread in degr BSPFIL%BFILED(19) = 1 40.31 ELSE CALL MSGERR (2, & 'instead of DSPR,DEGR b.c.file has: ' // HEDLIN(1:10)) 40.03 BSPFIL%BFILED(19) = 1 40.31 ENDIF ENDIF ! check Unit and Exception value READ (NDSD, '(A)') HEDLIN IF (IQUANT.EQ.3 .AND. EQCSTR(HEDLIN,'DEGR')) THEN IF (BSPFIL%BFILED(19).NE.1) THEN 40.31 CALL MSGERR (2, 'incompatible options in boundary file') BSPFIL%BFILED(19) = 1 40.31 ENDIF ENDIF IF (IQUANT.EQ.1) THEN 40.41 READ (NDSD, *) EXCV BSPFIL%BFILED(11) = NINT(EXCV) ELSE READ (NDSD, '(A)') HEDLIN END IF ENDDO NHEDF = NHEDF + 2 + 3*NQUANT ENDIF IF (ITEST.GE.60) WRITE (PRTEST,28) NQUANT 40.00 28 FORMAT (I6, ' quantities') ELSE CALL MSGERR (3, 'unsupported boundary data file') ENDIF ! ALLOCATE(BSPFIL%BSPLOC(NBOUNC)) 40.31 DO IBC = 1, NBOUNC BSPFIL%BSPLOC(IBC) = NBSPEC + IBC 40.31 ENDDO NBSPEC = NBSPEC + NBOUNC ! ! store file reading parameters in array BFILED ! BSPFIL%BFILED(1) = ISTATF 40.31 BSPFIL%BFILED(2) = -999999999 40.31 BSPFIL%BFILED(3) = -999999999 40.31 BSPFIL%BFILED(4) = NDSL 40.31 BSPFIL%BFILED(5) = NDSD 40.31 BSPFIL%BFILED(6) = IOPTT 40.31 CALL COPYCH (BTYPE, 'T', BSPFIL%BFILED(7), 1, IERR) 40.31 BSPFIL%BFILED(8) = NBOUNC 40.31 BSPFIL%BFILED(9) = DORDER 40.31 BSPFIL%BFILED(10) = NANG 40.31 BSPFIL%BFILED(12) = NFRE 40.31 ! ordering of data on file BSPFIL%BFILED(13) = 0 40.31 ! number of heading lines: per file, per time, per spectrum BSPFIL%BFILED(14) = NHEDF 40.31 BSPFIL%BFILED(15) = NHEDT 40.31 BSPFIL%BFILED(16) = NHEDS 40.31 ! IF (ITEST.GE.80) WRITE(PRINTF,81) NBFILS, NBSPEC, & (BSPFIL%BFILED(II), II=1,16) 40.31 81 FORMAT (' array BFILED: ', 2I4, 2(/,8I10)) ! RETURN ! end of subroutine BCFILE END !********************************************************************* ! * SUBROUTINE BCWAMN (FBCNAM, BCTYPE, BSPFIL, LBGP, 40.31 & XCGRID, YCGRID, KGRPNT, XYTST) 40.31 ! * !********************************************************************* USE OCPCOMM2 40.41 USE OCPCOMM4 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 USE SWCOMM4 40.41 USE M_BNDSPEC 40.31 USE M_PARALL IMPLICIT NONE 40.13 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 30.73: Nico Booij ! 30.90: IJsbrand Haagsma (Equivalence version) ! 34.01: Jeroen Adema ! 40.03: Nico Booij ! 40.13: N. Booij ! 40.31: Marcel Zijlema ! 40.41: Marcel Zijlema ! ! 1. Updates ! ! 30.73, Jan. 98: new subroutine, based on older version by Weimin Luo ! 30.90, Oct. 98: Introduced EQUIVALENCE POOL-arrays ! 34.01, Feb. 99: Introducing STPNOW ! 40.03, Nov. 99: THD (first spectral direction in radians) added in ! expression for RBSDIR (directions of boundary spectrum) ! 40.03, Aug. 00: correction WAM nest with spherical SWAN ! 40.13, May 01: order of boundary points in WAM nesting file differed ! from order assumed in SWAN ! 40.31, Nov. 03: removing POOL-mechanism, reconsideration of this ! subroutine ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. PURPOSE ! ! reads file data for WAM nesting boundary condition ! ! 3. METHOD ! ! ! 4. Argument variables ! ! FBCNAM char inp filename of boundary data file ! BCTYPE char inp if value is "NEST": nesting b.c. ! XCGRID real inp x-coordinate of computational grid points ! YCGRID real inp y-coordinate of computational grid points ! KGRPNT int inp indirect addresses of grid points ! XYTST int inp ix, iy of test points ! ! 7. Common blocks used ! !OMP INCLUDE 'tpcomm.inc' 40.41 ! ! 5. SUBROUTINES CALLING ! ! SWBOUN ! ! 6. SUBROUTINES USED ! ! Ocean Pack command reading routines ! LOGICAL :: STPNOW 34.01 ! 7. ERROR MESSAGES ! ! --- ! ! 8. REMARKS ! ! data concerning boundary files are stored in array BFILED ! there is a subarray for each file; it contains: ! 1. status; 0: stationary, 1: nonstat, -1: exhausted ! 2. time of boundary values read one before last ! 3. time of boundary values read last ! 4. NDSL: unit ref. num. of file containing filenames ! 5. NDSD: unit ref. num. of file containing data ! 6. time coding option for reading time from data file ! 8. number of locations for which spectra are in the file ! 9. order of reading directional information ! 10. number of spectral directions of spectra on file ! 12. number of spectral frequencies ! 14. number of heading lines per file ! 15. number of heading lines per time step ! 16. number of heading lines per spectrum ! 17. =1: energy dens., =2: variance density ! ! ! 9. STRUCTURE ! ! ----------------------------------------------------------------- ! Open file containing filnames ! Read data file name ! Open boundary condition data file ! Read number of b.points, frequencies and directions ! Generate spectral directions from file into array BSPDIR ! Generate spectral frequencies from file into array BSPFRQ ! For all boundary spectra do ! read location from data file ! transform into local cartesian or spherical coordinates 40.13 ! ----------------------------------------------------------------- ! Determine spatial step size in WAM nesting file 40.13 ! For all spatial points in WAM file do 40.13 ! For all other spatial points in WAM file do 40.13 ! If the two points are neighbours 40.13 ! Then For all computational grid points on boundary do ! if point is located between nest file grid points ! calculate interpolation coefficients ! and put these into array BGRIDP ! ----------------------------------------------------------------- ! Write file characteristics into array BFILED ! ----------------------------------------------------------------- ! ! 10. SOURCE TEXT ! INTEGER KGRPNT(MXC,MYC), XYTST(*) REAL XCGRID(MXC,MYC), YCGRID(MXC,MYC) TYPE(BSPCDAT) :: BSPFIL 40.31 LOGICAL :: LBGP 40.31 CHARACTER FBCNAM *(*), BCTYPE *(*) ! ! local variables ! INTEGER ISTATF, NDSL, NDSD, IOSTAT, IERR, NBOUNC, NANG, NFRE, & IBOUNC, IX1, IY1, IX2, IY2, IXP, IYP, IP, MIP, INDXGR, & DORDER, IOPTT ! ISTATF if >0 file contains nonstationary data ! NDSL unit ref num of namelist file ! NDSD unit ref num of data file ! IOSTAT io status ! IERR error status ! NBOUNC number of boundary locations ! NANG number of directions on file ! NFRE number of frequencies on file ! IBOUNC counter of boundary spectra ! IX1 ! IY1 ! IX2 ! IY2 ! IXP ! IYP ! IP ! MIP ! INDXGR counter of boundary grid points ! DORDER if <0 order of reading directions is reversed ! IOPTT time reading option INTEGER :: IIPT1=0, IIPT2=0 ! local and overall number of interpolated boundary grid points INTEGER :: NHEDF ! number of heading lines at begin of file 40.13 INTEGER :: NHEDT ! number of heading lines per time step 40.13 INTEGER :: NHEDS ! number of heading lines 40.13 INTEGER :: IBC, IDW, ISW, IFRE, II, ISIDE, IHD ! counters INTEGER :: IBNC1, IBNC2 ! counters of nesting points INTEGER :: IBSP1, IBSP2 ! counters of nesting points REAL XP, YP, XP1, YP1, XP2, YP2, RR, RX, RY, RL2, & XANG, XFRE, THD, FR1, CO, XBOU, XDELC, & XLON, XLAT, XDATE, EMEAN, THQ, FMEAN, USNEW, THWNEW ! XP problem coordinate of a comp. grid point on the boundary ! YP problem coordinate of a comp. grid point on the boundary ! XP1 problem coordinate of a boundary location ! YP1 problem coordinate of a boundary location ! XP2 problem coordinate of a boundary location ! YP2 problem coordinate of a boundary location ! RR ! RX vector connecting two boundary locations ! RY vector connecting two boundary locations ! RL2 length **2 of vector connecting two boundary locations DOUBLE PRECISION, ALLOCATABLE :: XPWAM(:), YPWAM(:) REAL, ALLOCATABLE :: SPAUX(:) 40.31 ! locations of nesting points 40.13 DOUBLE PRECISION :: DXWAM, DYWAM ! spatial step sizes in nesting file 40.13 DOUBLE PRECISION :: DXTEST, DYTEST ! distance between two nesting points 40.13 REAL :: DISXY ! dim.less distance 40.13 REAL :: PHI ! direction of vector (RX,RY) 40.13 REAL :: DPHI ! difference in direction 40.13 REAL :: EPS ! tolerance 40.13 REAL :: W2 ! interpolation coefficient 40.13 CHARACTER (LEN=4) :: BTYPE ! type of boundary cond. CHARACTER (LEN=12) :: CDATE ! date-time CHARACTER (LEN=80) :: HEDLIN ! heading line DOUBLE PRECISION :: DDATE, XLON0, XLAT0 ! DDATE date-time ! XLON0 longitude of origin of computational grid ! XLAT0 latitude of origin of computational grid TYPE(BGPDAT), POINTER :: BGPTMP 40.31 ! subroutines used LOGICAL :: KEYWIS INTEGER, SAVE :: IENT = 0 40.13 CALL STRACE (IENT, 'BCWAMN') ! ISTATF = 1 IOPTT = 6 NDSL = 0 ! open file with list of names CALL FOR (NDSL, FILENM,'OF',IOSTAT) IF (STPNOW()) RETURN 34.01 READ (NDSL,'(A36)') FILENM CALL INKEYW ('REQ', ' ') IF (KEYWIS('FRE')) THEN BTYPE = 'WAMF' ELSE IF (KEYWIS('UNF')) THEN CALL INKEYW ('REQ', ' ') IF (KEYWIS('WK')) THEN BTYPE = 'WAMW' ELSE CALL IGNORE ('CRAY') BTYPE = 'WAMC' ENDIF ENDIF ! open WAM data file NDSD=0 IOSTAT = 0 IF (BTYPE.EQ.'WAMF') THEN CALL FOR(NDSD,FILENM,'OF',IOSTAT) IF (STPNOW()) RETURN 34.01 ELSE CALL FOR(NDSD,FILENM,'OU',IOSTAT) IF (STPNOW()) RETURN 34.01 ENDIF ! ! --- initialize array BFILED of BSPFIL 40.31 BSPFIL%BFILED = 0 40.31 ! ! read spherical coordinates of point corresponding to ! [xpc], [ypc] in Cartesian coordinates ! not necessary if SWAN uses spherical coordinates 40.03 ! if not given first point in data file is assumed CALL INDBLE('XGC',XLON0,'STA',-999.D0) 40.01 CALL INDBLE('YGC',XLAT0,'STA',-999.D0) 40.01 ! ! start reading from the data file ! ! read resolution information from WAM input ! IF (BTYPE.EQ.'WAMF') THEN READ (NDSD,*) XANG, XFRE, THD, FR1, CO, XBOU, XDELC ELSE ! Cray and workstation version READ (NDSD) XANG, XFRE, THD, FR1, CO, XBOU, XDELC ENDIF ! ! number of WAM boundary points NBOUNC = NINT(XBOU) ! number of direction of WAM spectrum NANG = NINT(XANG) DORDER = -1 ! number of frequencies of WAM spectrum NFRE = NINT(XFRE) ! number of heading lines: per file, per time, per spectrum NHEDF = 1 NHEDT = 0 NHEDS = 1 IF (ITEST.GE.80) THEN WRITE(PRINTF,*) ' Number of frequencies in WAM:',NFRE WRITE(PRINTF,*) ' Number of directions in WAM:',NANG WRITE(PRINTF,*) ' Lowest frequency in WAM:',FR1 WRITE(PRINTF,*) ' fi/fi-1 in WAM:',CO ENDIF ! ! convert WAM wave directions to SWAN convention ! apparently WAM uses direction TO which waves propagate !! ! ALLOCATE(BSPFIL%BSPDIR(NANG)) 40.31 DO IDW = NANG,1,-1 BSPFIL%BSPDIR(NANG-IDW+1) = DNORTH*DEGRAD - 40.31 30.90 & THD - REAL(IDW-1)*PI2/REAL(NANG) 40.03 ENDDO IF (ITEST.GE.50) WRITE (PRTEST,132) NANG, & (BSPFIL%BSPDIR(IDW)*180./PI, IDW=1,NANG) 40.31 30.90 132 FORMAT (' WAMNEST dirs ', I3, (/, 20F6.0)) ! ! calculate WAM angular frequency array ! ALLOCATE(BSPFIL%BSPFRQ(NFRE)) 40.31 BSPFIL%BSPFRQ(1) = PI2*FR1 40.31 30.90 DO ISW = 2, NFRE BSPFIL%BSPFRQ(ISW) = CO * BSPFIL%BSPFRQ(ISW-1) 40.31 30.90 ENDDO IF (ITEST.GE.50) WRITE (PRTEST,133) NFRE, & (BSPFIL%BSPFRQ(ISW)*180./PI, ISW=1,NFRE) 40.31 30.90 133 FORMAT (' WAMNEST freqs ', I3, (/, 20F6.2)) IF (NBOUNC.EQ.1) CALL MSGERR (3, & 'WAM nest does not work with only one nesting point') 40.13 ! allocate arrays XPWAM and YPWAM ALLOCATE (XPWAM(1:NBOUNC), YPWAM(1:NBOUNC)) 40.13 ALLOCATE (SPAUX(NANG*NFRE)) 40.31 ! read geographical locations and determine DXWAM and DYWAM 40.13 IIPT2 = 0 40.03 DXWAM = 180. DYWAM = 180. DO IBOUNC = 1, NBOUNC IF (BTYPE.EQ.'WAMF') THEN ! read boundary point coordinates from file READ(NDSD,*) XLON, XLAT, DDATE, EMEAN, & THQ, FMEAN, USNEW, THWNEW IF (IBOUNC.EQ.1 .AND. ITEST.GE.80) WRITE (PRTEST, *) 40.13 & ' WAMNEST starting time ', DDATE 40.13 ! read spectral densities but ignore them for the moment DO IFRE=1,NFRE READ(NDSD,*) (SPAUX(II), II=1,NANG) 40.31 30.90 ENDDO ELSE IF (BTYPE.EQ.'WAMC') THEN ! read boundary point coordinates from file READ(NDSD) XLON, XLAT, XDATE, EMEAN, & THQ, FMEAN, USNEW, THWNEW IF (IBOUNC.EQ.1 .AND. ITEST.GE.80) WRITE (PRTEST, *) 40.13 & ' WAMNEST starting time ', XDATE 40.13 ! read spectral densities but ignore them for the moment READ(NDSD) (SPAUX(II), II=1,NANG*NFRE) 40.31 30.90 ELSE ! read boundary point coordinates from file READ(NDSD) XLON, XLAT, CDATE, EMEAN, THQ, FMEAN, USNEW, THWNEW IF (IBOUNC.EQ.1 .AND. ITEST.GE.80) WRITE (PRTEST, *) 40.13 & ' WAMNEST starting time ', CDATE 40.13 ! read spectral densities but ignore them for the moment READ(NDSD) (SPAUX(II), II=1,NANG*NFRE) 40.31 30.90 ENDIF IF (ITEST.GE.50) WRITE (PRINTF, 178) IBOUNC, XLON, XLAT 40.13 178 FORMAT (' boundary spectrum ', I3, ' at ', 2F12.4) 40.13 XPWAM(IBOUNC) = XLON 40.13 YPWAM(IBOUNC) = XLAT 40.13 ! determine DXWAM and DYWAM 40.13 IF (IBOUNC.GT.1) THEN 40.13 IF (ABS(XPWAM(IBOUNC)-XPWAM(IBOUNC-1)).GT.1.E-6) 40.13 & DXWAM = MIN (DXWAM, ABS(XPWAM(IBOUNC)-XPWAM(IBOUNC-1))) 40.13 IF (ABS(YPWAM(IBOUNC)-YPWAM(IBOUNC-1)).GT.1.E-6) 40.13 & DYWAM = MIN (DYWAM, ABS(YPWAM(IBOUNC)-YPWAM(IBOUNC-1))) 40.13 ENDIF IF (KSPHER.EQ.0) THEN 33.09 ! determine lower left corner of WAM nesting grid if not given by the user IF (IBOUNC.EQ.1) THEN IF (XLON0.LT.-900.) THEN XLON0 = XLON XLAT0 = XLAT ENDIF ENDIF ENDIF ENDDO 40.13 IF (ITEST.GE.50) WRITE (PRINTF, 182) DXWAM, DYWAM 40.13 182 FORMAT (' WAM step sizes: ', 2F12.4) 40.13 EPS = 0.01 * MIN(DXWAM,DYWAM) 40.13 ! determine interpolation coefficients for all couples of 40.13 ! neighbouring WAM nest points 40.13 DO IBNC1 = 1, NBOUNC 160 DO IBNC2 = IBNC1+1, NBOUNC 40.13 DXTEST = ABS(XPWAM(IBNC1)-XPWAM(IBNC2)) 40.13 DYTEST = ABS(YPWAM(IBNC1)-YPWAM(IBNC2)) 40.13 IF ((DXTEST.LT.EPS .AND. ABS(DYTEST-DYWAM).LT.EPS) .OR. 40.13 & (DYTEST.LT.EPS .AND. ABS(DXTEST-DXWAM).LT.EPS)) THEN 40.13 ! points IBNC1 and IBNC2 are neighbours 40.13 IF (KSPHER.EQ.0) THEN 33.09 ! transform to local Cartesian coordinates XP1 = XPC + LENDEG*COS(PI*XLAT0/180.)*(XPWAM(IBNC1)-XLON0) 40.13 YP1 = YPC + LENDEG*(YPWAM(IBNC1)-XLAT0) 40.13 XP2 = XPC + LENDEG*COS(PI*XLAT0/180.)*(XPWAM(IBNC2)-XLON0) 40.13 YP2 = YPC + LENDEG*(YPWAM(IBNC2)-XLAT0) 40.13 ELSE XP1 = XPWAM(IBNC1) - XOFFS 40.13 YP1 = YPWAM(IBNC1) - YOFFS 40.13 XP2 = XPWAM(IBNC2) - XOFFS 40.13 YP2 = YPWAM(IBNC2) - YOFFS 40.13 ENDIF ! Determine interpolation coefficients IBSP1 = NBSPEC+IBNC1 40.13 IBSP2 = NBSPEC+IBNC2 40.13 IIPT1 = 0 40.03 RX = XP2 - XP1 RY = YP2 - YP1 RL2 = RX**2 + RY**2 IF (RL2.GT.0.) THEN RX = RX/RL2 RY = RY/RL2 ! check whether direction of (RX,RY) corresponds to ALPC + k * 90 degr PHI = ATAN2(RY,RX) DPHI = MOD(PHI-ALPC+1.25*PI,0.5*PI)-0.25*PI IF (ABS(DPHI) .LT. 0.1) THEN ! loop over boundary of comp. grid, select points between ! (XP1,YP1) and (XP2,YP2) DO ISIDE = 1, 4 IF (ISIDE.EQ.1) THEN IX1 = 1 IY1 = 1 IX2 = MXC IY2 = 1 MIP = MXC ELSE IF (ISIDE.EQ.2) THEN IX1 = MXC IY1 = 1 IX2 = MXC IY2 = MYC MIP = MYC ELSE IF (ISIDE.EQ.3) THEN IX1 = MXC IY1 = MYC IX2 = 1 IY2 = MYC MIP = MXC ELSE IF (ISIDE.EQ.4) THEN IX1 = 1 IY1 = MYC IX2 = 1 IY2 = 1 MIP = MYC ENDIF DO IP = 1, MIP-1 RR = REAL(IP-1) / REAL(MIP-1) IXP = IX1 + NINT(RR*REAL(IX2-IX1)) IYP = IY1 + NINT(RR*REAL(IY2-IY1)) INDXGR = KGRPNT(IXP,IYP) IF (INDXGR.GT.1) THEN XP = XCGRID(IXP,IYP) YP = YCGRID(IXP,IYP) ! DISXY is relative distance from (XP,YP) to line ! (XP1,YP1) to (XP2,YP2) DISXY = ABS(RX*(YP-YP1)-RY*(XP-XP1)) IF (DISXY.LT.0.1) THEN ! W2 is relative length of projection on line ! (XP1,YP1) to (XP2,YP2) W2 = RX*(XP-XP1)+RY*(YP-YP1) IF (W2.GT.-0.001 .AND. W2.LT.1.001) THEN IF (W2.LT.0.01) W2 = 0. IF (W2.GT.0.99) W2 = 1. IF (ITEST.GE.80) WRITE (PRTEST, 223) IXP, IYP, & XP+XOFFS, YP+YOFFS, W2, IBNC1, IBNC2 40.13 223 FORMAT (' B.pnt', 2I5, 2F14.4, F6.3, & ' from ', 2I3) 40.13 NBGRPT = NBGRPT + 1 IIPT1 = IIPT1 + 1 40.03 IIPT2 = IIPT2 + 1 40.03 ALLOCATE(BGPTMP) 40.31 BGPTMP%BGP(1) = INDXGR 40.31 ! next item indicates type of boundary condition BGPTMP%BGP(2) = 1 40.31 BGPTMP%BGP(3) = NINT(1000. * W2) 40.31 BGPTMP%BGP(4) = IBSP2 40.31 BGPTMP%BGP(5) = NINT(1000. * (1.-W2)) 40.31 BGPTMP%BGP(6) = IBSP1 40.31 NULLIFY(BGPTMP%NEXTBGP) 40.31 IF ( .NOT.LBGP ) THEN 40.31 FBGP = BGPTMP 40.31 CUBGP => FBGP 40.31 LBGP = .TRUE. 40.31 ELSE 40.31 CUBGP%NEXTBGP => BGPTMP 40.31 CUBGP => BGPTMP 40.31 END IF 40.31 ! test output if point is a test point IF (NPTST.GT.0) THEN DO IPTST = 1, NPTST IF (IXP.EQ.XYTST(2*IPTST-1)+MXF-1 .AND. & IYP.EQ.XYTST(2*IPTST )+MYF-1) & WRITE (PRTEST, 223) IXP, IYP, & XP+XOFFS, YP+YOFFS, W2, IBSP2, IBSP1 ENDDO ENDIF ENDIF ENDIF ENDIF ENDDO ENDDO ENDIF ENDIF IF (IIPT1.EQ.0) THEN WRITE (PRINTF, 218) XP1+XOFFS, YP1+YOFFS, & XP2+XOFFS, YP2+YOFFS 218 FORMAT (' Warning: no grid points on interval from ', 40.03 & 2F14.4, ' to ', 2F14.4) ENDIF ENDIF ENDDO 40.13 ! first nesting point may have two valid neighbours 40.13 IF (IBNC1.EQ.1 .AND. IBNC2.EQ.2) GOTO 160 40.13 ENDDO IF (IIPT2.EQ.0) CALL MSGERR (2, & 'no grid points on nested boundary') 40.03 IF (ITEST.GE.60) WRITE (PRTEST,16) NBOUNC 16 FORMAT (I6, ' boundary locations') ! deallocate arrays XPWAM, YPWAM and SPAUX DEALLOCATE (XPWAM, YPWAM, SPAUX) 40.31 40.13 ALLOCATE(BSPFIL%BSPLOC(NBOUNC)) 40.31 DO IBC = 1, NBOUNC BSPFIL%BSPLOC(IBC) = NBSPEC + IBC 40.31 ENDDO NBSPEC = NBSPEC + NBOUNC ! ! store file reading parameters in array BFILED ! BSPFIL%BFILED(1) = ISTATF 40.31 BSPFIL%BFILED(2) = -999999999 40.31 BSPFIL%BFILED(3) = -999999999 40.31 BSPFIL%BFILED(4) = NDSL 40.31 BSPFIL%BFILED(5) = NDSD 40.31 BSPFIL%BFILED(6) = IOPTT 40.31 CALL COPYCH (BTYPE, 'T', BSPFIL%BFILED(7), 1, IERR) 40.31 BSPFIL%BFILED(8) = NBOUNC 40.31 BSPFIL%BFILED(9) = DORDER 40.31 BSPFIL%BFILED(10) = NANG 40.31 BSPFIL%BFILED(11) = 0 40.31 BSPFIL%BFILED(12) = NFRE 40.31 ! ordering of data on file BSPFIL%BFILED(13) = 0 40.31 ! number of heading lines: per file, per time, per spectrum BSPFIL%BFILED(14) = NHEDF 40.31 BSPFIL%BFILED(15) = NHEDT 40.31 BSPFIL%BFILED(16) = NHEDS 40.31 ! quantity on file is variance density BSPFIL%BFILED(17) = 2 40.31 ! IF (ITEST.GE.80) WRITE(PRINTF,81) NBFILS, NBSPEC, & (BSPFIL%BFILED(II), II=1,16) 40.31 81 FORMAT (' array BFILED: ', 2I4, 2(/,8I10)) ! ! Rewind input file for proper start REWIND (NDSD) ! read heading line IF (BTYPE.EQ.'WAMF') THEN DO IHD = 1, BSPFIL%BFILED(14) 40.31 READ (NDSD, '(A)') HEDLIN IF (ITEST.GE.80) WRITE (PRINTF, 212) HEDLIN 212 FORMAT (' heading line: ', A) ENDDO ELSE DO IHD = 1, BSPFIL%BFILED(14) 40.31 READ (NDSD) ENDDO ENDIF RETURN END SUBROUTINE BCWAMN !********************************************************************* ! * SUBROUTINE BCWW3N (FBCNAM, BCTYPE, BSPFIL, LBGP, 40.31 & XCGRID, YCGRID, KGRPNT, 40.31 & XYTST, KGRBND, DONALL) 40.31 ! * !********************************************************************* ! USE OCPCOMM2 40.41 USE OCPCOMM4 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 USE SWCOMM4 40.41 USE M_BNDSPEC 40.31 ! IMPLICIT NONE ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 40.05 : Ekaterini E. Kriezi ! 40.13 : Nico Booij ! 40.31 : Marcel Zijlema ! 40.41 : Marcel Zijlema ! ! 1. Updates ! ! 40.05, Aug. 00: new subroutine ! 40.13, Jan. 01: remove declarations of unused variables ! 40.31, Nov. 03: removing POOL-mechanism, reconsideration this ! subroutine ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. Purpose ! ! reads file data for WaveWatch III boundary condition ! ! 3. Method ! ! open boundaries files ! read from ASSCII files the points where the energy desity is given and ! interpolate them to the grid points of the SWAN computational grid ! ! 4. Argument variables ! INTEGER, INTENT(IN) :: KGRPNT(MXC,MYC) ! indirect addresses of computational grid points INTEGER, INTENT(IN) :: KGRBND(*) ! array of boundary grid points INTEGER, INTENT(IN) :: XYTST(*) ! array of (ix,iy) of test points REAL, INTENT(IN) :: XCGRID(MXC,MYC), YCGRID(MXC,MYC) ! coordinates of computational grid points ! FBCNAM char inp filename of boundary data file ! BCTYPE char inp boundary condition type, is 'WW3N' in this case ! CHARACTER FBCNAM *(*), BCTYPE *(*) ! ! DONALL : logic arguments declare if the boundary is open or close ! LOGICAL :: DONALL TYPE(BSPCDAT) :: BSPFIL 40.31 LOGICAL :: LBGP 40.31 ! ! 5. Parameter variables ! ! -- ! ! 6. Local variables ! ! IENT number of entries into this subroutine ! IBC spectrum counter ! INTEGER :: IENT, IBC ! REAL, ALLOCATABLE :: FRQ_ARRAY(:), DIR_ARRAY(:) ! ! IHD counter of heading lines ! WWDATE date in boundary file ! WWTIME time in boundary file ! INTEGER :: WWDATE, WWTIME ! ! ISTATF if >0 file contains nonstationary data ! NDSL unit ref num of namelist file ! NDSD unit ref num of data file ! IOSTAT io status ! IERR error status ! NBOUNC number of boundary locations ! NANG number of directions on file ! NFRE number of frequencies on file ! NBPT number of WW3 boundary points ! DORDER if <0 order of reading directions is reversed ! IOPTT time reading option ! IBOUNC counter of boundary points ! IGRBND counter of boundary grid(swan grid) points ! II counter ! NHEDF number of heading lines per file ! NHEDT number of heading lines per time step ! NHEDS number of heading lines per spectrum ! INTEGER :: ISTATF, NDSL, NDSD, IOSTAT, IERR INTEGER :: NBOUNC, NANG, NFRE INTEGER :: IBOUNC INTEGER :: DORDER, IOPTT INTEGER :: NHEDF, NHEDT, NHEDS, NBGRPT_PREV INTEGER :: IFRE, IANG,II, IIPT2 ! ! DUM_A real number used for reading a file but not used in any calculation ! XLON longitude ! XLAT latitude ! XP2 problem coordinate of a boundary location ! YP2 problem coordinate of a boundary location ! DIRRD1 ! NBGRPT_PREV is the prevous number of NBGRPT ! IIPT2 counter use for the chekinf if there are grid points on nested boundary ! REAL :: DUM_A, XLON, XLAT,XP2,YP2,DIRRD1 ! CHARACTER (LEN=4) :: BTYPE ! type of boundary cond. CHARACTER (LEN=24) :: HEDLINT ! WW3 version CHARACTER (LEN=30) :: GNAME ! name of test case readed from b. file CHARACTER (LEN=12) :: PTNME ! name of b. point ! XLON0 longitude of origin of computational grid ! XLAT0 latitude of origin of computational grid ! DOUBLE PRECISION :: XLON0, XLAT0 ! ! FBCNAM char inp filename of boundary data file ! BCTYPE char inp if value is "NEST": nesting b.c. ! XCGRID real inp x-coordinate of computational grid points ! YCGRID real inp y-coordinate of computational grid points ! KGRPNT int inp indirect addresses of grid points ! XYTST int inp ix, iy of test points ! ! 8. Subroutines used ! ! Ocean Pack command reading routines ! SWBCPT, STPNOW ! LOGICAL :: STPNOW ! ! 9. Subroutines calling ! ! SWREAD ! ! 10. Error messages ! ! --- ! ! 11. Remarks ! ! data concerning boundary files are stored in array BFILED ! there is a subarray for each file; it contains: ! 1. status; 0: stationary, 1: nonstat, -1: exhausted ! 2. time of boundary values read one before last ! 3. time of boundary values read last ! 4. NDSL: unit ref. num. of file containing filenames ! 5. NDSD: unit ref. num. of file containing data ! 6. time coding option for reading time from data file ! 8. number of locations for which spectra are in the file ! 9. order of reading directional information ! 10. number of spectral directions of spectra on file ! 12. number of spectral frequencies ! 14. number of heading lines per file ! 15. number of heading lines per time step ! 16. number of heading lines per spectrum ! 17. =1: energy dens., =2: variance density, =3 variance energy density (k) ! 18. =1: Cartesian direction, =2: Nautical dir. ! 19. =1: direction spread in degr, =2: Power of Cos. ! 20. depth of boundary points ! ! 12. Structure ! ! ----------------------------------------------------------------- ! Open boundary condition data file ! Read type of file from first line of file ! If the headline is WAVEWATCH III SPECTRA ! ! then b.c. type is WW3N ! ----------------------------------------------------------- ! Read spectral directions from file and write them into ! array BSPDIR ! Read spectral frequencies from fileand write them into ! array BSPFRQ ! For all boubdaries points do ! read location from data file ! transform into local cartesian coordinates (if nesesery) ! Then calculate data on grid points ! ! ----------------------------------------------------------------- ! Put file characteristics into array BFILED ! ----------------------------------------------------------------- ! ! 13. Source text ! SAVE IENT DATA IENT /0/ CALL STRACE (IENT, 'BCWW3N') ! ! NDSL unit ref number for namelist files NDSL = 0 ISTATF = 1 ! number of heading lines: per file, per time, per spectrum NHEDF = 0 NHEDT = 0 NHEDS = 0 DORDER = -1 IOPTT = 1 IIPT2 = 0 ! ! open data file NDSD unit ref number for data files NDSD = 0 IOSTAT = 0 ! CALL FOR (NDSD, FILENM , 'OF', IOSTAT) IF (STPNOW()) RETURN ! ! --- initialize array BFILED of BSPFIL 40.31 BSPFIL%BFILED = 0 40.31 ! ! start reading from the boundary data file ! HEDLINT = 'WAVEWATCH III SPECTRA' = header of the file ! read from boundary file the header , number of frequencies NFR, ! number of direction NANG, number of boundaries points NBOUNC, ! Name of the points GNAME ! READ (NDSD, 1944) HEDLINT,NFRE, NANG, NBOUNC, GNAME ! IF (HEDLINT(2:22) .NE. 'WAVEWATCH III SPECTRA') & CALL MSGERR (3, 'file is not a WW3 spectral file') IF (NBOUNC.LT.2) CALL MSGERR & (3, 'SWAN need at least 2 boundary points for nesting') ! ! ISTATF related with stationary and non stationary mode ! BTYPE = 'WW3N' ! ! read frequencies from WW3 boundary file ! ALLOCATE (FRQ_ARRAY(1:NFRE), DIR_ARRAY(1:NANG) ) ! ! read frequency ! FRQ_ARRAY(IFRE) = SIG(IK)/(2*PI) ! READ (NDSD,1945) (FRQ_ARRAY(IFRE) ,IFRE=1,NFRE) ! ALLOCATE(BSPFIL%BSPFRQ(NFRE)) 40.31 DO IFRE = 1, NFRE BSPFIL%BSPFRQ(IFRE) = FRQ_ARRAY(IFRE)*2*PI 40.31 ENDDO ! IF (ITEST.GE.60) THEN WRITE(PRTEST,*) ' HEDLINT ',' NFRE ',' NANG ',' NBOUNC ', & ' GNAME' WRITE(PRTEST,1944) HEDLINT,NFRE, NANG, NBOUNC, GNAME WRITE (PRTEST,*) 'Frequencies read from boundary file ', FILENM WRITE (PRTEST,*) (FRQ_ARRAY(IFRE),IFRE = 1,NFRE) ENDIF ! ! read direction from WW3 boundary file ! DIR_ARRAY(IANG) = MOD(2.5*PI-TH(ITH),TPI) ! there are in radians but is not in right order related to SWAN ! READ (NDSD,1946) (DIR_ARRAY(IANG),IANG=1,NANG) ! ! put values in right order. The value of the DIR_ARRAY(i) should be ! smaller that the DIR_ARRAY(i-1) ! in the opposite situation make DIR_ARRAY(i) = DIR_ARRAY(i) - 2*PI ! DIR_ARRAY(:) = PI*DNORTH/180 - DIR_ARRAY(:) 40.15 ALLOCATE(BSPFIL%BSPDIR(NANG)) 40.31 DO IANG = 1, NANG IF (IANG.EQ.1) THEN BSPFIL%BSPDIR(1) = DIR_ARRAY(IANG) 40.31 DIRRD1 = BSPFIL%BSPDIR(1) 40.31 ELSE IF (DIR_ARRAY(IANG).LT.DIRRD1) THEN BSPFIL%BSPDIR(IANG) = 2*PI+DIR_ARRAY(IANG) 40.31 DIRRD1 = BSPFIL%BSPDIR(IANG) 40.31 ELSE BSPFIL%BSPDIR(IANG) = DIR_ARRAY(IANG) 40.31 DIRRD1 = BSPFIL%BSPDIR(IANG) 40.31 ENDIF ENDIF ENDDO IF(ITEST.GE.60) THEN WRITE (PRTEST,*) 'Directions read from boundary file ', & FILENM WRITE (PRTEST,1946) (DIR_ARRAY(IANG),IANG = 1,NANG) ENDIF ! ! Time READ (NDSD, 900) WWDATE,WWTIME ! ! Read from boundary file info about the boundary points(b.p): name of b. p., ! geographical location of b.p., depth, wind u-velocity and direction at the b.p. ! current velocity and direction at the b.p. ! ! If DONALL = .TRUE. boundary data correspond to an open boundary otherwise ! it is continue the interpolation of the grid point between the last and the ! first point ! DO IBOUNC = 1, NBOUNC IERR = 0 ! ! latitude = XLAT ! longitude = XLON ! A real which is not used in the computation ! READ (NDSD,901) PTNME, XLAT, XLON, DUM_A, DUM_A, & DUM_A, DUM_A, DUM_A ! Pass over the lines where the energy spectra is written in the boundary file. ! The energy spectra is going to be read later, in the subroutine RESPEC ! READ (NDSD,902) ((DUM_A, IFRE = 1,NFRE),IANG = 1,NANG) ! IF (ITEST.GE.80) THEN WRITE (PRTEST, *) ' B. spectrum WW3 ', IBOUNC, XLON, & XLAT, IERR ENDIF ! ! in case of nesting coordinates on file are used to determine interpolation ! coefficients ! IF (KSPHER.EQ.0) THEN ! ! if SWAN uses Cartesian coordinates, then transform the spherical coordinates ! of the boundary point to local Cartesian coordinates ! IF (IBOUNC.EQ.1) THEN CALL INDBLE('XGC',XLON0,'REQ',-999.D0) CALL INDBLE('YGC',XLAT0,'REQ',-999.D0) IF (XLON0.LT.-900.) THEN XLON0 = (XOFFS -XPC)/LENDEG XLAT0 = (YOFFS-YPC)/LENDEG ENDIF ENDIF ! XP2 = XPC + LENDEG*COS(PI*XLAT0/180.)*(XLON-XLON0) YP2 = YPC + LENDEG*(XLAT-XLAT0) ! ELSE XP2 = XLON-XOFFS YP2 = XLAT-YOFFS ENDIF ! ! --- interpolate the boundaries points to the grid points of ! the SWAN computational grid ! NBGRPT_PREV = NBGRPT CALL SWBCPT ( LBGP, XCGRID, YCGRID, 40.41 40.31 & KGRPNT, XYTST, KGRBND,XP2,YP2,IBOUNC, & NBOUNC, DONALL ) ! check if the grid points are on nested boundary. ! if not, stop the calculation and give an error message IF (NBGRPT.NE.NBGRPT_PREV) THEN IIPT2 = IIPT2+1 ENDIF ENDDO ! IF (IIPT2.EQ.0) CALL MSGERR (2, & 'no grid points on nested boundary') ! IF (ITEST.GE.60) WRITE (PRTEST,16) NBOUNC 16 FORMAT (I6, ' boundary locations') ! ! quantity on file is energy density BSPFIL%BFILED(17) = 1 40.31 ! NHEDF = NHEDF + CEILING(NFRE/8.) + CEILING(NANG/7.)+1 NHEDS = 2 ! NHEDT: calculated in the RBFILE subroutine for each time step ! ALLOCATE(BSPFIL%BSPLOC(NBOUNC)) 40.31 DO IBC = 1, NBOUNC BSPFIL%BSPLOC(IBC) = NBSPEC + IBC 40.31 ENDDO ! NBSPEC = NBSPEC + NBOUNC ! ! store file reading parameters in array BFILED ! BSPFIL%BFILED(1) = ISTATF 40.31 BSPFIL%BFILED(2) = -999999999 40.31 BSPFIL%BFILED(3) = -999999999 40.31 BSPFIL%BFILED(4) = NDSL 40.31 BSPFIL%BFILED(5) = NDSD 40.31 BSPFIL%BFILED(6) = IOPTT 40.31 CALL COPYCH (BTYPE, 'T', BSPFIL%BFILED(7), 1, IERR) 40.31 BSPFIL%BFILED(8) = NBOUNC 40.31 BSPFIL%BFILED(9) = DORDER 40.31 BSPFIL%BFILED(10) = NANG 40.31 BSPFIL%BFILED(11) = 0 40.31 BSPFIL%BFILED(12) = NFRE 40.31 ! ordering of data on file BSPFIL%BFILED(13) = 0 40.31 ! number of heading lines: per file, per time, per spectrum BSPFIL%BFILED(14) = NHEDF 40.31 BSPFIL%BFILED(15) = NHEDT 40.31 BSPFIL%BFILED(16) = NHEDS 40.31 ! quantity on file is energy density (k) BSPFIL%BFILED(17) = 3 40.31 ! IF (ITEST.GE.80) WRITE(PRINTF,81) NBFILS, NBSPEC, & (BSPFIL%BFILED(II), II=1,16) 40.31 ! ! Rewind input file for proper start REWIND (NDSD) 81 FORMAT (' array BFILED: ', 2I4, 2(/,8I10)) 900 FORMAT (I8.8,I7.6) 901 FORMAT (A12,2F7.2,F10.1,2(F7.2,F6.1)) 902 FORMAT (7E11.3) 1944 FORMAT (A23,1X,3I6,1X,A33) 1945 FORMAT (8E10.3) 1946 FORMAT (7E11.3) DEALLOCATE (FRQ_ARRAY,DIR_ARRAY) RETURN END SUBROUTINE BCWW3N !*********************************************************************** ! SUBROUTINE SWBCPT ( LBGP, XCGRID, YCGRID, 40.41 40.31 & KGRPNT, XYTST, KGRBND,XP2,YP2,IBOUNC, & NBOUNC,DONALL ) ! !************************************************************************ ! USE OCPCOMM4 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 USE SWCOMM4 40.41 USE M_BNDSPEC USE M_PARALL ! IMPLICIT NONE ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 30.73, 40.13: Nico Booij ! 40.05: Ekaterini E. Kriezi ! 40.31: Tim Campbell and John Cazes ! 40.31: Marcel Zijlema ! 40.41: Marcel Zijlema ! ! 1. Updates ! ! 40.05, Aug. 00: remove the code to a separate subroutine ! 40.13, Jan. 01: remove declarations of unused variables ! 40.31, Jul. 03: initializations XP0, XP1, YP0, YP1 ! 40.31, Nov. 03: removing POOL mechanism ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. Purpose ! ! interpolation of own boundary points to nesting boundary points ! ! 3. Method ! ! calculate interpolation coefficients between the nesting boundary grid and ! the SWAN computational grid ! ! 4. Argument variables ! ! NBOUNC: max number of boundaries points ! IBOUNC: counter ! ! INTEGER, INTENT(IN) :: KGRPNT(MXC,MYC) ! indirect addresses of computational grid points INTEGER, INTENT(IN) :: KGRBND(*) ! array of boundary grid points INTEGER, INTENT(IN) :: XYTST(*) ! array of (ix,iy) of test points ! INTEGER, INTENT(IN) :: IBOUNC,NBOUNC ! REAL, INTENT(IN) :: XCGRID(MXC,MYC), YCGRID(MXC,MYC) ! coordinates of computational grid points ! ! DONALL : logic arguments declare if the nesting boundary is open or close ! it is defined by the users ! LOGICAL, INTENT(INOUT) :: DONALL LOGICAL :: LBGP 40.31 ! ! 5. Parameter variables ! ! -- ! ! 6. Local variables ! ! XP1 problem coordinate of a boundary location ! YP1 problem coordinate of a boundary location ! XP2 problem coordinate of a boundary location ! YP2 problem coordinate of a boundary location ! DISXY distance in (x,y)-space ! W2 is relative length of projection on line ! IIPT1 counter checking the grid points related to the nesting boundary ! INTEGER, SAVE :: IBSP0 = 1, IBSP1 = 1, IENT = 0 INTEGER :: IBSP2,IGRBND,INDXGR INTEGER :: IXP,IYP,IIPT1 ! REAL, SAVE :: XP0=0., XP1=0., YP0=0., YP1=0. 40.31 REAL :: XP2, YP2, XP, YP, RX, RY, RL2 REAL :: DISXY, W2 TYPE(BGPDAT), POINTER :: BGPTMP 40.31 ! ! 7. Common blocks used ! !OMP INCLUDE 'tpcomm.inc' 40.41 ! ! 8. Subroutines Used ! ! ! 9. Subroutines calling ! ! BCFILE : open and read Swan nesting files ! BCWW3N : open and read WW3 nesting files ! ! 10. Error messages ! ! --- ! ! 11. Remarks ! ! 12. Structure ! ! For all computational grid points on boundary do ! if point is located between nest file grid points ! calculate interpolation coefficients ! if DONALL is TRUE ! the nesting boundary remain open ! else DONALL is FALSE (default case) ! boundary is close, it do interpolation between the last and the first point ! put interpolation coefficients into array BGRIDP ! ! 13. Source text ! CALL STRACE (IENT, 'SWBCPT') 40.41 ! IBSP2 = NBSPEC+IBOUNC IIPT1 = 0 ! 201 IF (IBOUNC.EQ.1) THEN XP0 = XP2 YP0 = YP2 IBSP0 = IBSP2 IIPT1 = 1 40.41 ELSE RX = XP2 - XP1 RY = YP2 - YP1 RL2 = RX**2 + RY**2 IF (RL2.GT.0.) THEN RX = RX/RL2 RY = RY/RL2 ! ! loop over boundary of comp. grid, select points between ! (XP1,YP1) and (XP2,YP2) ! KGRBND grid addresses on boundary points, NGRBND number of grid points ! on computational grid boundary ! DO IGRBND = 1, NGRBND IXP = KGRBND(2*IGRBND-1) IYP = KGRBND(2*IGRBND) ! IF (IXP.GT.0 .AND.IYP.GT.0) THEN ! ! (IXP,IYP) is a valid boundary point ! INDXGR = KGRPNT(IXP,IYP) XP = XCGRID(IXP,IYP) YP = YCGRID(IXP,IYP) ! ! DISXY is relative distance from (XP,YP) to line ! (XP1,YP1) to (XP2,YP2) with respect to the ! length of that line ! DISXY = ABS(RX*(YP-YP1)-RY*(XP-XP1)) ! IF (DISXY.LT.0.1) THEN ! W2 is relative length of projection on line ! (XP1,YP1) to (XP2,YP2) W2 = RX*(XP-XP1)+RY*(YP-YP1) IF (W2.GT.-0.001 .AND. W2.LT.1.001) THEN IF (W2.LT.0.01) W2 = 0. IF (W2.GT.0.99) W2 = 1. IF (ITEST.GE.80) WRITE (PRTEST, *) ' B.pnt', & IXP, IYP, XP, YP, W2, IBSP2, 1.-W2, IBSP1 40.41 NBGRPT = NBGRPT + 1 IIPT1 = IIPT1 + 1 ! ALLOCATE(BGPTMP) 40.31 BGPTMP%BGP(1) = INDXGR 40.31 ! next item indicates type of boundary condition BGPTMP%BGP(2) = 1 40.31 BGPTMP%BGP(3) = NINT(1000. * W2) 40.31 BGPTMP%BGP(4) = IBSP2 40.31 BGPTMP%BGP(5) = NINT(1000. * (1.-W2)) 40.31 BGPTMP%BGP(6) = IBSP1 40.31 NULLIFY(BGPTMP%NEXTBGP) 40.31 IF ( .NOT.LBGP ) THEN 40.31 FBGP = BGPTMP 40.31 CUBGP => FBGP 40.31 LBGP = .TRUE. 40.31 ELSE 40.31 CUBGP%NEXTBGP => BGPTMP 40.31 CUBGP => BGPTMP 40.31 END IF 40.31 ! ! test output if point is a test point ! IF (NPTST.GT.0) THEN DO IPTST = 1, NPTST IF (IXP.EQ.XYTST(2*IPTST-1)+MXF-1 .AND. & IYP.EQ.XYTST(2*IPTST )+MYF-1) & WRITE (PRTEST, 223) & IXP-1,IYP-1, XP+XOFFS, YP+YOFFS, W2, IBSP2, & IBSP1 223 FORMAT (' B.pnt', 2I5, 2F9.0, F6.3, 2I3) ENDDO ENDIF ENDIF ENDIF ENDIF ENDDO ENDIF ENDIF ! IF (IIPT1.EQ.0) THEN WRITE (PRINTF, 218) XP1+XOFFS, YP1+YOFFS, & XP2+XOFFS, YP2+YOFFS 218 FORMAT (' Warning: no grid points on interval from ', 2F12.4, & ' to ', 2F12.4) ENDIF ! XP1 = XP2 YP1 = YP2 IBSP1 = IBSP2 ! IF (IBOUNC.EQ.NBOUNC) THEN IF (.NOT. DONALL) THEN ! ! process grid points between last and first boundary point ! DONALL = .TRUE. XP2 = XP0 YP2 = YP0 IBSP2 = IBSP0 GOTO 201 ENDIF ENDIF RETURN ! END SUBROUTINE SWBCPT ! !********************************************************************* ! * LOGICAL FUNCTION BOUNPT (IX,IY,KGRPNT) ! * !********************************************************************* ! USE SWCOMM3 40.41 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 40.00, 40.03: Nico Booij ! ! 1. UPDATE ! ! Feb. 1998, ver. 40.00: new subroutine ! 40.03, Sep. 00: inconsistency with manual corrected ! ! 2. PURPOSE ! ! determine whether a grid point is a point where a boundary condition ! can be applied ! ! 3. METHOD ! ! ! 4. PARAMETERLIST ! ! IX, IY int inp grid point indices ! KGRPNT int inp indirect addresses of grid points ! ! 5. SUBROUTINES CALLING ! ! BCFILE ! ! 6. SUBROUTINES USED ! ! --- ! ! 7. ERROR MESSAGES ! ! --- ! ! 8. REMARKS ! ! KGRPNT(IX,IY)=1 means that (IX,IY) is not an active grid point ! ! 9. STRUCTURE ! ! ----------------------------------------------------------------- ! Make BOUNPT = False ! If the grid point is not active ! Then return ! ----------------------------------------------------------------- ! If grid point is on the outer boundary ! Then make BOUNPT = True ! return ! ----------------------------------------------------------------- ! If a neighbouring grid point is inactive ! Then make BOUNPT = True ! return ! ----------------------------------------------------------------- ! ! 10. SOURCE TEXT ! INTEGER IX,IY,KGRPNT(MXC,MYC) SAVE IENT DATA IENT/0/ CALL STRACE (IENT, 'BOUNPT') ! BOUNPT = .FALSE. ! IF (IX.LE.0) RETURN IF (IY.LE.0) RETURN IF (IX.GT.MXC) RETURN IF (IY.GT.MYC) RETURN ! ! If the grid point is not active ! Then return ! IF (KGRPNT(IX,IY).LE.1) RETURN ! ! If grid point is on the outer boundary ! Then make BOUNPT = True ! return ! IF (IX.EQ.1) THEN BOUNPT = .TRUE. RETURN ENDIF IF (IX.EQ.MXC) THEN BOUNPT = .TRUE. RETURN ENDIF IF (IY.EQ.1) THEN BOUNPT = .TRUE. RETURN ENDIF IF (IY.EQ.MYC) THEN BOUNPT = .TRUE. RETURN ENDIF ! ! If a neighbouring grid point is inactive ! Then make BOUNPT = True ! return ! IF (KGRPNT(IX-1,IY).LE.1) THEN BOUNPT = .TRUE. RETURN ENDIF IF (KGRPNT(IX+1,IY).LE.1) THEN BOUNPT = .TRUE. RETURN ENDIF IF (KGRPNT(IX,IY-1).LE.1) THEN BOUNPT = .TRUE. RETURN ENDIF IF (KGRPNT(IX,IY+1).LE.1) THEN BOUNPT = .TRUE. RETURN ENDIF RETURN END !********************************************************************* ! * SUBROUTINE RETSTP (MPTST, XYTST, KGRPNT, KGRBND, XCGRID, YCGRID, & SPCSIG, SPCDIR) 40.31 ! * !********************************************************************* ! USE OCPCOMM4 40.41 USE SWCOMM1 40.41 USE SWCOMM2 40.41 USE SWCOMM3 40.41 USE SWCOMM4 40.41 USE OUTP_DATA 40.31 USE M_PARALL 40.31 ! ! ! --|-----------------------------------------------------------|-- ! | Delft University of Technology | ! | Faculty of Civil Engineering | ! | Environmental Fluid Mechanics Section | ! | P.O. Box 5048, 2600 GA Delft, The Netherlands | ! | | ! | Programmers: R.C. Ris, N. Booij, | ! | IJ.G. Haagsma, A.T.M.M. Kieftenburg, | ! | M. Zijlema, E.E. Kriezi, | ! | R. Padilla-Hernandez, L.H. Holthuijsen | ! --|-----------------------------------------------------------|-- ! ! ! SWAN (Simulating WAves Nearshore); a third generation wave model ! Copyright (C) 2004-2005 Delft University of Technology ! ! This program is free software; you can redistribute it and/or ! modify it under the terms of the GNU General Public License as ! published by the Free Software Foundation; either version 2 of ! the License, or (at your option) any later version. ! ! This program is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! A copy of the GNU General Public License is available at ! http://www.gnu.org/copyleft/gpl.html#SEC3 ! or by writing to the Free Software Foundation, Inc., ! 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA ! ! ! 0. Authors ! ! 30.82: IJsbrand Haagsma ! 40.00, 40.13: Nico Booij ! 34.01: Jeroen Adema ! 40.04: Annette Kieftenburg ! 40.30: Marcel Zijlema ! 40.31: Marcel Zijlema ! 40.41: Marcel Zijlema ! ! 1. Updates ! ! 40.00, Apr. 98: New subroutine ! 30.82, Oct. 98: Updated description of several variables ! 34.01, Feb. 99: Introducing STPNOW ! 40.03, Dec. 99: XOFFS and YOFFS introduced in write statements ! May 00: ITMOPT written to heading of file instead of 1 ! 40.04, Aug. 00: Added error message if testpoints are defined ! before bottom grid is read ! 40.13, Jan. 01: two output strings corrected ! May 01: two incorrect units changed from m2/2 to m2/s ! 40.30, Mar. 03: introduction distributed-memory approach using MPI ! 40.31, Dec. 03: removing POOL-mechanism ! 40.41, Oct. 04: common blocks replaced by modules, include files removed ! ! 2. PURPOSE ! ! read test points, generate output point set 'TESTPNTS', ! read source term filenames ! ! 3. METHOD ! ! ! 4. Argument variables ! ! i SPCDIR: (*,1); spectral directions (radians) 30.82 ! (*,2); cosine of spectral directions 30.82 ! (*,3); sine of spectral directions 30.82 ! (*,4); cosine^2 of spectral directions 30.82 ! (*,5); cosine*sine of spectral directions 30.82 ! (*,6); sine^2 of spectral directions 30.82 ! i SPCSIG: Relative frequencies in computational domain in sigma-space 30.82 ! i XCGRID: Coordinates of computational grid in x-direction 30.82 ! i YCGRID: Coordinates of computational grid in y-direction 30.82 ! REAL SPCDIR(MDC,6) 30.82 REAL SPCSIG(MSC) 30.82 REAL XCGRID(MXC,MYC), YCGRID(MXC,MYC) 30.82 ! ! i MPTST : Maximum number of test points 30.82 ! i XYTST : Grid point indices of test points 30.82 ! INTEGER MPTST 30.82 INTEGER XYTST(2*MPTST) 30.82 ! ! 5. SUBROUTINES CALLING ! ! SWREAD ! ! 6. SUBROUTINES USED ! LOGICAL STPNOW 34.01 ! ! 7. Common blocks used ! !OMP INCLUDE 'tpcomm.inc' 40.41 ! ! 8. REMARKS ! ! ! 9. STRUCTURE ! ! ----------------------------------------------------------------- ! Repeat ! read (i,j) identifying a test point ! store values in array XYTST ! ----------------------------------------------------------------- ! Generate output point set 'TESTPNTS' ! write coordinates of test points into array OUTDA ! ----------------------------------------------------------------- ! If 1D output of source terms is requested ! Then open file ! write general data into the file ! ----------------------------------------------------------------- ! If 2D output of source terms is requested ! Then open file ! write general data into the file ! ----------------------------------------------------------------- ! ! 10. SOURCE TEXT ! INTEGER KGRPNT(MXC,MYC), KGRBND(*) 40.31 LOGICAL KEYWIS, LOCGRI 40.00 TYPE(OPSDAT), POINTER :: OPSTMP 40.31 SAVE IENT DATA IENT /0/ CALL STRACE (IENT, 'RETSTP') ! CALL INKEYW ('STA','IJ') 40.00 IF (LEDS(1).GE.2) THEN 40.04 IF (KEYWIS('XY')) THEN LOCGRI = .TRUE. ELSE IF (KEYWIS('IJ')) THEN LOCGRI = .FALSE. ELSE CALL WRNKEY ENDIF ELSE 40.04 CALL MSGERR(3,'Testpoints must be defined after bottom is read') 40.04 ENDIF 40.04 ! 10 IF (LOCGRI) THEN CALL READXY ('X','Y',XP,YP, 'REP', -1.E10, -1.E10) 40.03 IF (XP.LT.-.9E10) THEN LXDMP = 0 GOTO 60 ELSE CALL CVMESH (XP, YP, XC, YC, KGRPNT, XCGRID, YCGRID, KGRBND) 40.00 IF (XC.LT.0.) THEN 40.31 IF (XP.GE.XCGMIN .AND. XP.LE.XCGMAX .AND. 40.31 & YP.GE.YCGMIN .AND. YP.LE.YCGMAX ) THEN 40.31 GOTO 50 40.31 ELSE 40.31 GOTO 40 40.31 END IF 40.31 END IF 40.31 LXDMP = NINT(XC) + MXF -1 40.31 LYDMP = NINT(YC) + MYF -1 40.31 IF (ITEST.GE.30) WRITE (PRTEST, 14) XP+XOFFS, YP+YOFFS, 40.03 & LXDMP, LYDMP 14 FORMAT (' test point ', 2F12.2, ' to grid point ', 2I4) ENDIF ELSE CALL ININTG ('I' , LXDMP, 'REP', -1) 40.03 IF (LXDMP .LT. 0) GOTO 60 CALL ININTG ('J' , LYDMP, 'REQ', 0) 40.03 ENDIF IF (LXDMP.GE.0 .AND. LXDMP.LE.MXCGL-1 .AND. 40.31 & LYDMP.GE.0 .AND. LYDMP.LE.MYCGL-1) THEN 40.31 LXDMP = LXDMP - MXF + 1 40.31 LYDMP = LYDMP - MYF + 1 40.31 IF (LXDMP.GE.0 .AND. LXDMP.LE.MXC-1 .AND. & LYDMP.GE.0 .AND. LYDMP.LE.MYC-1) THEN IF (KGRPNT(LXDMP+1,LYDMP+1) .GT. 1) THEN NPTST = NPTST + 1 XYTST(2*NPTST-1) = LXDMP+1 XYTST(2*NPTST) = LYDMP+1 GOTO 50 ENDIF ELSE 40.31 GOTO 50 40.31 ENDIF 40.31 ENDIF 40 CALL MSGERR (2, 'test point outside comp. grid') IF (LOCGRI) THEN WRITE (PRINTF, *) XP+XOFFS, YP+YOFFS 40.03 ELSE WRITE (PRINTF, *) LXDMP, LYDMP ENDIF 50 IF (NPTST.LE.MPTST) GOTO 10 CALL MSGERR (2, 'Too many test points') ! ! generate output point set 'TESTPNTS' ! 60 ALLOCATE(OPSTMP) 40.31 OPSTMP%PSNAME = 'TESTPNTS' 40.31 OPSTMP%PSTYPE = 'P' 40.31 OPSTMP%MIP = NPTST 40.31 ALLOCATE(OPSTMP%XP(NPTST)) 40.31 ALLOCATE(OPSTMP%YP(NPTST)) 40.31 DO 140 IPTST = 1, NPTST LXDMP = XYTST(2*IPTST-1) 40.00 LYDMP = XYTST(2*IPTST) 40.00 OPSTMP%XP(IPTST) = XCGRID(LXDMP,LYDMP) 40.31 OPSTMP%YP(IPTST) = YCGRID(LXDMP,LYDMP) 40.31 140 CONTINUE NULLIFY(OPSTMP%NEXTOPS) 40.31 IF ( .NOT.LOPS ) THEN 40.31 FOPS = OPSTMP 40.31 COPS => FOPS 40.31 LOPS = .TRUE. 40.31 ELSE 40.31 COPS%NEXTOPS => OPSTMP 40.31 COPS => OPSTMP 40.31 END IF 40.31 ! ! open output file for test output of wave parameters 40.00 ! CALL INKEYW ('STA', ' ') 40.00 IF (KEYWIS('PAR')) THEN 40.00 CALL INCSTR ('FNAME', FILENM, 'STA', 'SWSRCPA') 40.00 ! --- append node number to FILENM in case of 40.30 ! parallel computing 40.30 IF ( PARLL ) THEN 40.30 ILPOS = INDEX ( FILENM, ' ' )-1 40.30 WRITE(FILENM(ILPOS+1:ILPOS+4),99) INODE 40.30 99 FORMAT('-',I3.3) 40.30 END IF 40.30 IERR = 0 40.00 CALL FOR (IFPAR, FILENM, 'UF', IERR) 40.00 IF (STPNOW()) RETURN 34.01 WRITE (IFPAR, 101) 1 40.00 101 FORMAT ('SWAN', I4, T41, & 'Swan standard spectral file, version') 40.00 WRITE (IFPAR, 111) VERTXT 40.03 111 FORMAT ('$ Data produced by SWAN version ', A) 40.03 WRITE (IFPAR, 113) PROJID, PROJNR 40.03 113 FORMAT ('$ Project: ', A, '; run number: ', A) IF (NSTATM.EQ.1) THEN WRITE (IFPAR, 102) 'TIME', 'time-dependent data' 102 FORMAT (A, T41, A) 40.00 WRITE (IFPAR, 103) ITMOPT, 'time coding option' 40.03 103 FORMAT (I6, T41, A) 40.00 ELSE WRITE (IFPAR, 102) 'ITER', 'iteration-dependent data' WRITE (IFPAR, 103) 0 ENDIF IF (KSPHER.EQ.0) THEN WRITE (IFPAR, 102) 'LOCATIONS', 'locations in x-y-space' ELSE WRITE (IFPAR, 102) 'LONLAT', & 'locations in longitude, latitude' ENDIF WRITE (IFPAR, 103) NPTST, 'number of locations' DO 110 IPTST = 1, NPTST LXDMP = XYTST(2*IPTST-1) 40.00 LYDMP = XYTST(2*IPTST) 40.00 WRITE (IFPAR, 106) XCGRID(LXDMP,LYDMP)+XOFFS, & YCGRID(LXDMP,LYDMP)+YOFFS 40.00 106 FORMAT (2(1X,F12.2)) 110 CONTINUE WRITE (IFPAR, 132) 9 40.55 40.00 132 FORMAT ('QUANT', /, I6, T41, 'number of quantities in table') 40.00 WRITE (IFPAR, 102) OVSNAM(10), OVLNAM(10) 40.00 WRITE (IFPAR, 102) OVUNIT(10), 'unit' 40.00 WRITE (IFPAR, 104) OVEXCV(10), 'exception value' 40.00 104 FORMAT (F14.6, T41, A) 40.00 WRITE (IFPAR, 102) OVSNAM(28), OVLNAM(28) 40.41 40.00 WRITE (IFPAR, 102) OVUNIT(28), 'unit' 40.41 40.00 WRITE (IFPAR, 104) OVEXCV(28), 'exception value' 40.41 40.00 WRITE (IFPAR, 102) 'Swind', 'wind source term (of var. dens.)' 40.00 WRITE (IFPAR, 102) 'm2/s', 'unit' 40.00 WRITE (IFPAR, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFPAR, 102) 'Swcap', 'whitecapping dissipation' 40.00 WRITE (IFPAR, 102) 'm2/s', 'unit' 40.00 WRITE (IFPAR, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFPAR, 102) 'Sfric', 'bottom friction dissipation' 40.00 WRITE (IFPAR, 102) 'm2/s', 'unit' 40.00 WRITE (IFPAR, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFPAR, 102) 'Smud', 'fluid mud dissipation' 40.55 WRITE (IFPAR, 102) 'm2/s', 'unit' 40.55 WRITE (IFPAR, 104) OVEXCV(7),'exception value' 40.55 WRITE (IFPAR, 102) 'Ssurf', 'surf breaking dissipation' 40.00 WRITE (IFPAR, 102) 'm2/s', 'unit' 40.00 WRITE (IFPAR, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFPAR, 102) 'Snl3', 'total absolute 3-wave interaction' 40.13 WRITE (IFPAR, 102) 'm2/s', 'unit' 40.13 WRITE (IFPAR, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFPAR, 102) 'Snl4', 'total absolute 4-wave interaction' 40.13 WRITE (IFPAR, 102) 'm2/s', 'unit' 40.13 WRITE (IFPAR, 104) OVEXCV(7),'exception value' 40.00 ENDIF ! ! open output file for source terms if requested ! CALL INKEYW ('STA', ' ') 40.00 IF (KEYWIS('S1D')) THEN 40.00 CALL INCSTR ('FNAME', FILENM, 'STA', 'SWSRC1D') 40.00 ! --- append node number to FILENM in case of 40.30 ! parallel computing 40.30 IF ( PARLL ) THEN 40.30 ILPOS = INDEX ( FILENM, ' ' )-1 40.30 WRITE(FILENM(ILPOS+1:ILPOS+4),99) INODE 40.30 END IF 40.30 IERR = 0 40.00 CALL FOR (IFS1D, FILENM, 'UF', IERR) 40.00 IF (STPNOW()) RETURN 34.01 WRITE (IFS1D, 101) 1 WRITE (IFS1D, 111) VERTXT 40.03 WRITE (IFS1D, 113) PROJID, PROJNR 40.03 IF (NSTATM.EQ.1) THEN WRITE (IFS1D, 102) 'TIME', 'time-dependent data' WRITE (IFS1D, 103) ITMOPT, 'time coding option' 40.03 ELSE WRITE (IFS1D, 102) 'ITER', 'iteration-dependent data' WRITE (IFS1D, 103) 0 ENDIF IF (KSPHER.EQ.0) THEN WRITE (IFS1D, 102) 'LOCATIONS', 'locations in x-y-space' ELSE WRITE (IFS1D, 102) 'LONLAT', & 'locations in longitude, latitude' ENDIF WRITE (IFS1D, 103) NPTST, 'number of locations' DO 210 IPTST = 1, NPTST 40.00 LXDMP = XYTST(2*IPTST-1) 40.00 LYDMP = XYTST(2*IPTST) 40.00 WRITE (IFS1D, 106) XCGRID(LXDMP,LYDMP)+XOFFS, 40.00 & YCGRID(LXDMP,LYDMP)+YOFFS 40.00 210 CONTINUE 40.00 IF (ICUR.GT.0) THEN WRITE (IFS1D, 102) 'RFREQ', 'relative frequencies in Hz' 40.00 ELSE WRITE (IFS1D, 102) 'AFREQ', 'absolute frequencies in Hz' 40.00 ENDIF WRITE (IFS1D, 103) MSC, 'number of frequencies' 40.00 DO 220 IS = 1, MSC 40.00 WRITE (IFS1D, 214) SPCSIG(IS)/PI2 40.00 214 FORMAT (F10.4) 40.00 220 CONTINUE 40.00 WRITE (IFS1D, 132) 8 40.55 WRITE (IFS1D, 102) 'VaDens', 'variance densities' 40.00 WRITE (IFS1D, 102) 'm2/Hz', 'unit' 40.00 WRITE (IFS1D, 104) 0., 'exception value' 40.00 WRITE (IFS1D, 102) 'Swind', 'wind source term' 40.00 WRITE (IFS1D, 102) 'm2', 'unit' 40.00 WRITE (IFS1D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS1D, 102) 'Swcap', 'whitecapping dissipation' 40.00 WRITE (IFS1D, 102) 'm2', 'unit' 40.00 WRITE (IFS1D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS1D, 102) 'Sfric', 'bottom friction dissipation' 40.00 WRITE (IFS1D, 102) 'm2', 'unit' 40.00 WRITE (IFS1D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS1D, 102) 'Smud', 'fluid mud dissipation' 40.55 WRITE (IFS1D, 102) 'm2', 'unit' 40.55 WRITE (IFS1D, 104) OVEXCV(7),'exception value' 40.55 WRITE (IFS1D, 102) 'Ssurf', 'surf breaking dissipation' 40.00 WRITE (IFS1D, 102) 'm2', 'unit' 40.00 WRITE (IFS1D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS1D, 102) 'Snl3', 'triad interactions' 40.00 WRITE (IFS1D, 102) 'm2', 'unit' 40.00 WRITE (IFS1D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS1D, 102) 'Snl4', 'quadruplet interactions' 40.00 WRITE (IFS1D, 102) 'm2', 'unit' 40.00 WRITE (IFS1D, 104) OVEXCV(7),'exception value' 40.00 ENDIF ! 2D source terms CALL INKEYW ('STA', ' ') 40.00 IF (KEYWIS('S2D')) THEN 40.00 CALL INCSTR ('FNAME', FILENM, 'STA', 'SWSRC2D') 40.00 ! --- append node number to FILENM in case of 40.30 ! parallel computing 40.30 IF ( PARLL ) THEN 40.30 ILPOS = INDEX ( FILENM, ' ' )-1 40.30 WRITE(FILENM(ILPOS+1:ILPOS+4),99) INODE 40.30 END IF 40.30 IERR = 0 40.00 CALL FOR (IFS2D, FILENM, 'UF', IERR) 40.00 IF (STPNOW()) RETURN 34.01 WRITE (IFS2D, 101) 1 WRITE (IFS2D, 111) VERTXT 40.03 WRITE (IFS2D, 113) PROJID, PROJNR 40.03 IF (NSTATM.EQ.1) THEN WRITE (IFS2D, 102) 'TIME', 'time-dependent data' WRITE (IFS2D, 103) ITMOPT, 'time coding option' 40.03 ELSE WRITE (IFS2D, 102) 'ITER', 'iteration-dependent data' WRITE (IFS2D, 103) 0 ENDIF IF (KSPHER.EQ.0) THEN WRITE (IFS2D, 102) 'LOCATIONS', 'locations in x-y-space' ELSE WRITE (IFS2D, 102) 'LONLAT', & 'locations in longitude, latitude' ENDIF WRITE (IFS2D, 103) NPTST, 'number of locations' DO 310 IPTST = 1, NPTST 40.00 LXDMP = XYTST(2*IPTST-1) 40.00 LYDMP = XYTST(2*IPTST) 40.00 WRITE (IFS2D, 106) XCGRID(LXDMP,LYDMP)+XOFFS, 40.00 & YCGRID(LXDMP,LYDMP)+YOFFS 40.00 310 CONTINUE 40.00 IF (ICUR.GT.0) THEN WRITE (IFS2D, 102) 'RFREQ', 'relative frequencies in Hz' 40.00 ELSE WRITE (IFS2D, 102) 'AFREQ', 'absolute frequencies in Hz' 40.00 ENDIF WRITE (IFS2D, 103) MSC, 'number of frequencies' 40.00 DO 320 IS = 1, MSC 40.00 WRITE (IFS2D, 214) SPCSIG(IS)/PI2 40.00 320 CONTINUE 40.00 ! full 2-D spectrum WRITE (IFS2D, 102) 'CDIR', & 'spectral Cartesian directions in degr' 40.00 WRITE (IFS2D, 103) MDC, 'number of directions' DO 330 ID = 1, MDC 40.00 WRITE (IFS2D, 324) SPCDIR(ID,1)*180./PI 30.82 324 FORMAT (F10.4) 40.00 330 CONTINUE 40.00 WRITE (IFS2D, 132) 8 40.55 WRITE (IFS2D, 102) 'VaDens', 'variance densities' 40.00 WRITE (IFS2D, 102) 'm2/Hz/degr', 'unit' 40.00 WRITE (IFS2D, 104) 0., 'exception value' 40.00 WRITE (IFS2D, 102) 'Swind', 'wind source term' 40.00 WRITE (IFS2D, 102) 'm2/degr','unit' 40.00 WRITE (IFS2D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS2D, 102) 'Swcap', 'whitecapping dissipation' 40.00 WRITE (IFS2D, 102) 'm2/degr','unit' 40.00 WRITE (IFS2D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS2D, 102) 'Sfric', 'bottom friction dissipation' 40.00 WRITE (IFS2D, 102) 'm2/degr','unit' 40.00 WRITE (IFS2D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS2D, 102) 'Smud', 'fluid mud dissipation' 40.55 WRITE (IFS2D, 102) 'm2/degr','unit' 40.55 WRITE (IFS2D, 104) OVEXCV(7),'exception value' 40.55 WRITE (IFS2D, 102) 'Ssurf', 'surf breaking dissipation' 40.00 WRITE (IFS2D, 102) 'm2/degr','unit' 40.00 WRITE (IFS2D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS2D, 102) 'Snl3', 'triad interactions' 40.00 WRITE (IFS2D, 102) 'm2/degr','unit' 40.00 WRITE (IFS2D, 104) OVEXCV(7),'exception value' 40.00 WRITE (IFS2D, 102) 'Snl4', 'quadruplet interactions' 40.00 WRITE (IFS2D, 102) 'm2/degr','unit' 40.00 WRITE (IFS2D, 104) OVEXCV(7),'exception value' 40.00 ENDIF RETURN ! end of subroutine RETSTP END