DNDSR 0.2.1
Distributed Numeric Data Structure for CFV
Loading...
Searching...
No Matches
Quadrature.hpp
Go to the documentation of this file.
1#pragma once
2
3// =============================================================================
4// Quadrature.hpp - High-level quadrature interface
5// =============================================================================
6// This file provides the Quadrature class and related utilities for numerical
7// integration over finite elements. The actual quadrature data and low-level
8// access functions are in the Quadratures/ subdirectory.
9//
10// Structure:
11// - Quadratures/QuadratureBase.hpp: Constants and type aliases
12// - Quadratures/*.hpp: Individual parametric space quadrature data
13// - QuadratureHub.hpp: Unified interface for scheme selection and point access
14// - Quadrature.hpp: This file - Quadrature class and integration utilities
15// =============================================================================
16
17#include "QuadratureHub.hpp"
18
19namespace DNDS::Geom::Elem
20{
21 // =========================================================================
22 // SummationNoOp - Helper class for no-op accumulation
23 // =========================================================================
25 {
26 public:
28 {
29 }
30
32 {
33 return {};
34 }
35
37 {
38 return {};
39 }
40
42 {
43 return {};
44 }
45 };
46
47 // =========================================================================
48 namespace detail
49 {
50
51 // Precomputed Shape Function Buffers at Quadrature Points
52 // =========================================================================
53 static struct TNBufferAtQuadrature
54 {
55 std::array<std::array<std::vector<tD01Nj>, INT_ORDER_MAX + 1>, ElemType_NUM> buf;
56
57 TNBufferAtQuadrature()
58 {
59 for (t_index i = 1; i < ElemType_NUM; i++)
60 {
61 Element c_elem{ElemType(i)};
62 for (int order = 0; order <= INT_ORDER_MAX; order++)
63 {
64 auto int_scheme = GetQuadratureScheme(c_elem.GetParamSpace(), order);
65 buf.at(i).at(order).resize(int_scheme);
66 for (auto &m : buf.at(i).at(order))
67 m.resize(4, c_elem.GetNumNodes());
68 for (int iG = 0; iG < int_scheme; iG++)
69 {
70 tPoint pParam{0, 0, 0};
71 t_real w;
72 GetQuadraturePoint(c_elem.GetParamSpace(), int_scheme, iG, pParam, w);
73 c_elem.GetD01Nj(pParam, buf.at(i).at(order).at(iG));
74 }
75 }
76 }
77 }
78
79 } NBufferAtQuadrature{};
80
81 } // namespace detail
82
83 // =========================================================================
84 // Quadrature Class
85 // =========================================================================
86 // Main interface for numerical integration over elements.
87 //
88 // Usage:
89 // Element elem{Quad4};
90 // Quadrature quad(elem, 3); // Order 3 integration
91 //
92 // real sum = 0;
93 // quad.Integration(sum, [](real& acc, int iG, tPoint p, tD01Nj D01Nj) {
94 // acc = function_value(p);
95 // });
96 // =========================================================================
98 {
103
105 : elem(n_elem), int_order(n_int_order), ps(elem.GetParamSpace()),
107 {
108 }
109
110 /**
111 * @brief General integration method with shape function access
112 * @param f Callback: f(TAcc& acc, int iG, tPoint pParam, tD01Nj D01Nj)
113 */
114 template <class TAcc, class TFunc>
115 DNDS_DEVICE_CALLABLE void Integration(TAcc &buf, TFunc &&f)
116 {
117 for (t_index iG = 0; iG < int_scheme; iG++)
118 {
119 tPoint pParam{0, 0, 0};
120 t_real w;
121 GetQuadraturePoint(ps, int_scheme, iG, pParam, w);
122 TAcc acc;
123 f(acc, iG, pParam, detail::NBufferAtQuadrature.buf.at(elem.type).at(int_order).at(iG));
124 buf += acc * w;
125 }
126 }
127
128 /**
129 * @brief Simple integration without shape function access
130 * @param f Callback: f(TAcc& acc, int iG)
131 */
132 template <class TAcc, class TFunc>
133 DNDS_DEVICE_CALLABLE std::enable_if_t<std::is_invocable_v<TFunc, TAcc &, int>> IntegrationSimple(TAcc &buf, TFunc &&f)
134 {
135 for (t_index iG = 0; iG < int_scheme; iG++)
136 {
137 tPoint pParam{0, 0, 0};
138 t_real w;
139 GetQuadraturePoint(ps, int_scheme, iG, pParam, w);
140 TAcc acc;
141 f(acc, iG);
142 buf += acc * w;
143 }
144 }
145
146 /**
147 * @brief Simple integration with weight access
148 * @param f Callback: f(TAcc& acc, int iG, real w)
149 */
150 template <class TAcc, class TFunc>
151 DNDS_DEVICE_CALLABLE std::enable_if_t<std::is_invocable_v<TFunc, TAcc &, int, real>> IntegrationSimple(TAcc &buf, TFunc &&f)
152 {
153 for (t_index iG = 0; iG < int_scheme; iG++)
154 {
155 tPoint pParam{0, 0, 0};
156 t_real w;
157 GetQuadraturePoint(ps, int_scheme, iG, pParam, w);
158 TAcc acc;
159 f(acc, iG, w);
160 buf += acc * w;
161 }
162 }
163
164 DNDS_DEVICE_CALLABLE [[nodiscard]] auto GetQuadraturePointInfo(int iG) const
165 {
166 tPoint pParam{0, 0, 0};
167 t_real w;
168 GetQuadraturePoint(ps, int_scheme, iG, pParam, w);
169 return std::make_tuple(pParam, w);
170 }
171
172 DNDS_DEVICE_CALLABLE [[nodiscard]] real GetWeight(int iG) const
173 {
174 return std::get<1>(GetQuadraturePointInfo(iG));
175 }
176
177 DNDS_DEVICE_CALLABLE [[nodiscard]] t_index GetNumPoints() const { return int_scheme; }
178 };
179
180 // =========================================================================
181 // Quadrature Patch Generation (for visualization/post-processing)
182 // =========================================================================
183 // Generates triangle patches for visualizing quadrature point distributions.
184 // Used for debugging and visualization purposes.
185 // =========================================================================
187 {
188 std::vector<std::array<t_index, 3>> ret;
189 using tArr = std::array<t_index, 3>;
190 ret.reserve(40);
192 {
193 if (q.elem.GetOrder() == 1 || q.elem.GetOrder() == 2) // omitting the middle point
194 {
195 ret.emplace_back(tArr{1, -1, 0});
196 for (int i = 1; i < q.GetNumPoints(); i++)
197 ret.emplace_back(tArr{-i, -1 - i, 0});
198 ret.emplace_back(tArr{-q.GetNumPoints(), 2, 0});
199 }
200 }
202 {
203 if (q.elem.GetOrder() == 1 || q.elem.GetOrder() == 2) // omitting the middle point
204 {
205 auto nGL = static_cast<t_index>(std::round(std::sqrt(q.GetNumPoints())));
206 auto ijG = [=](int i, int j)
207 { return nGL * j + i; };
208 t_index iMax = nGL - 1;
209
210 if (q.elem.GetOrder() == 1)
211 {
212 ret.emplace_back(tArr{1, 2, -1 - ijG(0, 0)});
213 for (int i = 0; i < iMax; i++)
214 ret.emplace_back(tArr{-1 - ijG(i + 1, 0), -1 - ijG(i, 0), 2});
215 ret.emplace_back(tArr{3, 4, -1 - ijG(iMax, iMax)});
216 for (int i = 0; i < iMax; i++)
217 ret.emplace_back(tArr{-1 - ijG(i, iMax), -1 - ijG(i + 1, iMax), 4});
218 ret.emplace_back(tArr{2, 3, -1 - ijG(iMax, 0)});
219 for (int i = 0; i < iMax; i++)
220 ret.emplace_back(tArr{-1 - ijG(iMax, i + 1), -1 - ijG(iMax, i), 3});
221 ret.emplace_back(tArr{4, 1, -1 - ijG(0, iMax)});
222 for (int i = 0; i < iMax; i++)
223 ret.emplace_back(tArr{-1 - ijG(0, i), -1 - ijG(0, i + 1), 1});
224 }
225 else if (q.elem.GetOrder() == 2)
226 {
227 ret.emplace_back(tArr{1, 5, -1 - ijG(0, 0)});
228 ret.emplace_back(tArr{5, 2, -1 - ijG(iMax, 0)});
229 for (int i = 0; i < iMax; i++)
230 ret.emplace_back(tArr{-1 - ijG(i + 1, 0), -1 - ijG(i, 0), 5});
231 ret.emplace_back(tArr{3, 7, -1 - ijG(iMax, iMax)});
232 ret.emplace_back(tArr{7, 4, -1 - ijG(0, iMax)});
233 for (int i = 0; i < iMax; i++)
234 ret.emplace_back(tArr{-1 - ijG(i, iMax), -1 - ijG(i + 1, iMax), 7});
235 ret.emplace_back(tArr{2, 6, -1 - ijG(iMax, 0)});
236 ret.emplace_back(tArr{6, 3, -1 - ijG(iMax, iMax)});
237 for (int i = 0; i < iMax; i++)
238 ret.emplace_back(tArr{-1 - ijG(iMax, i + 1), -1 - ijG(iMax, i), 6});
239 ret.emplace_back(tArr{4, 8, -1 - ijG(0, iMax)});
240 ret.emplace_back(tArr{8, 1, -1 - ijG(0, 0)});
241 for (int i = 0; i < iMax; i++)
242 ret.emplace_back(tArr{-1 - ijG(0, i), -1 - ijG(0, i + 1), 8});
243 }
244
245 for (int i = 0; i < iMax; i++)
246 for (int j = 0; j < iMax; j++)
247 {
248 ret.emplace_back(tArr{-1 - ijG(i, j), -1 - ijG(i + 1, j), -1 - ijG(i, j + 1)});
249 ret.emplace_back(tArr{-1 - ijG(i + 1, j + 1), -1 - ijG(i, j + 1), -1 - ijG(i + 1, j)});
250 }
251 }
252 }
253 else if (q.elem.GetParamSpace() == ParamSpace::TriSpace)
254 {
255 if (q.GetNumPoints() == 1)
256 {
257 if (q.elem.GetOrder() == 1)
258 {
259 ret.emplace_back(tArr{1, 2, -1});
260 ret.emplace_back(tArr{2, 3, -1});
261 ret.emplace_back(tArr{3, 1, -1});
262 }
263 else if (q.elem.GetOrder() == 2)
264 {
265 ret.emplace_back(tArr{1, 4, -1});
266 ret.emplace_back(tArr{4, 2, -1});
267 ret.emplace_back(tArr{2, 5, -1});
268 ret.emplace_back(tArr{5, 3, -1});
269 ret.emplace_back(tArr{3, 6, -1});
270 ret.emplace_back(tArr{6, 1, -1});
271 }
272 }
273 else if (q.GetNumPoints() == 3)
274 {
275 if (q.elem.GetOrder() == 1)
276 {
277 ret.emplace_back(tArr{1, 2, -3});
278 ret.emplace_back(tArr{2, -1, -3});
279 ret.emplace_back(tArr{2, 3, -1});
280 ret.emplace_back(tArr{3, -2, -1});
281 ret.emplace_back(tArr{3, 1, -2});
282 ret.emplace_back(tArr{1, -3, -2});
283
284 ret.emplace_back(tArr{-1, -2, -3});
285 }
286 else if (q.elem.GetOrder() == 2)
287 {
288 ret.emplace_back(tArr{1, 4, -3});
289 ret.emplace_back(tArr{4, 2, -1});
290 ret.emplace_back(tArr{2, 5, -1});
291 ret.emplace_back(tArr{5, 3, -2});
292 ret.emplace_back(tArr{3, 6, -2});
293 ret.emplace_back(tArr{6, 1, -3});
294
295 ret.emplace_back(tArr{-1, -3, 4});
296 ret.emplace_back(tArr{-2, -1, 5});
297 ret.emplace_back(tArr{-3, -2, 6});
298
299 ret.emplace_back(tArr{-1, -2, -3});
300 }
301 }
302 else if (q.GetNumPoints() == 6)
303 {
304 if (q.elem.GetOrder() == 1)
305 {
306 ret.emplace_back(tArr{1, 2, -6});
307 ret.emplace_back(tArr{2, 3, -4});
308 ret.emplace_back(tArr{3, 1, -5});
309
310 ret.emplace_back(tArr{-6, -1, 1});
311 ret.emplace_back(tArr{-2, -6, 2});
312 ret.emplace_back(tArr{-4, -2, 2});
313 ret.emplace_back(tArr{-3, -4, 3});
314 ret.emplace_back(tArr{-5, -3, 3});
315 ret.emplace_back(tArr{-1, -5, 1});
316
317 ret.emplace_back(tArr{-4, -5, -6});
318 ret.emplace_back(tArr{-1, -6, -5});
319 ret.emplace_back(tArr{-6, -2, -4});
320 ret.emplace_back(tArr{-5, -4, -3});
321 }
322 else if (q.elem.GetOrder() == 2)
323 {
324 ret.emplace_back(tArr{1, 4, -6});
325 ret.emplace_back(tArr{4, 2, -6});
326 ret.emplace_back(tArr{2, 5, -4});
327 ret.emplace_back(tArr{5, 3, -4});
328 ret.emplace_back(tArr{3, 6, -5});
329 ret.emplace_back(tArr{6, 1, -5});
330
331 ret.emplace_back(tArr{-6, -1, 1});
332 ret.emplace_back(tArr{-2, -6, 2});
333 ret.emplace_back(tArr{-4, -2, 2});
334 ret.emplace_back(tArr{-3, -4, 3});
335 ret.emplace_back(tArr{-5, -3, 3});
336 ret.emplace_back(tArr{-1, -5, 1});
337
338 ret.emplace_back(tArr{-4, -5, -6});
339 ret.emplace_back(tArr{-1, -6, -5});
340 ret.emplace_back(tArr{-6, -2, -4});
341 ret.emplace_back(tArr{-5, -4, -3});
342 }
343 }
344 else if (q.GetNumPoints() == 7)
345 {
346 static const std::array<int, 7> idxTri6ToTri7{
347 -1, 5, 6, 7, 2, 3, 4};
348 if (q.elem.GetOrder() == 1)
349 {
350 ret.emplace_back(tArr{1, 2, -idxTri6ToTri7[6]});
351 ret.emplace_back(tArr{2, 3, -idxTri6ToTri7[4]});
352 ret.emplace_back(tArr{3, 1, -idxTri6ToTri7[5]});
353
354 ret.emplace_back(tArr{-idxTri6ToTri7[6], -idxTri6ToTri7[1], 1});
355 ret.emplace_back(tArr{-idxTri6ToTri7[2], -idxTri6ToTri7[6], 2});
356 ret.emplace_back(tArr{-idxTri6ToTri7[4], -idxTri6ToTri7[2], 2});
357 ret.emplace_back(tArr{-idxTri6ToTri7[3], -idxTri6ToTri7[4], 3});
358 ret.emplace_back(tArr{-idxTri6ToTri7[5], -idxTri6ToTri7[3], 3});
359 ret.emplace_back(tArr{-idxTri6ToTri7[1], -idxTri6ToTri7[5], 1});
360
361 ret.emplace_back(tArr{-idxTri6ToTri7[4], -idxTri6ToTri7[5], -1});
362 ret.emplace_back(tArr{-idxTri6ToTri7[5], -idxTri6ToTri7[6], -1});
363 ret.emplace_back(tArr{-idxTri6ToTri7[6], -idxTri6ToTri7[4], -1});
364
365 ret.emplace_back(tArr{-idxTri6ToTri7[1], -idxTri6ToTri7[6], -idxTri6ToTri7[5]});
366 ret.emplace_back(tArr{-idxTri6ToTri7[6], -idxTri6ToTri7[2], -idxTri6ToTri7[4]});
367 ret.emplace_back(tArr{-idxTri6ToTri7[5], -idxTri6ToTri7[4], -idxTri6ToTri7[3]});
368 }
369 else if (q.elem.GetOrder() == 2)
370 {
371 ret.emplace_back(tArr{1, 4, -idxTri6ToTri7[6]});
372 ret.emplace_back(tArr{4, 2, -idxTri6ToTri7[6]});
373 ret.emplace_back(tArr{2, 5, -idxTri6ToTri7[4]});
374 ret.emplace_back(tArr{5, 3, -idxTri6ToTri7[4]});
375 ret.emplace_back(tArr{3, 6, -idxTri6ToTri7[5]});
376 ret.emplace_back(tArr{6, 1, -idxTri6ToTri7[5]});
377
378 ret.emplace_back(tArr{-idxTri6ToTri7[6], -idxTri6ToTri7[1], 1});
379 ret.emplace_back(tArr{-idxTri6ToTri7[2], -idxTri6ToTri7[6], 2});
380 ret.emplace_back(tArr{-idxTri6ToTri7[4], -idxTri6ToTri7[2], 2});
381 ret.emplace_back(tArr{-idxTri6ToTri7[3], -idxTri6ToTri7[4], 3});
382 ret.emplace_back(tArr{-idxTri6ToTri7[5], -idxTri6ToTri7[3], 3});
383 ret.emplace_back(tArr{-idxTri6ToTri7[1], -idxTri6ToTri7[5], 1});
384
385 ret.emplace_back(tArr{-idxTri6ToTri7[4], -idxTri6ToTri7[5], -1});
386 ret.emplace_back(tArr{-idxTri6ToTri7[5], -idxTri6ToTri7[6], -1});
387 ret.emplace_back(tArr{-idxTri6ToTri7[6], -idxTri6ToTri7[4], -1});
388
389 ret.emplace_back(tArr{-idxTri6ToTri7[1], -idxTri6ToTri7[6], -idxTri6ToTri7[5]});
390 ret.emplace_back(tArr{-idxTri6ToTri7[6], -idxTri6ToTri7[2], -idxTri6ToTri7[4]});
391 ret.emplace_back(tArr{-idxTri6ToTri7[5], -idxTri6ToTri7[4], -idxTri6ToTri7[3]});
392 }
393 }
394 else if (q.GetNumPoints() == 12)
395 {
396 if (q.elem.GetOrder() == 1)
397 {
398 ret.emplace_back(tArr{1, 2, -4});
399 ret.emplace_back(tArr{2, 3, -5});
400 ret.emplace_back(tArr{3, 1, -6});
401
402 ret.emplace_back(tArr{-7, -4, 2});
403 ret.emplace_back(tArr{-9, -7, 2});
404 ret.emplace_back(tArr{-5, -9, 2});
405
406 ret.emplace_back(tArr{-11, -5, 3});
407 ret.emplace_back(tArr{-12, -11, 3});
408 ret.emplace_back(tArr{-6, -12, 3});
409
410 ret.emplace_back(tArr{-10, -6, 1});
411 ret.emplace_back(tArr{-8, -10, 1});
412 ret.emplace_back(tArr{-4, -8, 1});
413
414 ret.emplace_back(tArr{-4, -7, -8});
415 ret.emplace_back(tArr{-9, -5, -11});
416 ret.emplace_back(tArr{-10, -12, -6});
417 ret.emplace_back(tArr{-8, -7, -1});
418 ret.emplace_back(tArr{-2, -9, -11});
419 ret.emplace_back(tArr{-3, -12, -10});
420
421 ret.emplace_back(tArr{-7, -9, -2});
422 ret.emplace_back(tArr{-7, -2, -1});
423 ret.emplace_back(tArr{-11, -12, -3});
424 ret.emplace_back(tArr{-11, -3, -2});
425 ret.emplace_back(tArr{-10, -8, -1});
426 ret.emplace_back(tArr{-10, -1, -3});
427
428 ret.emplace_back(tArr{-1, -2, -3});
429 }
430 else if (q.elem.GetOrder() == 2)
431 {
432 ret.emplace_back(tArr{1, 4, -4});
433 ret.emplace_back(tArr{4, 2, -5});
434 ret.emplace_back(tArr{2, 5, -5});
435 ret.emplace_back(tArr{5, 3, -6});
436 ret.emplace_back(tArr{3, 6, -6});
437 ret.emplace_back(tArr{6, 1, -4});
438
439 ret.emplace_back(tArr{-7, -4, 4});
440 ret.emplace_back(tArr{-9, -7, 4});
441 ret.emplace_back(tArr{-5, -9, 4});
442
443 ret.emplace_back(tArr{-11, -5, 5});
444 ret.emplace_back(tArr{-12, -11, 5});
445 ret.emplace_back(tArr{-6, -12, 5});
446
447 ret.emplace_back(tArr{-10, -6, 6});
448 ret.emplace_back(tArr{-8, -10, 6});
449 ret.emplace_back(tArr{-4, -8, 6});
450
451 ret.emplace_back(tArr{-4, -7, -8});
452 ret.emplace_back(tArr{-9, -5, -11});
453 ret.emplace_back(tArr{-10, -12, -6});
454 ret.emplace_back(tArr{-8, -7, -1});
455 ret.emplace_back(tArr{-2, -9, -11});
456 ret.emplace_back(tArr{-3, -12, -10});
457
458 ret.emplace_back(tArr{-7, -9, -2});
459 ret.emplace_back(tArr{-7, -2, -1});
460 ret.emplace_back(tArr{-11, -12, -3});
461 ret.emplace_back(tArr{-11, -3, -2});
462 ret.emplace_back(tArr{-10, -8, -1});
463 ret.emplace_back(tArr{-10, -1, -3});
464
465 ret.emplace_back(tArr{-1, -2, -3});
466 }
467 }
468 else
469 DNDS_assert(false);
470 }
471 else
472 DNDS_assert(false);
473 return ret;
474 }
475
476 // =========================================================================
477 // Jacobian Determinant Helpers
478 // =========================================================================
479 template <int dim>
481 {
482 using namespace Geom;
483 real JDet{0};
484 tJacobi J = Elem::ShapeJacobianCoordD01Nj(coordsCell, DiNj);
485 if constexpr (dim == 2)
486 JDet = J(EigenAll, 0).cross(J(EigenAll, 1)).stableNorm();
487 else
488 JDet = J.fullPivLu().determinant();
489 return JDet;
490 }
491
492 inline real CellJacobianDet(int dim, const Geom::tSmallCoords &coordsCell, const Geom::Elem::tD01Nj &DiNj)
493 {
494 if (dim == 2)
495 return CellJacobianDet<2>(coordsCell, DiNj);
496 else
497 return CellJacobianDet<3>(coordsCell, DiNj);
498 }
499
500 template <int dim>
502 {
503 using namespace Geom;
504 real JDet{0};
505 tJacobi J = Elem::ShapeJacobianCoordD01Nj(coords, DiNj);
506 if constexpr (dim == 2)
507 JDet = J(EigenAll, 0).stableNorm();
508 else
509 JDet = J(EigenAll, 0).cross(J(EigenAll, 1)).stableNorm();
510 return JDet;
511 }
512
513 inline real FaceJacobianDet(int dim, const Geom::tSmallCoords &coords, const Geom::Elem::tD01Nj &DiNj)
514 {
515 if (dim == 2)
516 return FaceJacobianDet<2>(coords, DiNj);
517 else
518 return FaceJacobianDet<3>(coords, DiNj);
519 }
520
521} // namespace DNDS::Geom::Elem
#define DNDS_DEVICE_CALLABLE
Definition Defines.hpp:76
#define DNDS_assert(expr)
Debug-only assertion (compiled out when DNDS_NDEBUG is defined). Prints the expression + file/line + ...
Definition Errors.hpp:112
Eigen::Matrix< real, 3, 3 > m
SummationNoOp operator*(const SummationNoOp &R) const
friend SummationNoOp operator*(const SummationNoOp &L, t_real R)
void operator+=(const SummationNoOp &R)
friend SummationNoOp operator*(t_real L, const SummationNoOp &R)
real FaceJacobianDet(const Geom::tSmallCoords &coords, const Geom::Elem::tD01Nj &DiNj)
constexpr t_index GetQuadratureScheme(ParamSpace ps, int int_order)
void GetQuadraturePoint(ParamSpace ps, t_index scheme, int iG, TPoint &pParam, t_real &w)
Eigen::Matrix< t_real, 4, Eigen::Dynamic > tD01Nj
Definition Elements.hpp:203
real CellJacobianDet(const Geom::tSmallCoords &coordsCell, const Geom::Elem::tD01Nj &DiNj)
auto GetQuadPatches(Quadrature &q)
int32_t t_index
Definition Geometric.hpp:6
Eigen::Matrix3d tJacobi
Definition Geometric.hpp:10
Eigen::Vector3d tPoint
Definition Geometric.hpp:9
Eigen::Matrix< real, 3, Eigen::Dynamic > tSmallCoords
Definition Geometric.hpp:50
double t_real
Definition Geometric.hpp:8
double real
Canonical floating-point scalar used throughout DNDSR (double precision).
Definition Defines.hpp:110
DNDS_DEVICE_CALLABLE constexpr ParamSpace GetParamSpace() const
Definition Elements.hpp:212
DNDS_DEVICE_CALLABLE constexpr t_index GetOrder() const
Definition Elements.hpp:223
DNDS_DEVICE_CALLABLE std::enable_if_t< std::is_invocable_v< TFunc, TAcc &, int > > IntegrationSimple(TAcc &buf, TFunc &&f)
Simple integration without shape function access.
DNDS_DEVICE_CALLABLE t_index GetNumPoints() const
DNDS_DEVICE_CALLABLE Quadrature(Element n_elem=Element{UnknownElem}, int n_int_order=0)
DNDS_DEVICE_CALLABLE auto GetQuadraturePointInfo(int iG) const
DNDS_DEVICE_CALLABLE std::enable_if_t< std::is_invocable_v< TFunc, TAcc &, int, real > > IntegrationSimple(TAcc &buf, TFunc &&f)
Simple integration with weight access.
DNDS_DEVICE_CALLABLE real GetWeight(int iG) const
DNDS_DEVICE_CALLABLE void Integration(TAcc &buf, TFunc &&f)
General integration method with shape function access.
double order
Definition test_ODE.cpp:257