perale.f

Go to the documentation of this file.

!                   *****************
                    SUBROUTINE perale
!                   *****************
!
     &(pale,gamma,perpic,npale,pmin,pmax)
!
!***********************************************************************
! ARTEMIS   V6P3                                  21/08/2010
!***********************************************************************
!
!brief    DISCRETISES AN ENERGY SPECTRUM IN NPALE BANDS
!+                OF EQUAL ENERGY. THE RESULT IS A LIST OF
!+                PERIODS CORRESPONDING TO EACH BAND.
!
!history  J-M HERVOUET (LNH)
!+
!+
!+   LINKED TO BIEF 5.0
!
!history  D. AELBRECHT (LNH)
!+        02/06/1999
!+        V5P1
!+
!
!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
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| GAMMA          |-->| COEFFICIENT IN THE SPECTRUM FORMULA
!| NPALE          |-->| NUMBER OF DISCRETISATION PERIOD
!| PALE           |<--| PERIODS FOR SPECTRUM DISCRETISATION
!| PERPIC         |-->| PEAK PERIOD FOR SPECTRUM
!| PMAX           |-->| MAXIMUM FREQUENCY FOR SPECTRUM
!| PMIN           |-->| MINIMUM FREQUENCY FOR SPECTRUM
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE interface_artemis, ONLY: spe
      USE bief
!
      USE declarations_artemis, ONLY: fp,gam,delta
      USE declarations_special
      IMPLICIT NONE
!
      INTEGER NPALE,NPAS,I,K
!
      DOUBLE PRECISION PALE(npale)
      DOUBLE PRECISION PERPIC,GAMMA ,SUMB  ,SUMICI   ,DF    ,VAR
      DOUBLE PRECISION PMIN  ,PMAX  ,FMIN  ,FMAX
      DOUBLE PRECISION B     ,B1    ,B2
!
      INTRINSIC log,float
!
!-----------------------------------------------------------------------
!
! PEAK FREQUENCY
      fp   = 1.d0 / perpic
      fmin = 1.d0 / pmax
      fmax = 1.d0 / pmin
      IF (fmax.GE.99.d0) THEN
        fmax = 2.5d0 * fp
      ENDIF
!
! GAMMA IS IN THE COMMON STATEMENT OF FUNCTION SPE (CANNOT BE
! CALLED GAMMA BECAUSE IT IS AN ARGUMENT OF THIS SUBROUTINE)
      gam  = gamma
!
!-----------------------------------------------------------------------
!
      IF (gamma.GT.0.99d0 .AND. gamma.LT.1.01d0) THEN
!
!
!        PIERSON-MOSKOWITZ SPECTRUM
!        ----------------------------
!
        b1 = exp(-1.25d0 * (fp/fmax)**4)
        b2 = exp(-1.25d0 * (fp/fmin)**4)
        b  = b1 - b2
        DO i=1,npale
          pale(npale-i+1) = perpic *
     &    (-0.8d0*log( b2 + b*float(2*i-1)/float(2*npale) ))**(0.25d0)
        ENDDO
!
      ELSE
!
!
!        JONSWAP SPECTRUM
!        ------------------
!
!        THE FREQUENCIES LIMITING THE SPECTRUM TO THE LEFT AND RIGHT
!        ARE GIVEN BY KEYWORDS IN THE ARGUMENTS
!
        IF (fmax.LE.fp) THEN
          fmax = 2.5d0 * fp
          WRITE(lu,110) fmax
 110      FORMAT(/,1x,'(PERALE) : FMAX < FP ??? =>',1x,
     &          'CORRECTION : FMAX =',1x,f5.3,' HZ',/)
        ENDIF
!
!       NUMBER OF INTEGRATION INTERVALS FOR THE TRAPEZOIDS METHOD
!
        npas = 2000*npale
!
!       WIDTH OF AN INTEGRATION INTERVAL
!
        df = (fmax-fmin)/float(npas)
!
!       COEFFICIENT FOR THE FUNCTION OF THE SPECTRUM (COMPUTED HERE
!       SO THAT IT'S NOT RECOMPUTED WHEN CALLS SPE)
!
        delta = 0.0624d0 * fp**4 /
     &           ( 0.230d0 + 0.0336d0*gam - 0.185d0 / (1.9d0+gam) )
!
!       SURFACE OF THE SPECTRUM (TRAPEZOIDS METHOD)
!
        sumb = (spe(fmin) + spe(fmax))/2.d0
        DO i = 2,npas-1
          sumb = sumb + spe(fmin+float(i)*df)
        ENDDO
!
!       DIVIDES THE SPECTRUM INTO 2*NPALE BANDS OF EQUAL ENERGY
!
        sumb = sumb/float(2*npale)
!
!       IDENTIFIES THE FREQUENCIES EVERY (2*I-1)*SUMB (I=1,NPALE)
!
        sumici = spe(fmin)/2.d0
        i   = 1
        DO k=1,npas
          var = spe(fmin+df*float(k))
          sumici = sumici + var/2.d0
          IF (sumici.GT.sumb*float(2*i-1)) THEN
            pale(npale-i+1) = 1.d0 / ( fmin + df*(float(k)-0.5d0) )
            i = i + 1
            IF (i.GT.npale) RETURN
          ENDIF
            sumici = sumici + var/2.d0
        ENDDO
!
!
      ENDIF
!
!-----------------------------------------------------------------------
!
      RETURN
      END