/home/runner/work/kynema/kynema/kynema/src/math/gll_quadrature.hpp Source File

Kynema API: /home/runner/work/kynema/kynema/kynema/src/math/gll_quadrature.hpp Source File
Kynema API
A flexible multibody structural dynamics code for wind turbines
Loading...
Searching...
No Matches
gll_quadrature.hpp
Go to the documentation of this file.
1#pragma once
2
3#include <cmath>
4#include <stdexcept>
5#include <vector>
6
7namespace kynema::math {
8
9inline std::vector<double> GetGllLocations(size_t order) {
10 switch (order) {
11 case 1UL:
12 return std::vector{-1., 1.};
13 case 2UL:
14 return std::vector{-1., 0., 1.};
15 case 3UL:
16 return std::vector{-1., -std::sqrt(1. / 5.), std::sqrt(1. / 5.), 1.};
17 case 4UL:
18 return std::vector{-1., -std::sqrt(3. / 7.), 0., std::sqrt(3. / 7.), 1.};
19 case 5UL:
20 return std::vector{
21 -1.,
22 -std::sqrt((1. / 3.) + (2. * std::sqrt(7.) / 21.)),
23 -std::sqrt((1. / 3.) - (2. * std::sqrt(7.) / 21.)),
24 std::sqrt((1. / 3.) - (2. * std::sqrt(7.) / 21.)),
25 std::sqrt((1. / 3.) + (2. * std::sqrt(7.) / 21.)),
26 1.
27 };
28 case 6UL:
29 return std::vector{
30 -1.,
31 -std::sqrt((5. / 11.) + ((2. / 11.) * std::sqrt(5. / 3.))),
32 -std::sqrt((5. / 11.) - ((2. / 11.) * std::sqrt(5. / 3.))),
33 0.,
34 std::sqrt((5. / 11.) - ((2. / 11.) * std::sqrt(5. / 3.))),
35 std::sqrt((5. / 11.) + ((2. / 11.) * std::sqrt(5. / 3.))),
36 1.
37 };
38 case 7UL:
39 return std::vector{
40 -1.,
41 -0.87174014850960657,
42 -0.59170018143314229,
43 -0.20929921790247885,
44 0.20929921790247885,
45 0.59170018143314229,
46 0.87174014850960657,
47 1.
48 };
49 case 8UL:
50 return std::vector{-1., -0.89975799541146018, -0.67718627951073773, -0.36311746382617816,
51 0., 0.36311746382617816, 0.67718627951073773, 0.89975799541146018,
52 1.};
53 case 9UL:
54 return std::vector{
55 -1.,
56 -0.91953390816645886,
57 -0.73877386510550502,
58 -0.47792494981044453,
59 -0.16527895766638701,
60 0.16527895766638701,
61 0.47792494981044448,
62 0.73877386510550502,
63 0.91953390816645886,
64 1.
65 };
66 case 10UL:
67 return std::vector{-1.,
68 -0.93400143040805916,
69 -0.78448347366314441,
70 -0.56523532699620493,
71 -0.29575813558693936,
72 0.,
73 0.29575813558693936,
74 0.56523532699620493,
75 0.78448347366314441,
76 0.93400143040805916,
77 1.};
78 case 11UL:
79 return std::vector{
80 -1.,
81 -0.94489927222288217,
82 -0.81927932164400663,
83 -0.63287615303186073,
84 -0.39953094096534891,
85 -0.13655293285492756,
86 0.13655293285492756,
87 0.39953094096534891,
88 0.63287615303186073,
89 0.81927932164400663,
90 0.94489927222288217,
91 1.
92 };
93 case 12UL:
94 return std::vector{
95 -1.,
96 -0.95330984664216400,
97 -0.84634756465187233,
98 -0.68618846908175746,
99 -0.48290982109133623,
100 -0.24928693010623998,
101 0.,
102 0.24928693010623998,
103 0.48290982109133623,
104 0.68618846908175746,
105 0.84634756465187233,
106 0.95330984664216400,
107 1.
108 };
109 case 13UL:
110 return std::vector{
111 -1.,
112 -0.95993504526726092,
113 -0.86780105383034722,
114 -0.72886859909132617,
115 -0.55063940292864699,
116 -0.34272401334271285,
117 -0.11633186888370387,
118 0.11633186888370387,
119 0.34272401334271285,
120 0.55063940292864710,
121 0.72886859909132617,
122 0.86780105383034722,
123 0.95993504526726092,
124 1.
125 };
126 case 14UL:
127 return std::vector{
128 -1.,
129 -0.96524592650383856,
130 -0.88508204422297632,
131 -0.76351968995181518,
132 -0.60625320546984574,
133 -0.42063805471367249,
134 -0.21535395536379423,
135 0.,
136 0.21535395536379423,
137 0.42063805471367249,
138 0.60625320546984574,
139 0.76351968995181518,
140 0.88508204422297632,
141 0.96524592650383856,
142 1.
143 };
144 case 15UL:
145 return std::vector{
146 -1.,
147 -0.96956804627021798,
148 -0.89920053309347214,
149 -0.79200829186181509,
150 -0.65238870288249307,
151 -0.48605942188713763,
152 -0.29983046890076320,
153 -0.10132627352194944,
154 0.10132627352194944,
155 0.29983046890076320,
156 0.48605942188713763,
157 0.65238870288249307,
158 0.79200829186181509,
159 0.89920053309347214,
160 0.96956804627021798,
161 1.
162 };
163 default:
164 throw std::runtime_error("Supported orders are 1 - 15");
165 }
166}
167
168inline std::vector<double> GetGllWeights(size_t order) {
169 switch (order) {
170 case 1UL:
171 return std::vector{1., 1.};
172 case 2UL:
173 return std::vector{1. / 3., 4. / 3., 1. / 3.};
174 case 3UL:
175 return std::vector{1. / 6., 5. / 6., 5. / 6., 1. / 6.};
176 case 4UL:
177 return std::vector{1. / 10., 49. / 90., 32. / 45., 49. / 90., 1. / 10.};
178 case 5UL:
179 return std::vector{
180 1. / 15.,
181 (14. - std::sqrt(7.)) / 30.,
182 (14. + std::sqrt(7.)) / 30.,
183 (14. + std::sqrt(7.)) / 30.,
184 (14. - std::sqrt(7.)) / 30.,
185 1. / 15.
186 };
187 case 6UL:
188 return std::vector{1. / 21.,
189 (124. - 7. * std::sqrt(15.)) / 350.,
190 (124. + 7. * std::sqrt(15.)) / 350.,
191 256. / 525.,
192 (124. + 7. * std::sqrt(15.)) / 350.,
193 (124. - 7. * std::sqrt(15.)) / 350,
194 1. / 21.};
195 case 7UL:
196 return std::vector{3.5714285714285712E-002, 0.21070422714350615, 0.34112269248350441,
197 0.41245879465870372, 0.41245879465870372, 0.34112269248350441,
198 0.21070422714350615, 3.5714285714285712E-002};
199 case 8UL:
200 return std::vector{2.7777777777777776E-002, 0.16549536156080552,
201 0.27453871250016160, 0.34642851097304617,
202 0.37151927437641724, 0.34642851097304617,
203 0.27453871250016160, 0.16549536156080552,
204 2.7777777777777776E-002};
205 case 9UL:
206 return std::vector{2.2222222222222223E-002, 0.13330599085107006, 0.22488934206312636,
207 0.29204268367968378, 0.32753976118389744, 0.32753976118389744,
208 0.29204268367968378, 0.22488934206312636, 0.13330599085107006,
209 2.2222222222222223E-002};
210 case 10UL:
211 return std::vector{1.8181818181818181E-002, 0.10961227326699498, 0.18716988178030541,
212 0.24804810426402829, 0.28687912477900801, 0.30021759545569071,
213 0.28687912477900823, 0.24804810426402829, 0.18716988178030541,
214 0.10961227326699498, 1.8181818181818181E-002};
215 case 11UL:
216 return std::vector{1.5151515151515152E-002, 9.1684517413196109E-002,
217 0.15797470556437004, 0.21250841776102122,
218 0.25127560319920106, 0.27140524091069618,
219 0.27140524091069618, 0.25127560319920106,
220 0.21250841776102122, 0.15797470556437004,
221 9.1684517413196109E-002, 1.5151515151515152E-002};
222 case 12UL:
223 return std::vector{1.2820512820512820E-002, 7.7801686746818866E-002,
224 0.13498192668960840, 0.18364686520355006,
225 0.22076779356611012, 0.24401579030667625,
226 0.25193084933344673, 0.24401579030667625,
227 0.22076779356611012, 0.18364686520355006,
228 0.13498192668960840, 7.7801686746818866E-002,
229 1.2820512820512820E-002};
230 case 13UL:
231 return std::vector{1.0989010989010990E-002, 6.6837284497681185E-002, 0.11658665589871173,
232 0.16002185176295217, 0.19482614937341614, 0.21912625300977098,
233 0.23161279446845698, 0.23161279446845698, 0.21912625300977098,
234 0.19482614937341600, 0.16002185176295217, 0.11658665589871173,
235 6.6837284497681185E-002, 1.0989010989010990E-002};
236 case 14UL:
237 return std::vector{
238 9.5238095238095247E-003, 5.8029893028601086E-002, 0.10166007032571801,
239 0.14051169980242798, 0.17278964725360088, 0.19698723596461334,
240 0.21197358592682095, 0.21704811634881566, 0.21197358592682075,
241 0.19698723596461334, 0.17278964725360088, 0.14051169980242798,
242 0.10166007032571801, 5.8029893028601086E-002, 9.5238095238095247E-003
243 };
244 case 15UL:
245 return std::vector{8.3333333333333332E-003, 5.0850361005920039E-002,
246 8.9393697325930832E-002, 0.12425538213251400,
247 0.15402698080716443, 0.17749191339170411,
248 0.19369002382520362, 0.20195830817822993,
249 0.20195830817822993, 0.19369002382520362,
250 0.17749191339170411, 0.15402698080716443,
251 0.12425538213251400, 8.9393697325930832E-002,
252 5.0850361005920039E-002, 8.3333333333333332E-003};
253 default:
254 throw std::runtime_error("Supported orders are 1 - 15");
255 }
256}
257
258} // namespace kynema::math
Definition gll_quadrature.hpp:7
std::vector< double > GetGllLocations(size_t order)
Definition gll_quadrature.hpp:9
std::vector< double > GetGllWeights(size_t order)
Definition gll_quadrature.hpp:168