qwind2.f

Go to the documentation of this file.

!                   *****************
                    SUBROUTINE qwind2
!                   *****************
!
     &( tstot , tsder , f     , xk    , usold , usnew , twold , twnew ,
     &  nf    , ndire , npoin2, usn   , uso   )
!
!***********************************************************************
! TOMAWAC   V6P3                                   27/06/2011
!***********************************************************************
!
!brief    COMPUTES THE CONTRIBUTION OF THE WAVE GENERATION
!+               (BY WIND) SOURCE TERM BASED ON SNYDER ET AL. (1981),
!+                MODIFIED BY KOMEN ET AL. (1984) TO MAKE USE OF THE
!+                FRICTION VELOCITY U* INSTEAD OF THE U5 VELOCITY
!+               (MEASURED 5 METERS ABOVE) FOR THE WIND.
!+
!+            THIS GENERATION THEORY IS IDENTICAL TO THAT IN WAM-CYCLE 3.
!
!reference  SNYDER ET AL. (1981) :
!+                     "ARRAY MEASUREMENTS OF ATMOSPHERIC PRESSURE
!+                      FLUCTUATIONS ABOVE SURFACE GRAVITY WAVES".
!+                      JOURNAL OF FLUID MECH., VOL 102., PP 1-59.
!reference KOMEN ET AL.  (1984) :
!+                     "ON THE EXISTENCE OF A FULLY DEVELOPED
!+                      WINDSEA SPECTRUM". JPO, VOL 14, PP 1271-1285.
!
!history  M. BENOIT (EDF/DER/LNH)
!+        26/03/96
!+        V1P1
!+
!
!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
!
!history  J-M HERVOUET (EDF/LNHE)
!+        23/12/2012
!+        V6P3
!+   A first optimisation.
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| F              |-->| DIRECTIONAL SPECTRUM
!| NF             |-->| NUMBER OF FREQUENCIES
!| NDIRE          |-->| NUMBER OF DIRECTIONS
!| NPOIN2         |-->| NUMBER OF POINTS IN 2D MESH
!| 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
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE declarations_tomawac, ONLY : deupi,usdpi,roair,roeau,cimpli,
     &                teta, freq
!FROM TOMAWAC MODULE
! CIMPLI          IMPLICITATION COEFFICIENT FOR SOURCE TERMS
      USE interface_tomawac, ex_qwind2 => qwind2
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER, INTENT(IN)            :: NF,NDIRE,NPOIN2
      DOUBLE PRECISION, INTENT(IN)   :: TWOLD(npoin2),TWNEW(npoin2)
      DOUBLE PRECISION, INTENT(IN)   :: USNEW(npoin2),USOLD(npoin2)
      DOUBLE PRECISION, INTENT(INOUT):: USO(npoin2,ndire)
      DOUBLE PRECISION, INTENT(INOUT):: USN(npoin2,ndire)
      DOUBLE PRECISION, INTENT(INOUT):: TSTOT(npoin2,ndire,nf)
      DOUBLE PRECISION, INTENT(INOUT):: TSDER(npoin2,ndire,nf)
      DOUBLE PRECISION, INTENT(IN)   :: F(npoin2,ndire,nf),XK(npoin2,nf)
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER JP,JF,IP
      DOUBLE PRECISION C1,DIREC,CONST,COEPHAS,SURDEUPIFREQ,BETAO,BETAN
!
!-----------------------------------------------------------------------
!
      c1 = 0.25d0 * (roair/roeau) * deupi
!
!     ARRAYS DEPENDING ONLY ON POINTS AND DIRECTION
!     COULD BE OPTIMISED MORE BY DECOMPOSING THE COS...
!
      DO jp=1,ndire
        direc=teta(jp)
        DO ip=1,npoin2
          uso(ip,jp)=28.d0*usold(ip)*cos(direc-twold(ip))
          usn(ip,jp)=28.d0*usnew(ip)*cos(direc-twnew(ip))
        ENDDO
      ENDDO
!
!     LOOP ON THE DISCRETISED DIRECTIONS
!
      DO jf=1,nf
        const=c1*freq(jf)
        surdeupifreq=usdpi/freq(jf)
        DO jp=1,ndire
!
!         COMPUTES THE PROPORTIONALITY FACTORS BETA
!         AND TAKES THE SOURCE TERM INTO ACCOUNT
!
          DO ip=1,npoin2
            coephas = xk(ip,jf)*surdeupifreq
            betao=max(uso(ip,jp)*coephas-1.d0,0.d0)*const
            betan=max(usn(ip,jp)*coephas-1.d0,0.d0)*const
            tstot(ip,jp,jf) = tstot(ip,jp,jf)
     &                      + (betao+cimpli*(betan-betao))*f(ip,jp,jf)
            tsder(ip,jp,jf)=tsder(ip,jp,jf)+betan
          ENDDO
!
        ENDDO
!
      ENDDO
!
!-----------------------------------------------------------------------
!
      RETURN
      END