|
Packit |
67cb25 |
/* integration/qk31.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[16] = /* abscissae of the 31-point kronrod rule */
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
0.998002298693397060285172840152271,
|
|
Packit |
67cb25 |
0.987992518020485428489565718586613,
|
|
Packit |
67cb25 |
0.967739075679139134257347978784337,
|
|
Packit |
67cb25 |
0.937273392400705904307758947710209,
|
|
Packit |
67cb25 |
0.897264532344081900882509656454496,
|
|
Packit |
67cb25 |
0.848206583410427216200648320774217,
|
|
Packit |
67cb25 |
0.790418501442465932967649294817947,
|
|
Packit |
67cb25 |
0.724417731360170047416186054613938,
|
|
Packit |
67cb25 |
0.650996741297416970533735895313275,
|
|
Packit |
67cb25 |
0.570972172608538847537226737253911,
|
|
Packit |
67cb25 |
0.485081863640239680693655740232351,
|
|
Packit |
67cb25 |
0.394151347077563369897207370981045,
|
|
Packit |
67cb25 |
0.299180007153168812166780024266389,
|
|
Packit |
67cb25 |
0.201194093997434522300628303394596,
|
|
Packit |
67cb25 |
0.101142066918717499027074231447392,
|
|
Packit |
67cb25 |
0.000000000000000000000000000000000
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* xgk[1], xgk[3], ... abscissae of the 15-point gauss rule.
|
|
Packit |
67cb25 |
xgk[0], xgk[2], ... abscissae to optimally extend the 15-point gauss rule */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double wg[8] = /* weights of the 15-point gauss rule */
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
0.030753241996117268354628393577204,
|
|
Packit |
67cb25 |
0.070366047488108124709267416450667,
|
|
Packit |
67cb25 |
0.107159220467171935011869546685869,
|
|
Packit |
67cb25 |
0.139570677926154314447804794511028,
|
|
Packit |
67cb25 |
0.166269205816993933553200860481209,
|
|
Packit |
67cb25 |
0.186161000015562211026800561866423,
|
|
Packit |
67cb25 |
0.198431485327111576456118326443839,
|
|
Packit |
67cb25 |
0.202578241925561272880620199967519
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double wgk[16] = /* weights of the 31-point kronrod rule */
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
0.005377479872923348987792051430128,
|
|
Packit |
67cb25 |
0.015007947329316122538374763075807,
|
|
Packit |
67cb25 |
0.025460847326715320186874001019653,
|
|
Packit |
67cb25 |
0.035346360791375846222037948478360,
|
|
Packit |
67cb25 |
0.044589751324764876608227299373280,
|
|
Packit |
67cb25 |
0.053481524690928087265343147239430,
|
|
Packit |
67cb25 |
0.062009567800670640285139230960803,
|
|
Packit |
67cb25 |
0.069854121318728258709520077099147,
|
|
Packit |
67cb25 |
0.076849680757720378894432777482659,
|
|
Packit |
67cb25 |
0.083080502823133021038289247286104,
|
|
Packit |
67cb25 |
0.088564443056211770647275443693774,
|
|
Packit |
67cb25 |
0.093126598170825321225486872747346,
|
|
Packit |
67cb25 |
0.096642726983623678505179907627589,
|
|
Packit |
67cb25 |
0.099173598721791959332393173484603,
|
|
Packit |
67cb25 |
0.100769845523875595044946662617570,
|
|
Packit |
67cb25 |
0.101330007014791549017374792767493
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
void
|
|
Packit |
67cb25 |
gsl_integration_qk31 (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[16], fv2[16];
|
|
Packit |
67cb25 |
gsl_integration_qk (16, xgk, wg, wgk, fv1, fv2, f, a, b, result, abserr, resabs, resasc);
|
|
Packit |
67cb25 |
}
|