C***************************************************************
C       Subroutine VQUAD by boundary-element-method.com                  
C***************************************************************
C
C VQUAD.FOR Version 1, Sept 2015
C by Stephen Kirkup
C  School of Engineering
C  University of Central Lancashire
C  http://www.researchgate.net/profile/Stephen_Kirkup

C A test on the quadrature rule input to *2LC (and *ALC) routines in 
C the 2D and axisymmetric boundary element method.
C Source of thiS code:
C  http://www.boundary-element-methods.com/fortran/VQUAD.FOR
C Source of the user-guide:
C  http://www.boundary-element-methods.com/fortran/VQUAD_FOR.pdf
C Source of the test code:
C  http://www.boundary-element-methods.com/fortran/VQUAD_TESTS.FOR

      LOGICAL FUNCTION VQUAD(MAXN, N, W, X, XSTAR, TOLQUAD, IU)
      INTEGER MAXN
      INTEGER N
      REAL*8  W(MAXN)
      REAL*8  X(MAXN)
      REAL*8  XSTAR
      REAL*8  TOLQUAD
      INTEGER IU
      
      LOGICAL LERROR
      REAL*8  SUM
      INTEGER I

      VQUAD=.TRUE.
      
C Test 1: Verify tolQuad
C  Verify that tolQuad is positive
C   Tests 1(a)
      IF (TOLQUAD.LE.0.0D0) THEN
        VQUAD=.FALSE.
        WRITE(IU,*) 'Error(*LC)- Test 1(a) : tolerance (tolQuad) must'
        WRITE(IU,*) ' be positive'
      END IF
C  Verify that tolQuad is not too large       
C   Tests 1(b)
      IF (TOLQUAD.GT.0.001D0) THEN
        VQUAD=.FALSE.
        WRITE(IU,*) 'Error(*LC)- Test 1(b) : tolerance (tolQuad) is'
        WRITE(IU,*) ' too large'
      END IF

C Test 2: 
C  Test 2(a) Verify N is positive
      IF (N.LE.0) THEN
        VQUAD=.FALSE.
        WRITE(IU,*) 'Error(*LC)- Test 2(a) : number of quadrature '
        WRITE(IU,*) ' points (N) must be positive'
      END IF
C  Test 2(b) Verify MAXN is positive
      IF (MAXN.LE.0) THEN
        VQUAD=.FALSE.
        WRITE(IU,*) 'Error(*LC)- Test 2(b) : max number of quadrature '
        WRITE(IU,*) ' points (MAXN) must be positive'
      END IF
C  Test 2(c) Verify N<=MAXN is positive
      IF (N.GT.MAXN) THEN
        VQUAD=.FALSE.
        WRITE(IU,*) 'Error(*LC)- Test 2(c) : number of quadrature '
        WRITE(IU,*) ' points (N) must be at most the max (MAXN)'
      END IF          

C  Test 3: Verify that all of the abscissae are in [0,1]
      LERROR=.FALSE   .
      DO 100 I=1,N
        IF (X(I).LT.0.0D0) LERROR=.TRUE.
        IF (X(I).GT.1.0D0) LERROR=.TRUE.
100   CONTINUE
      IF (LERROR) THEN
        VQUAD=.FALSE.
        WRITE(IU,*) 'Error(*LC)- Test 3 : quadrature abscissae must'
        WRITE(IU,*) ' be on [0,1]'
      END IF           
        
C Test 4: Verify that all of the weights sum to one
      SUM=0.0D0
      DO 110 I=1,N
        SUM=SUM+W(I)
110   CONTINUE
      IF (ABS(SUM-1.0D0).GT.N*TOLQUAD) THEN
        VQUAD=.FALSE.
        WRITE(IU,*) 'Error(*LC)- Test 4 : quadrature weights must'
        WRITE(IU,*) ' sum to one'
      END IF          
           
C Test 5: Verify that the integral of f(x)=x is 0.5
      SUM=0.0D0
      DO 120 I=1,N
        SUM=SUM+W(I)*X(I)
120   CONTINUE
      IF (ABS(SUM-0.5D0).GT.N*TOLQUAD) THEN
        VQUAD=.FALSE.
        WRITE(IU,*) 'Error(*LC)- Test 5 : integration of f(x)=x failed'
      END IF          

        
C  The following tests are particularly for the case in which there is
C   a possible discontinuity at xStar
      IF (XSTAR.GT.0.0D0.AND.XSTAR.LE.1.0D0) THEN

C Test D1: Verify for f(x)=1 (x<xStar), f(x)=-1 (x>xStar)
        SUM=0.0D0
        DO 130 I=1,N
          IF (X(I).LE.XSTAR) THEN
            SUM=SUM+W(I)
          ELSE 
            SUM=SUM-W(I)
          END IF
130     CONTINUE
        IF (ABS(SUM-1.0D0+2.0D0*XSTAR).GT.N*TOLQUAD) THEN
          VQUAD=.FALSE.
          WRITE(IU,*) 'Error(*LC) - Test D1 : split quadrature rule'         
          WRITE(IU,*) ' must integrate function with discontinuity at'
          WRITE(IU,*) ' x_discontinutiy accurately'
        END IF

C Test D2: Verify for f(x)=1 (x<xStar), f(x)=-1 (x>xStar)
        SUM=0.0D0
        DO 140 I=1,N
          IF (X(I).LE.XSTAR) THEN
            SUM=SUM+W(I)
          ELSE 
            SUM=SUM-W(I)*0.5D0
          END IF
140     CONTINUE
        IF (ABS(SUM+0.5D0-1.5D0*XSTAR).GT.N*TOLQUAD) THEN
          VQUAD=.FALSE.
          WRITE(IU,*) 'Error(*LC) - Test D2 : split quadrature rule'         
          WRITE(IU,*) ' must integrate function with discontinuity at'
          WRITE(IU,*) ' x_discontinutiy accurately'
        END IF  

C Test D3: Verify for f(x)=x (x<xStar), f(x)=1 (x>xStar)
        SUM=0.0D0
        DO 150 I=1,N
          IF (X(I).LE.XSTAR) THEN
            SUM=SUM+W(I)*X(I)
          ELSE
            SUM=SUM+W(I)
          END IF
150     CONTINUE
        IF (ABS(SUM-1.0D0+XSTAR-0.5D0*XSTAR**2).GT.N*TOLQUAD) THEN
          VQUAD=.FALSE.
          WRITE(IU,*) 'Error(*LC) - Test D3 : split quadrature rule'         
          WRITE(IU,*) ' must integrate function with discontinuity at'
          WRITE(IU,*) ' x_discontinutiy accurately'
        END IF     
      END IF
      END