qgauss.f

Go to the documentation of this file.

!                   ***************
                    FUNCTION qgauss
!                   ***************
!
     &( b     , n     , a     , xm    )
!
!***********************************************************************
! TOMAWAC   V6P1                                   23/06/2011
!***********************************************************************
!
!brief    COMPUTES THE INTEGRAL (0 TO INFINITY) OF THE FUNCTION
!+                GIVEN BY 'FONCRO', USING GAUSS QUADRATURES.
!
!history  F. BECQ (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)
!+        23/06/2011
!+        V6P1
!+   Translation of French names of the variables in argument
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| A              |-->| PARAMETER A OF THE FUNCTION TO BE INTEGRATED
!| B              |-->| PARAMETER B OF THE FUNCTION TO BE INTEGRATED
!| N              |-->| EXPONENT N OF THE FUNCTION TO BE INTEGRATED
!| XM             |-->| PARAMETER M OF THE FUNCTION TO BE INTEGRATED
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE interface_tomawac, ex_qgauss => qgauss
      IMPLICIT NONE
!
!     VARIABLES IN ARGUMENT
!     """""""""""""""""""""
      INTEGER, INTENT(IN)    ::  N
      DOUBLE PRECISION QGAUSS
      DOUBLE PRECISION, INTENT(IN)    :: B , A , XM
!
!     LOCAL VARIABLES
!     """"""""""""""""""
      INTEGER  J     , I     , NFOIS
      DOUBLE PRECISION :: XB    , XR    , DX    , DA    , SS    , W(5)
      DOUBLE PRECISION :: A1    , A2    , A3    , Y2    , X(5)
      parameter( x = (/ .1488743389d0,.4333953941d0,.6794095682d0,
     &                   .8650633666d0,.9739065285d0 /) )
      parameter( w = (/ .2955242247d0,.2692667193d0,.2190863625d0,
     &                   .1494513491d0,.0666713443d0 /) )
!
      nfois = 1
!
      CALL bornes
     &( b     , n     , a     , xm    , a2    , a3    )
      qgauss = 0.d0
      da = (a3-a2)/dble(nfois)
!
      DO i=1,nfois
        a1 = a2
        a2 = a2+da
        xb = 0.5d0*(a1+a2)
        xr = 0.5d0*(a2-a1)
        ss = 0.d0
        DO j=1,5
          dx = xr*x(j)
          ss = ss + w(j)*(foncro(xb+dx,b,n,a,xm)
     &                   +foncro(xb-dx,b,n,a,xm))
        ENDDO
        y2 = xr*ss
        qgauss = qgauss + y2
      ENDDO
!
      RETURN
      END