Hex: src/gauss.h Source File

| Hex 2.2 Hydrogen-electron collision solver

gauss.h
Go to the documentation of this file.
1 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
2 // //
3 // / / / / __ \ \ / / //
4 // / /__ / / / _ \ \ \/ / //
5 // / ___ / | |/_/ / /\ \ //
6 // / / / / \_\ / / \ \ //
7 // //
8 // //
9 // Copyright (c) 2016, Jakub Benda, Charles University in Prague //
10 // //
11 // MIT License: //
12 // //
13 // Permission is hereby granted, free of charge, to any person obtaining a //
14 // copy of this software and associated documentation files (the "Software"), //
15 // to deal in the Software without restriction, including without limitation //
16 // the rights to use, copy, modify, merge, publish, distribute, sublicense, //
17 // and/or sell copies of the Software, and to permit persons to whom the //
18 // Software is furnished to do so, subject to the following conditions: //
19 // //
20 // The above copyright notice and this permission notice shall be included //
21 // in all copies or substantial portions of the Software. //
22 // //
23 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS //
24 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, //
25 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE //
26 // AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, //
27 // WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF //
28 // OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. //
29 // //
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31
32 #ifndef HEX_GAUSS
33 #define HEX_GAUSS
34
35 #include <map>
36 #include <vector>
37
38 #include "bspline.h"
39
40 class GaussLegendreData
41 {
42 public:
43
52 static void gauss_nodes_and_weights (int points, const Real*& vx, const Real*& vw);
53
59 static void precompute_nodes_and_weights (int points);
60
61 private:
62
63 // precomputed nodes and weights, common to all instances
64 static std::vector<std::pair<Real*,Real*>> data_;
65 };
66
79 class GaussLegendre : public GaussLegendreData
80 {
81 public:
82
90 void scaled_nodes_and_weights (int points, Complex x1, Complex x2, Complex* xs, Complex* ws) const;
91
115 template <class Functor, class... Data> Complex quad
116 (
117 Functor f, Bspline const & bspline,
118 int points, int iknot, Complex x1, Complex x2,
119 Data... data
120 ) const
121 {
122 // check boundaries
123 if (x1.real() < bspline.t(iknot).real() or bspline.t(iknot+1).real() < x1.real() or
124 x2.real() < bspline.t(iknot).real() or bspline.t(iknot+1).real() < x2.real())
125 {
126 HexException("[quad] Error: boundaries not for this iknot!");
127 }
128
129 // get evaluation points and weights
130 std::vector<Complex> xs(points), ws(points), fs(points);
131 scaled_nodes_and_weights(points, x1, x2, xs.data(), ws.data());
132
133 // evaluate the function
134 f(points, xs.data(), fs.data(), data...);
135
136 // sum the results
137 Complex result = 0.;
138 for (int ipt = 0; ipt < points; ipt++)
139 result += ws[ipt] * fs[ipt];
140
141 return result;
142 }
143
155 template <class ClassPtr, class Functor, class... Data> Complex quadMFP
156 (
157 ClassPtr ptr, Functor f, Bspline const & bspline,
158 int points, int iknot, Complex x1, Complex x2,
159 Data... data
160 ) const
161 {
162 // check boundaries
163 if (x1.real() < bspline.t(iknot).real() or bspline.t(iknot+1).real() < x1.real() or
164 x2.real() < bspline.t(iknot).real() or bspline.t(iknot+1).real() < x2.real())
165 {
166 HexException("[quad] Error: boundaries not for this iknot!");
167 }
168
169 // get evaluation points and weights
170 std::vector<Complex> xs(points), ws(points), fs(points);
171 scaled_nodes_and_weights(points, x1, x2, xs.data(), ws.data());
172
173 // evaluate the function
174 (ptr->*f)(points, xs.data(), fs.data(), data...);
175
176 // sum the results
177 Complex result = 0.;
178 for (int ipt = 0; ipt < points; ipt++)
179 result += ws[ipt] * fs[ipt];
180
181 return result;
182 }
183 };
184
185 #endif
GaussLegendreData::precompute_nodes_and_weights
static void precompute_nodes_and_weights(int points)
Precalculate nodes and weights so that the retrieval is fast and thread-safe.
Definition: gauss.cpp:49
bspline.h
Bspline
B-spline environment.
Definition: bspline.h:57
GaussLegendre
Gauss-Legendre quadrature.
Definition: gauss.h:79
GaussLegendreData::gauss_nodes_and_weights
static void gauss_nodes_and_weights(int points, const Real *&vx, const Real *&vw)
Retrieve Gauss-Legendre data.
Definition: gauss.cpp:80
Bspline::t
cArrayView t() const
B-spline knot sequence.
Definition: bspline.h:266
GaussLegendreData
Definition: gauss.h:40