qwind1.f

Go to the documentation of this file.

!                   *****************
                    SUBROUTINE qwind1
!                   *****************
!
     &( tstot , tsder , f     , xk    , usold , usnew , twold , twnew ,
     &  z0old , z0new , nf    , ndire , npoin2, told  , tnew  ,
     &  usn   , uso   , omnew , omold , betan , betao )
!
!***********************************************************************
! TOMAWAC   V6P1                                   27/06/2011
!***********************************************************************
!
!brief    COMPUTES THE CONTRIBUTION OF THE WAVE GENERATION
!+               (BY WIND) SOURCE TERM BASED ON JANSSEN (1989,1991).
!
!reference  JANSSEN P.A.E.M (1989) :
!+                     "WIND-INDUCED STRESS AND THE DRAG OF AIR
!+                      FLOW OVER SEA WAVES". JPO, VOL 19, PP 745-754.
!reference JANSSEN P.A.E.M (1991) :
!+                     "QUASI-LINEAR THEORY OF WIND-WAVE GENERATION
!+                      APPLIED TO WAVE FORECASTING". JPO, VOL 21, PP 1631-1642.
!
!history  P. THELLIER; M. BENOIT (EDF/DER/LNH)
!+        11/04/95
!+        V1P0
!+
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        13/07/2010
!+        V6P0
!+   Translation of French comments within the FORTRAN sources into
!+   English comments
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        21/08/2010
!+        V6P0
!+   Creation of DOXYGEN tags for automated documentation and
!+   cross-referencing of the FORTRAN sources
!
!history  G.MATTAROLO (EDF - LNHE)
!+        27/06/2011
!+        V6P1
!+   Translation of French names of the variables in argument
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| BETAM          |-->| WIND GENERATION COEFFICIENT
!| BETAN          |<--| WORK TABLE
!| BETAO          |<--| WORK TABLE
!| CPHAS          |<--| WORK TABLE
!| F              |-->| DIRECTIONAL SPECTRUM
!| NF             |-->| NUMBER OF FREQUENCIES
!| NDIRE          |-->| NUMBER OF DIRECTIONS
!| NPOIN2         |-->| NUMBER OF POINTS IN 2D MESH
!| OMNEW          |<--| WORK TABLE
!| OMOLD          |<--| WORK TABLE
!| TNEW           |<--| WORK TABLE
!| TOLD           |<--| WORK TABLE
!| TSDER          |<->| DERIVED PART OF THE SOURCE TERM CONTRIBUTION
!| TSTOT          |<->| TOTAL PART OF THE SOURCE TERM CONTRIBUTION
!| TWNEW          |-->| WIND DIRECTION AT TIME N+1
!| TWOLD          |-->| WIND DIRECTION AT TIME N
!| USN            |<--| WORK TABLE
!| USNEW          |-->| FRICTION VELOCITY AT TIME N+1
!| USO            |<--| WORK TABLE
!| USOLD          |-->| FRICTION VELOCITY AT TIME N
!| XK             |-->| DISCRETIZED WAVE NUMBER
!| XKAPPA         |-->| VON KARMAN CONSTANT
!| Z0NEW          |-->| SURFACE ROUGHNESS LENGTH AT TIME N+1
!| Z0OLD          |-->| SURFACE ROUGHNESS LENGTH AT TIME N
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE declarations_tomawac, ONLY : deupi, gravit, roair, roeau,
     &                    teta, freq, decal, xkappa,betam, cimpli

!
      USE interface_tomawac, ex_qwind1 => qwind1
      IMPLICIT NONE
!
!.....VARIABLES IN ARGUMENT
!     """"""""""""""""""""
      INTEGER, INTENT(IN)   ::  NF    , NDIRE , NPOIN2
      DOUBLE PRECISION, INTENT(IN)   :: F(npoin2,ndire,nf),XK(npoin2,nf)
      DOUBLE PRECISION, INTENT(IN)   :: TWOLD(npoin2), TWNEW(npoin2)
      DOUBLE PRECISION, INTENT(IN)   :: USOLD(npoin2), USNEW(npoin2)
      DOUBLE PRECISION, INTENT(IN)   :: Z0OLD(npoin2), Z0NEW(npoin2)
      DOUBLE PRECISION, INTENT(INOUT):: TNEW(npoin2,ndire)
      DOUBLE PRECISION, INTENT(INOUT):: TOLD(npoin2,ndire)
      DOUBLE PRECISION, INTENT(INOUT):: USO(npoin2),USN(npoin2)
      DOUBLE PRECISION, INTENT(INOUT):: OMNEW(npoin2),OMOLD(npoin2)
      DOUBLE PRECISION, INTENT(INOUT):: BETAN(npoin2), BETAO(npoin2)
      DOUBLE PRECISION, INTENT(INOUT):: TSTOT(npoin2,ndire,nf)
      DOUBLE PRECISION, INTENT(INOUT):: TSDER(npoin2,ndire,nf)
!.....VARIABLES FROM MODULE TOMAWAC
!     """"""""""""""""""""
!     DECAL           SHIFT GROWING CURVE DUE TO WIND
!     XKAPPA          VON KARMAN CONSTANT
!     BETAM           WIND GENERATION COEFFICIENT
!     CIMPLI          IMPLICITATION COEFFICIENT FOR SOURCE TERMS
!
!.....LOCAL VARIABLES
!     """""""""""""""""
      INTEGER  JP    , JF    , IP
      DOUBLE PRECISION C1,DIREC,CONST,DIMPLI,XX,ZLOGMU, CPHAS
!
      c1 = deupi * (roair/roeau) * (betam/xkappa**2)
      dimpli=1.d0-cimpli
!
!.....COMPUTES (1ST PASS) THE DIRECTIONAL DEPENDENCES
!     """"""""""""""""""""""""""""""""""""""""""""""
      DO jp=1,ndire
        direc=teta(jp)
        DO ip=1,npoin2
          told(ip,jp)=cos(direc-twold(ip))
          tnew(ip,jp)=cos(direc-twnew(ip))
        ENDDO
      ENDDO
!
!.....LOOP ON THE DISCRETISED FREQUENCIES
!     """"""""""""""""""""""""""""""""""""""""""""
      DO jf=1,nf
        const=c1*freq(jf)
!
!.......COMPUTES (1ST PASS) THE FREQUENCIES (OMEGA AND UETOILE/CPHASE)
!       """""""""""""""""""""""""""""""""""""""""""""""""
        DO ip=1,npoin2
!.......COMPUTES THE PHASE VELOCITY
          cphas = deupi * freq(jf) / xk(ip,jf)
          omold(ip) = gravit * z0old(ip) / cphas**2
          omnew(ip) = gravit * z0new(ip) / cphas**2
          uso(ip) = (usold(ip) / cphas) + decal
          usn(ip) = (usnew(ip) / cphas) + decal
        ENDDO
!
!.......LOOP ON THE DISCRETISED DIRECTIONS
!       """"""""""""""""""""""""""""""""""""""""""""
        DO jp=1,ndire
!
          DO ip=1,npoin2
            betao(ip)=0.d0
            betan(ip)=0.d0
          ENDDO
!
!.........COMPUTES THE SOURCE TERM
!         """"""""""""""""""""""
          DO ip=1,npoin2
            IF(told(ip,jp).GT.0.01d0) THEN
              xx = uso(ip) * told(ip,jp)
              zlogmu = log(omold(ip)) + xkappa/xx
              IF(zlogmu.LT.0.d0) THEN
                betao(ip) = const*omold(ip)*exp(xkappa/xx)*
     &                          zlogmu**4*xx**2
              ENDIF
            ENDIF
          ENDDO
          DO ip=1,npoin2
            IF(tnew(ip,jp).GT.0.01d0) THEN
              xx = usn(ip) * tnew(ip,jp)
              zlogmu = log(omnew(ip)) + xkappa/xx
              IF(zlogmu.LT.0.d0) THEN
                betan(ip) = const*omnew(ip)*exp(xkappa/xx)*
     &                          zlogmu**4*xx**2
              ENDIF
            ENDIF
          ENDDO
!
!.........TAKES THE SOURCE TERM INTO ACCOUNT
!         """"""""""""""""""""""""""""""""
          DO ip=1,npoin2
            tstot(ip,jp,jf) = tstot(ip,jp,jf)
     &              + (dimpli*betao(ip)+cimpli*betan(ip)) * f(ip,jp,jf)
            tsder(ip,jp,jf) = tsder(ip,jp,jf) + betan(ip)
          ENDDO
!
        ENDDO
      ENDDO
!
      RETURN
      END