Provided by: libquadrule-dev_0~20121001-2build1_amd64 bug

NAME

       quadrule ‐ C library of quadrature rules

SYNOPSIS

       #include <quadrule.h>

DESCRIPTION

       APPROXIMATING INTEGRALS

       The  following  functions  apply  an  existing  set of abscissas and weights to compute an
       approximation to an integral.

       double summer ( double func ( double x ), int order, double xtab[], double weight[] )

              Carries out a quadrature rule over a single interval.

       double sum_sub ( double func ( double x ), double a, double b, int nsub, int order, double
              xlo, double xhi, double xtab[], double weight[] )

              Carries out a composite quadrature rule.

              xlo  and  xhi  are  the  left  and  right  endpoints of the interval over which the
              quadrature rule was defined.

              a and b are the lower and upper limits of integration.

              nsub is the number of equal subintervals into which the finite interval (a,b) is to
              be subdivided for higher accuracy.  nsub must be at least 1.

       double bdf_sum ( double func ( double x ), int order, double xtab[], double weight[] )

              Carries  out  explicit  backward difference quadrature on [0,1].  Requires xtab and
              weight to be pre‐computed from bdfc_set or bdfp_set.

       double laguerre_sum ( double func ( double x ), double a, int order, double xtab[], double
              weight[] )

              Carries  out  Laguerre  quadrature over [ A, +oo ).  Requires xtab and weight to be
              pre‐computed from laguerre_compute, gen_laguerre_compute, or laguerre_set.

       void summer_gk ( double func ( double x ), int orderg, double weightg[], double  *resultg,
              int orderk , double xtabk[], double weightk[], double *resultk )

              Carries  out  Gauss‐Kronrod quadrature over a single interval.  Before calling this
              function, weightg should be pre‐computed from legendre_compute_dr or  legendre_set.
              Also, xtabk and weightk should be pre‐computed from kronrod_set.  The orders should
              follow the relation orderk = 2 * orderg + 1.  resultk will contain  the  result  of
              the  Gauss‐Kronrod sum.  resultg contains an intermediate sum, using only the Gauss
              quadrature.

       void sum_sub_gk ( double func ( double x ), double a, double  b,  int  nsub,  int  orderg,
              double  weightg[],  double  *resultg, int orderk, double xtabk[], double weightk[],
              double *resultk, double *error )

              Carries out a composite Gauss‐Kronrod rule.  Approximates the integral of func from
              a  to  b, by dividing the domain into nsub subdomains, and applying a Gauss‐Kronrod
              rule on each subdomain.

       The following are utility function(s), to  modify  existing  weights  /  abscissas  before
       computing an integral.

       void  rule_adjust  (  double  a,  double  b, double c, double d, int order, double xtab[],
              double weight[] )

              Maps a quadrature rule from [A,B] to [C,D]. xtab and weight are overwritten.

       COMPUTING QUADRATURE ABSCISSAS / WEIGHTS

       The following functions compute arrays of abscissas  and  weights  (xtab,  weight)  for  a
       particular  quadrature rule, for any given order.  The range of the abscissas is ( -1 <= x
       < 1 ) unless otherwise given.

       NOTE: xtab and weight must be allocated before calling these functions.

       void chebyshev1_compute ( int order, double xtab[], double weight[] )

              Computes a Gauss‐Chebyshev type 1 quadrature rule.

       void chebyshev2_compute ( int order, double xtab[], double weight[] )

              Computes a Gauss‐Chebyshev type 2 quadrature rule.

       void chebyshev3_compute ( int order, double xtab[], double weight[] )

              Computes a Gauss‐Chebyshev type 3 quadrature rule.

       void clenshaw_curtis_compute ( int order, double xtab[], double weight[] )

              Computes a Clenshaw Curtis quadrature rule.

       void fejer1_compute ( int order, double xtab[], double weight[] )

              Computes a Fejer type 1 quadrature rule.

       void fejer2_compute ( int order, double xtab[], double weight[] )

              Computes a Fejer type 2 quadrature rule.

       void gegenbauer_compute ( int order, double alpha, double xtab[], double weight[] )

              Computes a Gauss‐Gegenbauer quadrature rule.

       void gen_hermite_compute ( int order, double alpha, double xtab[], double weight[] )

              Computes a generalized Gauss‐Hermite rule for  the  interval  (  -infinity  <  x  <
              +infinity )

       void gen_laguerre_compute ( int order, double alpha, double xtab[], double weight[] )

              Computes a generalized Gauss‐Laguerre quadrature rule for the interval ( alpha <= x
              < infinity )

       void hermite_ek_compute ( int order, double xtab[], double weight[] )

              Computes a Gauss‐Hermite quadrature rule, using an algorithm by Elhay and  Kautsky.
              Interval is ( -infinity < x < infinity )

       void hermite_ss_compute ( int order, double xtab[], double weight[] )

              Computes  a  Gauss‐Hermite quadrature rule, using an algorithm by Arthur Stroud and
              Don Secrest. Interval is ( -infinity < x < infinity )

       void jacobi_compute ( int order, double alpha, double beta, double xtab[], double weight[]
              )

              Computes a Gauss‐Jacobi quadrature rule.

       void laguerre_compute ( int order, double xtab[], double weight[] )

              Computes a Gauss‐Laguerre quadrature rule for the inverval ( 0 <= x < infinity )

       void legendre_compute_dr ( int order, double xtab[], double weight[] )

              Gauss‐Legendre quadrature by Davis‐Rabinowitz method.

       void lobatto_compute ( int order, double xtab[], double weight[] )

              Computes a Lobatto quadrature rule.

       void nc_compute ( int order, double x_min, double x_max, double xtab[], double weight[] )

              Computes a Newton‐Cotes quadrature rule.

              ( x_min <= x <= x_max )

       void ncc_compute ( int order, double xtab[], double weight[] )

              Computes a Newton‐Cotes Closed quadrature rule.

       void nco_compute ( int order, double xtab[], double weight[] )

              Computes a Newton‐Cotes Open quadrature rule.

       void ncoh_compute ( int order, double xtab[], double weight[] )

              Computes a Newton‐Cotes "open half" quadrature rule.

              ( x_min <= x <= x_max )

       void radau_compute ( int order, double xtab[], double weight[] )

              Computes a Radau quadrature rule.

       PRE‐COMPUTED QUADRATURES

       The  following  functions  return abscissas/weights for various quadrature rules, but only
       for particular orders.

       NOTE: Check the valid range of order parameters for the particular function  before  using
       it.  Using an unsupported order will cause an abort.

       void bashforth_set ( int order, double xtab[], double weight [] )

              Sets  an  Adams‐Bashforth quadrature for order = 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 12,
              14, 16, 18, or 20.

       void bdfc_set ( int order, double xtab[], double weight[] )

              Sets weights for backward differentiation corrector quadrature, for 1 <=  order  <=
              19.

       void bdfp_set ( int order, double xtab[], double weight[] )

              Sets  weights  for backward differentiation predictor quadrature, for 1 <= order <=
              19.

       void cheb_set ( int order, double xtab[], double weight[] )

              Sets abscissas and weights for Chebyshev quadrature, for order = 1, 2, 3, 4, 5,  6,
              7, or 9.

       void clenshaw_curtis_set ( int order, double xtab[], double weight[] )

              Sets  a Clenshaw‐Curtis quadrature rule, for order = 1, 2, 3, 4, 5, 6, 7, 8, 9, 10,
              11, 12, 13, 14, 15, 16, 17, 33, 65, or 129.

       void fejer1_set ( int order, double xtab[], double weight[] )

              Sets abscissas and weights for Fejer type 1 quadrature, for order = 1, 2, 3, 4,  5,
              6, 7, 8, 9, or 10.

       void fejer2_set ( int order, double xtab[], double weight[] )

              Sets  abscissas and weights for Fejer type 2 quadrature, for order = 1, 2, 3, 4, 5,
              6, 7, 8, 9, or 10.

       void hermite_genz_keister_set ( int order, double xtab[], double weight[] )

              Sets a Hermite Genz‐Keister rule, for order = 1, 3, 7, 9, 17, 19, 31, 33, or 35.

       void hermite_set ( int order, double xtab[], double weight[] )

              Sets abscissas and weights for Hermite quadrature.  The order must be between 1 and
              20, or 31/32/33, 63/64/65, 127/128/129.

       void kronrod_set ( int order, double xtab[], double weight[] )

              Sets  abscissas and weights for Gauss‐Kronrod quadrature.  The order may be 15, 21,
              31 or 41, corresponding to Gauss‐Legendre rules of order 7, 10, 15 or 20.

       void laguerre_set ( int order, double xtab[], double weight[] )

              Sets abscissas and weights for Laguerre quadrature.  The order must  be  between  1
              and 20, or 31/32/33, 63/64/65, 127/128/129.

       void legendre_set ( int order, double xtab[], double weight[] )

              Sets  abscissas  and  weights  for  Gauss‐Legendre  quadrature.   The order must be
              between 1 and 33 or 63/64/65, 127/128/129, 255/256/257.

       void legendre_set_cos ( int order, double xtab[], double weight[] )

              Sets Gauss‐Legendre rules for COS(X)*F(X) on [-PI/2,PI/2].  The order must be 1, 2,
              4, 8 or 16.

       void legendre_set_cos2 ( int order, double xtab[], double weight[] )

              Sets  Gauss‐Legendre rules for COS(X)*F(X) on [0,PI/2].  The order must be  2, 4, 8
              or 16.

       void legendre_set_log ( int order, double xtab[], double weight[] )

              Sets a Gauss‐Legendre rule for - LOG(X) * F(X) on [0,1].  The order must be between
              1 through 8, or 16.

       void legendre_set_sqrtx_01 ( int order, double xtab[], double weight[] )

              Sets  Gauss‐Legendre  rules for SQRT(X)*F(X) on [0,1].  The order must be between 1
              and 16 or 63/64, 127/128.

       void legendre_set_sqrtx2_01 ( int order, double xtab[], double weight[] )

              Sets Gauss‐Legendre rules for F(X)/SQRT(X) on [0,1].  The order must be be tween  1
              and 16 or 63/64, 127/128.

       void legendre_set_x0_01 ( int order, double xtab[], double weight[] )

              Sets a Gauss‐Legendre rule for F(X) on [0,1].  The order must be between 1 and 8.

       void legendre_set_x1 ( int order, double xtab[], double weight[] )

              Sets  a  Gauss‐Legendre  rule  for  (  1 + X ) * F(X) on [-1,1].  The order must be
              between 1 and 9.

       void legendre_set_x1_01 ( int order, double xtab[], double weight[] )

              Sets a Gauss‐Legendre rule for X * F(X) on [0,1].  The order must be between 1  and
              8.

       void legendre_set_x2 ( int order, double xtab[], double weight[] )

              Sets  Gauss‐Legendre  rules  for  (  1  +  X )^2*F(X) on [-1,1].  The order must be
              between 1 and 9.

       void legendre_set_x2_01 ( int order, double xtab[], double weight[] )

              Sets a Gauss‐Legendre rule for X*X * F(X) on [0,1].  The order must  be  between  1
              and 8.

       void lobatto_set ( int order, double xtab[], double weight[] )

              Sets abscissas and weights for Lobatto quadrature.  The order must be between 2 and
              20.

       void moulton_set ( int order, double xtab[], double weight[] )

              Sets weights for Adams‐Moulton quadrature.  The order must be between 1 and  10  or
              12, 14, 16, 18 or 20.

       void ncc_set ( int order, double xtab[], double weight[] )

              Sets  abscissas  and weights for closed Newton‐Cotes quadrature.  The order must be
              between 1 and 21.

       void nco_set ( int order, double xtab[], double weight[] )

              Sets abscissas and weights for open Newton‐Cotes quadrature.   The  order  must  be
              between 1 and 7, or 9.

       void ncoh_set ( int order, double xtab[], double weight[] )

              Sets  abscissas  and weights for Newton‐Cotes "open half" rules.  The order must be
              between 1 and 17.

       void patterson_set ( int order, double xtab[], double weight[] )

              Sets abscissas and weights for Gauss‐Patterson quadrature.  The order must be 1, 3,
              7, 15, 31, 63, 127 or 255.

       void radau_set ( int order, double xtab[], double weight[] )

              Sets  abscissas  and weights for Radau quadrature.  The order must be between 1 and
              15.

AUTHOR

       John Burkardt

SEE ALSO

       The official quadrule web page
       http://people.sc.fsu.edu/~jburkardt/c_src/quadrule/quadrule.html

                                           15 Mar 2013                                quadrule(3)