1 ////////////////////////////////////////////////////////////////////////
2 //
3 // Copyright (C) 1993-2021 The Octave Project Developers
4 //
5 // See the file COPYRIGHT.md in the top-level directory of this
6 // distribution or <https://octave.org/copyright/>.
7 //
8 // This file is part of Octave.
9 //
10 // Octave is free software: you can redistribute it and/or modify it
11 // under the terms of the GNU General Public License as published by
12 // the Free Software Foundation, either version 3 of the License, or
13 // (at your option) any later version.
14 //
15 // Octave is distributed in the hope that it will be useful, but
16 // WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
18 // GNU General Public License for more details.
19 //
20 // You should have received a copy of the GNU General Public License
21 // along with Octave; see the file COPYING.  If not, see
22 // <https://www.gnu.org/licenses/>.
23 //
24 ////////////////////////////////////////////////////////////////////////
25 
26 #if ! defined (octave_Quad_h)
27 #define octave_Quad_h 1
28 
29 #include "octave-config.h"
30 
31 #include "dColVector.h"
32 #include "fColVector.h"
33 
34 typedef double (*integrand_fcn) (double x);
35 typedef float (*float_integrand_fcn) (float x);
36 
37 #include "Quad-opts.h"
38 
39 class
40 OCTAVE_API
41 Quad : public Quad_options
42 {
43 public:
44 
Quad(integrand_fcn fcn)45   Quad (integrand_fcn fcn)
46     : Quad_options (), f (fcn), ff () { }
47 
Quad(float_integrand_fcn fcn)48   Quad (float_integrand_fcn fcn)
49     : Quad_options (), f (), ff (fcn) { }
50 
51   virtual ~Quad (void) = default;
52 
integrate(void)53   virtual double integrate (void)
54   {
55     octave_idx_type ier, neval;
56     double abserr;
57     return do_integrate (ier, neval, abserr);
58   }
59 
float_integrate(void)60   virtual float float_integrate (void)
61   {
62     octave_idx_type ier, neval;
63     float abserr;
64     return do_integrate (ier, neval, abserr);
65   }
66 
integrate(octave_idx_type & ier)67   virtual double integrate (octave_idx_type& ier)
68   {
69     octave_idx_type neval;
70     double abserr;
71     return do_integrate (ier, neval, abserr);
72   }
73 
float_integrate(octave_idx_type & ier)74   virtual float float_integrate (octave_idx_type& ier)
75   {
76     octave_idx_type neval;
77     float abserr;
78     return do_integrate (ier, neval, abserr);
79   }
80 
integrate(octave_idx_type & ier,octave_idx_type & neval)81   virtual double integrate (octave_idx_type& ier, octave_idx_type& neval)
82   {
83     double abserr;
84     return do_integrate (ier, neval, abserr);
85   }
86 
float_integrate(octave_idx_type & ier,octave_idx_type & neval)87   virtual float float_integrate (octave_idx_type& ier, octave_idx_type& neval)
88   {
89     float abserr;
90     return do_integrate (ier, neval, abserr);
91   }
92 
integrate(octave_idx_type & ier,octave_idx_type & neval,double & abserr)93   virtual double integrate (octave_idx_type& ier, octave_idx_type& neval,
94                             double& abserr)
95   {
96     return do_integrate (ier, neval, abserr);
97   }
98 
float_integrate(octave_idx_type & ier,octave_idx_type & neval,float & abserr)99   virtual float float_integrate (octave_idx_type& ier, octave_idx_type& neval,
100                                  float& abserr)
101   {
102     return do_integrate (ier, neval, abserr);
103   }
104 
105   virtual double do_integrate (octave_idx_type& ier, octave_idx_type& neval,
106                                double& abserr) = 0;
107 
108   virtual float do_integrate (octave_idx_type& ier, octave_idx_type& neval,
109                               float& abserr) = 0;
110 
111 protected:
112 
113   integrand_fcn f;
114   float_integrand_fcn ff;
115 };
116 
117 class
118 OCTAVE_API
119 DefQuad : public Quad
120 {
121 public:
122 
DefQuad(integrand_fcn fcn)123   DefQuad (integrand_fcn fcn)
124     : Quad (fcn), lower_limit (0.0), upper_limit (1.0), singularities () { }
125 
DefQuad(integrand_fcn fcn,double ll,double ul)126   DefQuad (integrand_fcn fcn, double ll, double ul)
127     : Quad (fcn), lower_limit (ll), upper_limit (ul), singularities () { }
128 
DefQuad(integrand_fcn fcn,double ll,double ul,const ColumnVector & sing)129   DefQuad (integrand_fcn fcn, double ll, double ul,
130            const ColumnVector& sing)
131     : Quad (fcn), lower_limit (ll), upper_limit (ul),
132       singularities (sing) { }
133 
DefQuad(integrand_fcn fcn,const ColumnVector & sing)134   DefQuad (integrand_fcn fcn, const ColumnVector& sing)
135     : Quad (fcn), lower_limit (0.0), upper_limit (1.0),
136       singularities (sing) { }
137 
138   ~DefQuad (void) = default;
139 
140   double do_integrate (octave_idx_type& ier, octave_idx_type& neval,
141                        double& abserr);
142 
143   OCTAVE_NORETURN float do_integrate (octave_idx_type& ier,
144                                       octave_idx_type& neval, float& abserr);
145 
146 private:
147 
148   double lower_limit;
149   double upper_limit;
150 
151   ColumnVector singularities;
152 };
153 
154 class
155 OCTAVE_API
156 IndefQuad : public Quad
157 {
158 public:
159 
160   enum IntegralType { bound_to_inf, neg_inf_to_bound, doubly_infinite };
161 
IndefQuad(integrand_fcn fcn)162   IndefQuad (integrand_fcn fcn)
163     : Quad (fcn), bound (0.0), type (bound_to_inf) { }
164 
IndefQuad(integrand_fcn fcn,double b,IntegralType t)165   IndefQuad (integrand_fcn fcn, double b, IntegralType t)
166     : Quad (fcn), bound (b), type (t) { }
167 
168   ~IndefQuad (void) = default;
169 
170   double do_integrate (octave_idx_type& ier, octave_idx_type& neval,
171                        double& abserr);
172 
173   OCTAVE_NORETURN float do_integrate (octave_idx_type& ier,
174                                       octave_idx_type& neval, float& abserr);
175 
176 private:
177 
178   double bound;
179   IntegralType type;
180 };
181 
182 class
183 OCTAVE_API
184 FloatDefQuad : public Quad
185 {
186 public:
187 
FloatDefQuad(float_integrand_fcn fcn)188   FloatDefQuad (float_integrand_fcn fcn)
189     : Quad (fcn), lower_limit (0.0), upper_limit (1.0), singularities () { }
190 
FloatDefQuad(float_integrand_fcn fcn,float ll,float ul)191   FloatDefQuad (float_integrand_fcn fcn, float ll, float ul)
192     : Quad (fcn), lower_limit (ll), upper_limit (ul), singularities () { }
193 
FloatDefQuad(float_integrand_fcn fcn,float ll,float ul,const FloatColumnVector & sing)194   FloatDefQuad (float_integrand_fcn fcn, float ll, float ul,
195                 const FloatColumnVector& sing)
196     : Quad (fcn), lower_limit (ll), upper_limit (ul),
197       singularities (sing) { }
198 
FloatDefQuad(float_integrand_fcn fcn,const FloatColumnVector & sing)199   FloatDefQuad (float_integrand_fcn fcn, const FloatColumnVector& sing)
200     : Quad (fcn), lower_limit (0.0), upper_limit (1.0),
201       singularities (sing) { }
202 
203   ~FloatDefQuad (void) = default;
204 
205   OCTAVE_NORETURN double do_integrate (octave_idx_type& ier,
206                                        octave_idx_type& neval, double& abserr);
207 
208   float do_integrate (octave_idx_type& ier, octave_idx_type& neval,
209                       float& abserr);
210 
211 private:
212 
213   float lower_limit;
214   float upper_limit;
215 
216   FloatColumnVector singularities;
217 };
218 
219 class
220 OCTAVE_API
221 FloatIndefQuad : public Quad
222 {
223 public:
224 
225   enum IntegralType { bound_to_inf, neg_inf_to_bound, doubly_infinite };
226 
FloatIndefQuad(float_integrand_fcn fcn)227   FloatIndefQuad (float_integrand_fcn fcn)
228     : Quad (fcn), bound (0.0), type (bound_to_inf) { }
229 
FloatIndefQuad(float_integrand_fcn fcn,double b,IntegralType t)230   FloatIndefQuad (float_integrand_fcn fcn, double b, IntegralType t)
231     : Quad (fcn), bound (b), type (t) { }
232 
233   ~FloatIndefQuad (void) = default;
234 
235   OCTAVE_NORETURN double do_integrate (octave_idx_type& ier,
236                                        octave_idx_type& neval, double& abserr);
237 
238   float do_integrate (octave_idx_type& ier, octave_idx_type& neval,
239                       float& abserr);
240 
241 private:
242 
243   float bound;
244   IntegralType type;
245 };
246 
247 #endif
248