|
Packit |
67cb25 |
/* integration/qk21.c
|
|
Packit |
67cb25 |
*
|
|
Packit |
67cb25 |
* Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Brian Gough
|
|
Packit |
67cb25 |
*
|
|
Packit |
67cb25 |
* This program is free software; you can redistribute it and/or modify
|
|
Packit |
67cb25 |
* it under the terms of the GNU General Public License as published by
|
|
Packit |
67cb25 |
* the Free Software Foundation; either version 3 of the License, or (at
|
|
Packit |
67cb25 |
* your option) any later version.
|
|
Packit |
67cb25 |
*
|
|
Packit |
67cb25 |
* This program is distributed in the hope that it will be useful, but
|
|
Packit |
67cb25 |
* WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
Packit |
67cb25 |
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
|
|
Packit |
67cb25 |
* General Public License for more details.
|
|
Packit |
67cb25 |
*
|
|
Packit |
67cb25 |
* You should have received a copy of the GNU General Public License
|
|
Packit |
67cb25 |
* along with this program; if not, write to the Free Software
|
|
Packit |
67cb25 |
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
#include <config.h>
|
|
Packit |
67cb25 |
#include <gsl/gsl_integration.h>
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Gauss quadrature weights and kronrod quadrature abscissae and
|
|
Packit |
67cb25 |
weights as evaluated with 80 decimal digit arithmetic by
|
|
Packit |
67cb25 |
L. W. Fullerton, Bell Labs, Nov. 1981. */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double xgk[11] = /* abscissae of the 21-point kronrod rule */
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
0.995657163025808080735527280689003,
|
|
Packit |
67cb25 |
0.973906528517171720077964012084452,
|
|
Packit |
67cb25 |
0.930157491355708226001207180059508,
|
|
Packit |
67cb25 |
0.865063366688984510732096688423493,
|
|
Packit |
67cb25 |
0.780817726586416897063717578345042,
|
|
Packit |
67cb25 |
0.679409568299024406234327365114874,
|
|
Packit |
67cb25 |
0.562757134668604683339000099272694,
|
|
Packit |
67cb25 |
0.433395394129247190799265943165784,
|
|
Packit |
67cb25 |
0.294392862701460198131126603103866,
|
|
Packit |
67cb25 |
0.148874338981631210884826001129720,
|
|
Packit |
67cb25 |
0.000000000000000000000000000000000
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* xgk[1], xgk[3], ... abscissae of the 10-point gauss rule.
|
|
Packit |
67cb25 |
xgk[0], xgk[2], ... abscissae to optimally extend the 10-point gauss rule */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double wg[5] = /* weights of the 10-point gauss rule */
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
0.066671344308688137593568809893332,
|
|
Packit |
67cb25 |
0.149451349150580593145776339657697,
|
|
Packit |
67cb25 |
0.219086362515982043995534934228163,
|
|
Packit |
67cb25 |
0.269266719309996355091226921569469,
|
|
Packit |
67cb25 |
0.295524224714752870173892994651338
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double wgk[11] = /* weights of the 21-point kronrod rule */
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
0.011694638867371874278064396062192,
|
|
Packit |
67cb25 |
0.032558162307964727478818972459390,
|
|
Packit |
67cb25 |
0.054755896574351996031381300244580,
|
|
Packit |
67cb25 |
0.075039674810919952767043140916190,
|
|
Packit |
67cb25 |
0.093125454583697605535065465083366,
|
|
Packit |
67cb25 |
0.109387158802297641899210590325805,
|
|
Packit |
67cb25 |
0.123491976262065851077958109831074,
|
|
Packit |
67cb25 |
0.134709217311473325928054001771707,
|
|
Packit |
67cb25 |
0.142775938577060080797094273138717,
|
|
Packit |
67cb25 |
0.147739104901338491374841515972068,
|
|
Packit |
67cb25 |
0.149445554002916905664936468389821
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
void
|
|
Packit |
67cb25 |
gsl_integration_qk21 (const gsl_function * f, double a, double b,
|
|
Packit |
67cb25 |
double *result, double *abserr,
|
|
Packit |
67cb25 |
double *resabs, double *resasc)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
double fv1[11], fv2[11];
|
|
Packit |
67cb25 |
gsl_integration_qk (11, xgk, wg, wgk, fv1, fv2, f, a, b, result, abserr, resabs, resasc);
|
|
Packit |
67cb25 |
}
|