V3FIT
integral.f
1  FUNCTION integral (n, x, y1, y2)
2  USE stel_kinds
3  IMPLICIT NONE
4 !-----------------------------------------------
5 ! D u m m y A r g u m e n t s
6 !-----------------------------------------------
7  INTEGER, INTENT(IN) :: n
8  REAL(rprec), INTENT(IN), DIMENSION(n) :: y1, y2 , x
9  REAL(rprec), PARAMETER :: one = 1
10  REAL(rprec) :: h, integral
11 !------------------------------------------------------------------------------
12 ! WARNING: The integral formula assumes values on the half MESH,
13 ! that is, the independent variable is given from x(3/2) to x(M-1/2),
14 ! but the limits of the integral are x(1)=-1 and x(M)=1. N = M-1 is the
15 ! number of points on the half mesh.
16 ! In the current version, fulfilment of this condition is ONLY checked
17 ! for the case in which the interval of integration is the interval
18 ! of orthogonality of Legendre polynomias, i.e., [-1,1].
19 ! Integration formula is a 2-nd order HALF-MESH formula from
20 ! "Numerical Recipes", W.Press et al, Chapter 4, page 110.
21 !------------------------------------------------------------------------------
22  IF (n < 10) stop 'Too few points to carry out integration.!'
23  IF (x(n) < x(1)) stop ' B < A in INTEGRAL!'
24  IF (x(n) == one .or. x(1) == -one) stop 'HALF MESH INTEGRAL!'
25 
26  h = (x(n)-x(1))/(n-1)
27  integral = h* sum(y1*y2)
28 
29  END FUNCTION integral