inte_gauss56_cern.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2020, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 /** \file inte_gauss56_cern.h
24  \brief File defining \ref o2scl::inte_gauss56_cern
25 */
26 #ifndef O2SCL_INTE_GAUSS56_CERN_H
27 #define O2SCL_INTE_GAUSS56_CERN_H
28 
29 #ifdef O2SCL_LD_TYPES
30 #include <boost/multiprecision/cpp_dec_float.hpp>
31 #endif
32 
33 #include <o2scl/misc.h>
34 #include <o2scl/inte.h>
35 #include <o2scl/funct.h>
36 
37 #ifndef DOXYGEN_NO_O2NS
38 namespace o2scl {
39 #endif
40 
41  /** \brief Integration weights and abcissas for
42  \ref o2scl::inte_gauss56_cern in double precision
43  */
45 
46  public:
47 
48  /** \brief Fifth order integration abscissas for
49  \ref o2scl::inte_gauss56_cern in double precision
50  */
51  double x5[5];
52 
53  /** \brief Fifth order integration weights for
54  \ref o2scl::inte_gauss56_cern in double precision
55  */
56  double w5[5];
57 
58  /** \brief Sixth order integration abscissas for
59  \ref o2scl::inte_gauss56_cern in double precision
60  */
61  double x6[6];
62 
63  /** \brief Sixth order integration weights for
64  \ref o2scl::inte_gauss56_cern in double precision
65  */
66  double w6[6];
67 
69  x5[0]=4.6910077030668004e-02;
70  x5[1]=2.3076534494715846e-01;
71  x5[2]=5.0000000000000000e-01;
72  x5[3]=7.6923465505284154e-01;
73  x5[4]=9.5308992296933200e-01;
74 
75  w5[0]=1.1846344252809454e-01;
76  w5[1]=2.3931433524968324e-01;
77  w5[2]=2.8444444444444444e-01;
78  w5[3]=2.3931433524968324e-01;
79  w5[4]=1.1846344252809454e-01;
80 
81  x6[0]=3.3765242898423989e-02;
82  x6[1]=1.6939530676686775e-01;
83  x6[2]=3.8069040695840155e-01;
84  x6[3]=6.1930959304159845e-01;
85  x6[4]=8.3060469323313225e-01;
86  x6[5]=9.6623475710157601e-01;
87 
88  w6[0]=8.5662246189585178e-02;
89  w6[1]=1.8038078652406930e-01;
90  w6[2]=2.3395696728634552e-01;
91  w6[3]=2.3395696728634552e-01;
92  w6[4]=1.8038078652406930e-01;
93  w6[5]=8.5662246189585178e-02;
94  }
95  };
96 
97  /** \brief Integration weights and abcissas for
98  \ref o2scl::inte_gauss56_cern in long double precision
99 
100  \note The long double type doesn't work uniformly across systems
101  and so the accuracy when using these coefficients varies.
102  */
104 
105  public:
106 
107  /** \brief Fifth order integration abscissas for
108  \ref o2scl::inte_gauss56_cern in long double precision
109  */
110  long double x5[5];
111 
112  /** \brief Fifth order integration weights for
113  \ref o2scl::inte_gauss56_cern in long double precision
114  */
115  long double w5[5];
116 
117  /** \brief Sixth order integration abscissas for
118  \ref o2scl::inte_gauss56_cern in long double precision
119  */
120  long double x6[6];
121 
122  /** \brief Sixth order integration weights for
123  \ref o2scl::inte_gauss56_cern in long double precision
124  */
125  long double w6[6];
126 
128 
129  x5[0]=0.04691007703066800360118656085030352L;
130  x5[1]=0.23076534494715845448184278964989560L;
131  x5[2]=0.5L;
132  x5[3]=0.76923465505284154551815721035010440L;
133  x5[4]=0.95308992296933199639881343914969648L;
134 
135  w5[0]=0.11846344252809454375713202035995868L;
136  w5[1]=0.23931433524968323402064575741781910L;
137  w5[2]=0.28444444444444444444444444444444444L;
138  w5[3]=0.23931433524968323402064575741781910L;
139  w5[4]=0.11846344252809454375713202035995868L;
140 
141  x6[0]=0.03376524289842398609384922275300270L;
142  x6[1]=0.16939530676686774316930020249004733L;
143  x6[2]=0.38069040695840154568474913915964403L;
144  x6[3]=0.61930959304159845431525086084035597L;
145  x6[4]=0.83060469323313225683069979750995267L;
146  x6[5]=0.96623475710157601390615077724699730L;
147 
148  w6[0]=0.08566224618958517252014807108636645L;
149  w6[1]=0.18038078652406930378491675691885806L;
150  w6[2]=0.23395696728634552369493517199477550L;
151  w6[3]=0.23395696728634552369493517199477550L;
152  w6[4]=0.18038078652406930378491675691885806L;
153  w6[5]=0.08566224618958517252014807108636645L;
154 
155  }
156 
157  };
158 
159 #if defined(O2SCL_LD_TYPES) || defined(DOXYGEN)
160 
161  /** \brief Integration weights and abcissas for
162  \ref o2scl::inte_gauss56_cern in long double precision
163 
164  \note Experimental, and only included if
165  O2SCL_LD_TYPES is defined during the library configuration.
166 
167  \comment
168  Weights and abcissas originally generated using cpp_dec_float_100
169  numbers by AWS using code in ~/wcs/int5/sbox on 10/7/19.
170  \endcomment
171  */
173 
174  public:
175 
176  /** \brief Fifth order integration abscissas for
177  \ref o2scl::inte_gauss56_cern in cpp_dec_float_50 precision
178  */
179  boost::multiprecision::cpp_dec_float_50 x5[5];
180 
181  /** \brief Fifth order integration weights for
182  \ref o2scl::inte_gauss56_cern in cpp_dec_float_50 precision
183  */
184  boost::multiprecision::cpp_dec_float_50 w5[5];
185 
186  /** \brief Sixth order integration abscissas for
187  \ref o2scl::inte_gauss56_cern in cpp_dec_float_50 precision
188  */
189  boost::multiprecision::cpp_dec_float_50 x6[6];
190 
191  /** \brief Sixth order integration weights for
192  \ref o2scl::inte_gauss56_cern in cpp_dec_float_50 precision
193  */
194  boost::multiprecision::cpp_dec_float_50 w6[6];
195 
197 
198  x5[0]=boost::multiprecision::cpp_dec_float_50
199  ("4.69100770306680036011865608503035174371740446187346e-02");
200  w5[0]=boost::multiprecision::cpp_dec_float_50
201  ("1.18463442528094543757132020359958681321630001106207e-01");
202  x5[1]=boost::multiprecision::cpp_dec_float_50
203  ("2.30765344947158454481842789649895597516356696547220e-01");
204  w5[1]=boost::multiprecision::cpp_dec_float_50
205  ("2.39314335249683234020645757417819096456147776671571e-01");
206  x5[2]=boost::multiprecision::cpp_dec_float_50
207  ("5.00000000000000000000000000000000000000000000000000e-01");
208  w5[2]=boost::multiprecision::cpp_dec_float_50
209  ("2.84444444444444444444444444444444444444444444444444e-01");
210  x5[3]=boost::multiprecision::cpp_dec_float_50
211  ("7.69234655052841545518157210350104402483643303452780e-01");
212  w5[3]=boost::multiprecision::cpp_dec_float_50
213  ("2.39314335249683234020645757417819096456147776671571e-01");
214  x5[4]=boost::multiprecision::cpp_dec_float_50
215  ("9.53089922969331996398813439149696482562825955381265e-01");
216  w5[4]=boost::multiprecision::cpp_dec_float_50
217  ("1.18463442528094543757132020359958681321630001106207e-01");
218 
219  x6[0]=boost::multiprecision::cpp_dec_float_50
220  ("3.37652428984239860938492227530026954326171311438551e-02");
221  w6[0]=boost::multiprecision::cpp_dec_float_50
222  ("8.56622461895851725201480710863664467634112507420220e-02");
223  x6[1]=boost::multiprecision::cpp_dec_float_50
224  ("1.69395306766867743169300202490047326496775717802415e-01");
225  w6[1]=boost::multiprecision::cpp_dec_float_50
226  ("1.80380786524069303784916756918858055830760946373373e-01");
227  x6[2]=boost::multiprecision::cpp_dec_float_50
228  ("3.80690406958401545684749139159644032290694684929989e-01");
229  w6[2]=boost::multiprecision::cpp_dec_float_50
230  ("2.33956967286345523694935171994775497405827802884605e-01");
231  x6[3]=boost::multiprecision::cpp_dec_float_50
232  ("6.19309593041598454315250860840355967709305315070011e-01");
233  w6[3]=boost::multiprecision::cpp_dec_float_50
234  ("2.33956967286345523694935171994775497405827802884605e-01");
235  x6[4]=boost::multiprecision::cpp_dec_float_50
236  ("8.30604693233132256830699797509952673503224282197585e-01");
237  w6[4]=boost::multiprecision::cpp_dec_float_50
238  ("1.80380786524069303784916756918858055830760946373373e-01");
239  x6[5]=boost::multiprecision::cpp_dec_float_50
240  ("9.66234757101576013906150777246997304567382868856145e-01");
241  w6[5]=boost::multiprecision::cpp_dec_float_50
242  ("8.56622461895851725201480710863664467634112507420220e-02");
243 
244  }
245 
246  };
247 
248 #endif
249 
250  /** \brief 5,6-point Gaussian quadrature (CERNLIB)
251 
252  If \f$ I_5 \f$ is the 5-point approximation, and \f$ I_6 \f$ is the
253  6-point approximation to the integral, then integ_err() returns
254  the result \f$ \frac{1}{2}(I_5+I_6) \f$ with uncertainty
255  \f$ |I_5-I_6| \f$.
256 
257  This class is based on the CERNLIB routines RGS56P and
258  DGS56P which are documented at
259  http://wwwasdoc.web.cern.ch/wwwasdoc/shortwrupsdir/d106/top.html
260  */
261  template<class func_t=funct, class fp_t=double,
262  class weights_t=inte_gauss56_coeffs_double>
263  class inte_gauss56_cern : public inte<func_t,fp_t> {
264 
265  protected:
266 
267  const fp_t *w5, *x5, *w6, *x6;
268  weights_t wgts;
269 
270  public:
271 
273  w5=&(wgts.w5[0]);
274  x5=&(wgts.x5[0]);
275  w6=&(wgts.w6[0]);
276  x6=&(wgts.x6[0]);
277  }
278 
279  /** \brief Integrate function \c func from \c a to \c b
280  giving result \c res and error \c err
281 
282  This function always returns \ref success.
283  */
284  virtual int integ_err(func_t &func, fp_t a, fp_t b,
285  fp_t &res, fp_t &err) {
286 
287  fp_t rang=b-a, e5=0.0, e6=0.0, ytmp;
288 
289  for(int i=0;i<5;i++) {
290  ytmp=func(a+rang*x5[i]);
291  e5+=w5[i]*ytmp;
292  ytmp=func(a+rang*x6[i]);
293  e6+=w6[i]*ytmp;
294  }
295  ytmp=func(a+rang*x6[5]);
296  e6+=w6[5]*ytmp;
297  res=(e6+e5)*rang/2.0;
298  err=o2scl::o2abs(e5-e6)*rang;
299 
300  return success;
301  }
302 
303  };
304 
305 #ifndef DOXYGEN_NO_O2NS
306 }
307 #endif
308 
309 #endif
o2scl::inte_gauss56_coeffs_cpp_dec_float_50
Integration weights and abcissas for o2scl::inte_gauss56_cern in long double precision.
Definition: inte_gauss56_cern.h:172
o2scl
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
Definition: anneal.h:42
o2scl::inte_gauss56_coeffs_double::x6
double x6[6]
Sixth order integration abscissas for o2scl::inte_gauss56_cern in double precision.
Definition: inte_gauss56_cern.h:61
o2scl::inte_gauss56_cern
5,6-point Gaussian quadrature (CERNLIB)
Definition: inte_gauss56_cern.h:263
o2scl::inte_gauss56_coeffs_cpp_dec_float_50::x6
boost::multiprecision::cpp_dec_float_50 x6[6]
Sixth order integration abscissas for o2scl::inte_gauss56_cern in cpp_dec_float_50 precision.
Definition: inte_gauss56_cern.h:189
o2scl::inte_gauss56_coeffs_long_double::x6
long double x6[6]
Sixth order integration abscissas for o2scl::inte_gauss56_cern in long double precision.
Definition: inte_gauss56_cern.h:120
o2scl::success
@ success
Success.
Definition: err_hnd.h:47
o2scl::inte_gauss56_coeffs_double::w6
double w6[6]
Sixth order integration weights for o2scl::inte_gauss56_cern in double precision.
Definition: inte_gauss56_cern.h:66
o2scl::inte_gauss56_coeffs_double
Integration weights and abcissas for o2scl::inte_gauss56_cern in double precision.
Definition: inte_gauss56_cern.h:44
o2scl::inte_gauss56_coeffs_long_double::w6
long double w6[6]
Sixth order integration weights for o2scl::inte_gauss56_cern in long double precision.
Definition: inte_gauss56_cern.h:125
o2scl::inte_gauss56_coeffs_cpp_dec_float_50::x5
boost::multiprecision::cpp_dec_float_50 x5[5]
Fifth order integration abscissas for o2scl::inte_gauss56_cern in cpp_dec_float_50 precision.
Definition: inte_gauss56_cern.h:179
o2scl::inte_gauss56_cern::integ_err
virtual int integ_err(func_t &func, fp_t a, fp_t b, fp_t &res, fp_t &err)
Integrate function func from a to b giving result res and error err.
Definition: inte_gauss56_cern.h:284
o2scl::inte_gauss56_coeffs_long_double::w5
long double w5[5]
Fifth order integration weights for o2scl::inte_gauss56_cern in long double precision.
Definition: inte_gauss56_cern.h:115
o2scl::inte_gauss56_coeffs_double::x5
double x5[5]
Fifth order integration abscissas for o2scl::inte_gauss56_cern in double precision.
Definition: inte_gauss56_cern.h:51
o2scl::inte
Base integration class [abstract base].
Definition: inte.h:51
o2scl::funct
std::function< double(double)> funct
One-dimensional function typedef in src/base/funct.h.
Definition: funct.h:48
o2scl::o2abs
float o2abs(const float x)
Absolute value for single precision numbers.
o2scl::inte_gauss56_coeffs_cpp_dec_float_50::w5
boost::multiprecision::cpp_dec_float_50 w5[5]
Fifth order integration weights for o2scl::inte_gauss56_cern in cpp_dec_float_50 precision.
Definition: inte_gauss56_cern.h:184
o2scl::inte_gauss56_coeffs_cpp_dec_float_50::w6
boost::multiprecision::cpp_dec_float_50 w6[6]
Sixth order integration weights for o2scl::inte_gauss56_cern in cpp_dec_float_50 precision.
Definition: inte_gauss56_cern.h:194
o2scl::inte_gauss56_coeffs_long_double::x5
long double x5[5]
Fifth order integration abscissas for o2scl::inte_gauss56_cern in long double precision.
Definition: inte_gauss56_cern.h:110
o2scl::inte_gauss56_coeffs_long_double
Integration weights and abcissas for o2scl::inte_gauss56_cern in long double precision.
Definition: inte_gauss56_cern.h:103
o2scl::inte_gauss56_coeffs_double::w5
double w5[5]
Fifth order integration weights for o2scl::inte_gauss56_cern in double precision.
Definition: inte_gauss56_cern.h:56

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).