gquad6

   ==== Six Point Composite Gauss Formula ====
   ====   With Weight Factors Included    ====

  area = gquad6(fun,xlow,xhigh,mparts)
  This function determines the area under an externally
  defined function fun(x) between limits xlow and xhigh. The
  numerical integration is performed using a composite gauss
  integration rule.  The whole interval is divided into mparts
  subintervals and the integration over each subinterval
  is done with a six point Gauss formula which involves base
  points bp and weight factors wf.  The normalized interval
  of integration for the bp and wf constants is -1 to +1.  the
  algorithm is structured in terms of a parameter mquad = 6 which
  can be changed along with bp and wf to accommodate a different
  order formula.  The composite algorithm is described by the
  following summation relation
  x=b                     j=n k=m
  integral( f(x)*dx ) = d1*sum sum( wf(j)*fun(a1+d*k+d1*bp(j)) )
  x=a                     j=1 k=1
        where d = (b-a)/m, d1 = d/2, a1 = a-d1,
              m = mparts, and n = nquad.