MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
quadinterpolator.cpp
Go to the documentation of this file.
1// Copyright (c) 2010-2025, Lawrence Livermore National Security, LLC. Produced
2// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3// LICENSE and NOTICE for details. LLNL-CODE-806117.
4//
5// This file is part of the MFEM library. For more information and source code
6// availability visit https://mfem.org.
7//
8// MFEM is free software; you can redistribute it and/or modify it under the
9// terms of the BSD-3 license. We welcome feedback and contributions, see file
10// CONTRIBUTING.md for details.
11
12#include "quadinterpolator.hpp"
13#include "qinterp/grad.hpp"
14#include "qinterp/eval.hpp"
15#include "qspace.hpp"
16#include "../general/forall.hpp"
17#include "../linalg/dtensor.hpp"
18#include "../linalg/kernels.hpp"
19
20namespace mfem
21{
22
23namespace internal
24{
25namespace quadrature_interpolator
26{
27void InitEvalByNodesKernels();
28void InitEvalByVDimKernels();
29void InitEvalKernels();
30void InitDetKernels();
31template <bool P> void InitGradByNodesKernels();
32template <bool P> void InitGradByVDimKernels();
33void InitTensorEvalHDivKernels();
34struct Kernels
35{
36 Kernels()
37 {
38 using namespace internal::quadrature_interpolator;
39
40 InitEvalByNodesKernels();
41 InitEvalByVDimKernels();
42 // Non-phys grad kernels
43 InitGradByNodesKernels<false>();
44 InitGradByVDimKernels<false>();
45 // Phys grad kernels
46 InitGradByNodesKernels<true>();
47 InitGradByVDimKernels<true>();
48 // Determinants
49 InitDetKernels();
50 // Non-tensor
51 InitEvalKernels();
52 // Tensor (quad,hex) H(div)
53 InitTensorEvalHDivKernels();
54 }
55};
56}
57}
58
59QuadratureInterpolator::QuadratureInterpolator(const FiniteElementSpace &fes,
60 const IntegrationRule &ir):
61
62 fespace(&fes),
63 qspace(nullptr),
64 IntRule(&ir),
65 q_layout(QVectorLayout::byNODES),
66 use_tensor_products(UsesTensorBasis(fes))
67{
68 static internal::quadrature_interpolator::Kernels kernels;
69
70 d_buffer.UseDevice(true);
71 if (fespace->GetNE() == 0) { return; }
72 const FiniteElement *fe = fespace->GetTypicalFE();
73 MFEM_VERIFY(fe->GetMapType() == FiniteElement::MapType::VALUE ||
74 fe->GetMapType() == FiniteElement::MapType::H_DIV,
75 "Only elements with MapType VALUE and H_DIV are supported!");
76}
77
78QuadratureInterpolator::QuadratureInterpolator(const FiniteElementSpace &fes,
79 const QuadratureSpace &qs):
80
81 fespace(&fes),
82 qspace(&qs),
83 IntRule(nullptr),
84 q_layout(QVectorLayout::byNODES),
85 use_tensor_products(UsesTensorBasis(fes))
86{
87 d_buffer.UseDevice(true);
88 if (fespace->GetNE() == 0) { return; }
89 const FiniteElement *fe = fespace->GetTypicalFE();
90 MFEM_VERIFY(fe->GetMapType() == FiniteElement::MapType::VALUE ||
91 fe->GetMapType() == FiniteElement::MapType::H_DIV,
92 "Only elements with MapType VALUE and H_DIV are supported!");
93}
94
95namespace internal
96{
97
98namespace quadrature_interpolator
99{
100
101// Compute kernel for 1D quadrature interpolation:
102// * non-tensor product version,
103// * assumes 'e_vec' is using ElementDofOrdering::NATIVE,
104// * assumes 'maps.mode == FULL'.
105static void Eval1D(const int NE,
106 const int vdim,
107 const QVectorLayout q_layout,
108 const GeometricFactors *geom,
109 const DofToQuad &maps,
110 const Vector &e_vec,
111 Vector &q_val,
112 Vector &q_der,
113 Vector &q_det,
114 const int eval_flags)
115{
116 using QI = QuadratureInterpolator;
117
118 const int nd = maps.ndof;
119 const int nq = maps.nqpt;
120 MFEM_ASSERT(maps.mode == DofToQuad::FULL, "internal error");
121 MFEM_ASSERT(!geom || geom->mesh->SpaceDimension() == 1, "");
122 MFEM_VERIFY(vdim == 1 || !(eval_flags & QI::DETERMINANTS), "");
123 MFEM_VERIFY(bool(geom) == bool(eval_flags & QI::PHYSICAL_DERIVATIVES),
124 "'geom' must be given (non-null) only when evaluating physical"
125 " derivatives");
126 const auto B = Reshape(maps.B.Read(), nq, nd);
127 const auto G = Reshape(maps.G.Read(), nq, nd);
128 const auto J = Reshape(geom ? geom->J.Read() : nullptr, nq, NE);
129 const auto E = Reshape(e_vec.Read(), nd, vdim, NE);
130 auto val = q_layout == QVectorLayout::byNODES ?
131 Reshape(q_val.Write(), nq, vdim, NE):
132 Reshape(q_val.Write(), vdim, nq, NE);
133 auto der = q_layout == QVectorLayout::byNODES ?
134 Reshape(q_der.Write(), nq, vdim, NE):
135 Reshape(q_der.Write(), vdim, nq, NE);
136 auto det = Reshape(q_det.Write(), nq, NE);
137 mfem::forall(NE, [=] MFEM_HOST_DEVICE (int e)
138 {
139 for (int q = 0; q < nq; ++q)
140 {
141 if (eval_flags & QI::VALUES)
142 {
143 for (int c = 0; c < vdim; c++)
144 {
145 real_t q_val = 0.0;
146 for (int d = 0; d < nd; ++d)
147 {
148 q_val += B(q,d)*E(d,c,e);
149 }
150 if (q_layout == QVectorLayout::byVDIM) { val(c,q,e) = q_val; }
151 if (q_layout == QVectorLayout::byNODES) { val(q,c,e) = q_val; }
152 }
153 }
154 if ((eval_flags & QI::DERIVATIVES) ||
155 (eval_flags & QI::PHYSICAL_DERIVATIVES) ||
156 (eval_flags & QI::DETERMINANTS))
157 {
158 for (int c = 0; c < vdim; c++)
159 {
160 real_t q_d = 0.0;
161 for (int d = 0; d < nd; ++d)
162 {
163 q_d += G(q,d)*E(d,c,e);
164 }
165 if (eval_flags & QI::PHYSICAL_DERIVATIVES)
166 {
167 q_d /= J(q,e);
168 }
169 if (eval_flags & QI::DERIVATIVES || eval_flags & QI::PHYSICAL_DERIVATIVES)
170 {
171 if (q_layout == QVectorLayout::byVDIM) { der(c,q,e) = q_d; }
172 if (q_layout == QVectorLayout::byNODES) { der(q,c,e) = q_d; }
173 }
174 if (vdim == 1 && (eval_flags & QI::DETERMINANTS))
175 {
176 det(q,e) = q_d;
177 }
178 }
179 }
180 }
181 });
182}
183
184// Template compute kernel for 2D quadrature interpolation:
185// * non-tensor product version,
186// * assumes 'e_vec' is using ElementDofOrdering::NATIVE,
187// * assumes 'maps.mode == FULL'.
188template<const int T_VDIM, const int T_ND, const int T_NQ>
189static void Eval2D(const int NE,
190 const int vdim,
191 const QVectorLayout q_layout,
192 const GeometricFactors *geom,
193 const DofToQuad &maps,
194 const Vector &e_vec,
195 Vector &q_val,
196 Vector &q_der,
197 Vector &q_det,
198 const int eval_flags)
199{
200 using QI = QuadratureInterpolator;
201
202 const int nd = maps.ndof;
203 const int nq = maps.nqpt;
204 const int ND = T_ND ? T_ND : nd;
205 const int NQ = T_NQ ? T_NQ : nq;
206 const int NMAX = NQ > ND ? NQ : ND;
207 const int VDIM = T_VDIM ? T_VDIM : vdim;
208 MFEM_ASSERT(maps.mode == DofToQuad::FULL, "internal error");
209 MFEM_ASSERT(!geom || geom->mesh->SpaceDimension() == 2, "");
210 MFEM_VERIFY(ND <= QI::MAX_ND2D, "");
211 MFEM_VERIFY(NQ <= QI::MAX_NQ2D, "");
212 MFEM_VERIFY(bool(geom) == bool(eval_flags & QI::PHYSICAL_DERIVATIVES),
213 "'geom' must be given (non-null) only when evaluating physical"
214 " derivatives");
215 const auto B = Reshape(maps.B.Read(), NQ, ND);
216 const auto G = Reshape(maps.G.Read(), NQ, 2, ND);
217 const auto J = Reshape(geom ? geom->J.Read() : nullptr, NQ, 2, 2, NE);
218 const auto E = Reshape(e_vec.Read(), ND, VDIM, NE);
219 auto val = q_layout == QVectorLayout::byNODES ?
220 Reshape(q_val.Write(), NQ, VDIM, NE):
221 Reshape(q_val.Write(), VDIM, NQ, NE);
222 auto der = q_layout == QVectorLayout::byNODES ?
223 Reshape(q_der.Write(), NQ, VDIM, 2, NE):
224 Reshape(q_der.Write(), VDIM, 2, NQ, NE);
225 auto det = Reshape(q_det.Write(), NQ, NE);
226 mfem::forall_2D(NE, NMAX, 1, [=] MFEM_HOST_DEVICE (int e)
227 {
228 const int ND = T_ND ? T_ND : nd;
229 const int NQ = T_NQ ? T_NQ : nq;
230 const int VDIM = T_VDIM ? T_VDIM : vdim;
231 constexpr int max_ND = T_ND ? T_ND : QI::MAX_ND2D;
232 constexpr int max_VDIM = T_VDIM ? T_VDIM : QI::MAX_VDIM2D;
233 MFEM_SHARED real_t s_E[max_VDIM*max_ND];
234 MFEM_FOREACH_THREAD(d, x, ND)
235 {
236 for (int c = 0; c < VDIM; c++)
237 {
238 s_E[c+d*VDIM] = E(d,c,e);
239 }
240 }
241 MFEM_SYNC_THREAD;
242
243 MFEM_FOREACH_THREAD(q, x, NQ)
244 {
245 if (eval_flags & QI::VALUES)
246 {
247 real_t ed[max_VDIM];
248 for (int c = 0; c < VDIM; c++) { ed[c] = 0.0; }
249 for (int d = 0; d < ND; ++d)
250 {
251 const real_t b = B(q,d);
252 for (int c = 0; c < VDIM; c++) { ed[c] += b*s_E[c+d*VDIM]; }
253 }
254 for (int c = 0; c < VDIM; c++)
255 {
256 if (q_layout == QVectorLayout::byVDIM) { val(c,q,e) = ed[c]; }
257 if (q_layout == QVectorLayout::byNODES) { val(q,c,e) = ed[c]; }
258 }
259 }
260 if ((eval_flags & QI::DERIVATIVES) ||
261 (eval_flags & QI::PHYSICAL_DERIVATIVES) ||
262 (eval_flags & QI::DETERMINANTS))
263 {
264 // use MAX_VDIM2D to avoid "subscript out of range" warnings
265 real_t D[QI::MAX_VDIM2D*2];
266 for (int i = 0; i < 2*VDIM; i++) { D[i] = 0.0; }
267 for (int d = 0; d < ND; ++d)
268 {
269 const real_t wx = G(q,0,d);
270 const real_t wy = G(q,1,d);
271 for (int c = 0; c < VDIM; c++)
272 {
273 real_t s_e = s_E[c+d*VDIM];
274 D[c+VDIM*0] += s_e * wx;
275 D[c+VDIM*1] += s_e * wy;
276 }
277 }
278 if (eval_flags & QI::DERIVATIVES)
279 {
280 for (int c = 0; c < VDIM; c++)
281 {
282 if (q_layout == QVectorLayout::byVDIM)
283 {
284 der(c,0,q,e) = D[c+VDIM*0];
285 der(c,1,q,e) = D[c+VDIM*1];
286 }
287 if (q_layout == QVectorLayout::byNODES)
288 {
289 der(q,c,0,e) = D[c+VDIM*0];
290 der(q,c,1,e) = D[c+VDIM*1];
291 }
292 }
293 }
294 if (eval_flags & QI::PHYSICAL_DERIVATIVES)
295 {
296 real_t Jloc[4], Jinv[4];
297 Jloc[0] = J(q,0,0,e);
298 Jloc[1] = J(q,1,0,e);
299 Jloc[2] = J(q,0,1,e);
300 Jloc[3] = J(q,1,1,e);
301 kernels::CalcInverse<2>(Jloc, Jinv);
302 for (int c = 0; c < VDIM; c++)
303 {
304 const real_t u = D[c+VDIM*0];
305 const real_t v = D[c+VDIM*1];
306 const real_t JiU = Jinv[0]*u + Jinv[1]*v;
307 const real_t JiV = Jinv[2]*u + Jinv[3]*v;
308 if (q_layout == QVectorLayout::byVDIM)
309 {
310 der(c,0,q,e) = JiU;
311 der(c,1,q,e) = JiV;
312 }
313 if (q_layout == QVectorLayout::byNODES)
314 {
315 der(q,c,0,e) = JiU;
316 der(q,c,1,e) = JiV;
317 }
318 }
319 }
320 if (eval_flags & QI::DETERMINANTS)
321 {
322 if (VDIM == 2) { det(q,e) = kernels::Det<2>(D); }
323 else
324 {
325 DeviceTensor<2> j(D, 3, 2);
326 const double E = j(0,0)*j(0,0) + j(1,0)*j(1,0) + j(2,0)*j(2,0);
327 const double F = j(0,0)*j(0,1) + j(1,0)*j(1,1) + j(2,0)*j(2,1);
328 const double G = j(0,1)*j(0,1) + j(1,1)*j(1,1) + j(2,1)*j(2,1);
329 det(q,e) = sqrt(E*G - F*F);
330 }
331 }
332 }
333 }
334 });
335}
336
337// Template compute kernel for 3D quadrature interpolation:
338// * non-tensor product version,
339// * assumes 'e_vec' is using ElementDofOrdering::NATIVE,
340// * assumes 'maps.mode == FULL'.
341template<const int T_VDIM, const int T_ND, const int T_NQ>
342static void Eval3D(const int NE,
343 const int vdim,
344 const QVectorLayout q_layout,
345 const GeometricFactors *geom,
346 const DofToQuad &maps,
347 const Vector &e_vec,
348 Vector &q_val,
349 Vector &q_der,
350 Vector &q_det,
351 const int eval_flags)
352{
353 using QI = QuadratureInterpolator;
354
355 const int nd = maps.ndof;
356 const int nq = maps.nqpt;
357 const int ND = T_ND ? T_ND : nd;
358 const int NQ = T_NQ ? T_NQ : nq;
359 const int NMAX = NQ > ND ? NQ : ND;
360 const int VDIM = T_VDIM ? T_VDIM : vdim;
361 MFEM_ASSERT(maps.mode == DofToQuad::FULL, "internal error");
362 MFEM_ASSERT(!geom || geom->mesh->SpaceDimension() == 3, "");
363 MFEM_VERIFY(ND <= QI::MAX_ND3D, "");
364 MFEM_VERIFY(NQ <= QI::MAX_NQ3D, "");
365 MFEM_VERIFY(VDIM == 3 || !(eval_flags & QI::DETERMINANTS), "");
366 MFEM_VERIFY(bool(geom) == bool(eval_flags & QI::PHYSICAL_DERIVATIVES),
367 "'geom' must be given (non-null) only when evaluating physical"
368 " derivatives");
369 const auto B = Reshape(maps.B.Read(), NQ, ND);
370 const auto G = Reshape(maps.G.Read(), NQ, 3, ND);
371 const auto J = Reshape(geom ? geom->J.Read() : nullptr, NQ, 3, 3, NE);
372 const auto E = Reshape(e_vec.Read(), ND, VDIM, NE);
373 auto val = q_layout == QVectorLayout::byNODES ?
374 Reshape(q_val.Write(), NQ, VDIM, NE):
375 Reshape(q_val.Write(), VDIM, NQ, NE);
376 auto der = q_layout == QVectorLayout::byNODES ?
377 Reshape(q_der.Write(), NQ, VDIM, 3, NE):
378 Reshape(q_der.Write(), VDIM, 3, NQ, NE);
379 auto det = Reshape(q_det.Write(), NQ, NE);
380 mfem::forall_2D(NE, NMAX, 1, [=] MFEM_HOST_DEVICE (int e)
381 {
382 const int ND = T_ND ? T_ND : nd;
383 const int NQ = T_NQ ? T_NQ : nq;
384 const int VDIM = T_VDIM ? T_VDIM : vdim;
385 constexpr int max_ND = T_ND ? T_ND : QI::MAX_ND3D;
386 constexpr int max_VDIM = T_VDIM ? T_VDIM : QI::MAX_VDIM3D;
387 MFEM_SHARED real_t s_E[max_VDIM*max_ND];
388 MFEM_FOREACH_THREAD(d, x, ND)
389 {
390 for (int c = 0; c < VDIM; c++)
391 {
392 s_E[c+d*VDIM] = E(d,c,e);
393 }
394 }
395 MFEM_SYNC_THREAD;
396
397 MFEM_FOREACH_THREAD(q, x, NQ)
398 {
399 if (eval_flags & QI::VALUES)
400 {
401 real_t ed[max_VDIM];
402 for (int c = 0; c < VDIM; c++) { ed[c] = 0.0; }
403 for (int d = 0; d < ND; ++d)
404 {
405 const real_t b = B(q,d);
406 for (int c = 0; c < VDIM; c++) { ed[c] += b*s_E[c+d*VDIM]; }
407 }
408 for (int c = 0; c < VDIM; c++)
409 {
410 if (q_layout == QVectorLayout::byVDIM) { val(c,q,e) = ed[c]; }
411 if (q_layout == QVectorLayout::byNODES) { val(q,c,e) = ed[c]; }
412 }
413 }
414 if ((eval_flags & QI::DERIVATIVES) ||
415 (eval_flags & QI::PHYSICAL_DERIVATIVES) ||
416 (eval_flags & QI::DETERMINANTS))
417 {
418 // use MAX_VDIM3D to avoid "subscript out of range" warnings
419 real_t D[QI::MAX_VDIM3D*3];
420 for (int i = 0; i < 3*VDIM; i++) { D[i] = 0.0; }
421 for (int d = 0; d < ND; ++d)
422 {
423 const real_t wx = G(q,0,d);
424 const real_t wy = G(q,1,d);
425 const real_t wz = G(q,2,d);
426 for (int c = 0; c < VDIM; c++)
427 {
428 real_t s_e = s_E[c+d*VDIM];
429 D[c+VDIM*0] += s_e * wx;
430 D[c+VDIM*1] += s_e * wy;
431 D[c+VDIM*2] += s_e * wz;
432 }
433 }
434 if (eval_flags & QI::DERIVATIVES)
435 {
436 for (int c = 0; c < VDIM; c++)
437 {
438 if (q_layout == QVectorLayout::byVDIM)
439 {
440 der(c,0,q,e) = D[c+VDIM*0];
441 der(c,1,q,e) = D[c+VDIM*1];
442 der(c,2,q,e) = D[c+VDIM*2];
443 }
444 if (q_layout == QVectorLayout::byNODES)
445 {
446 der(q,c,0,e) = D[c+VDIM*0];
447 der(q,c,1,e) = D[c+VDIM*1];
448 der(q,c,2,e) = D[c+VDIM*2];
449 }
450 }
451 }
452 if (eval_flags & QI::PHYSICAL_DERIVATIVES)
453 {
454 real_t Jloc[9], Jinv[9];
455 for (int col = 0; col < 3; col++)
456 {
457 for (int row = 0; row < 3; row++)
458 {
459 Jloc[row+3*col] = J(q,row,col,e);
460 }
461 }
462 kernels::CalcInverse<3>(Jloc, Jinv);
463 for (int c = 0; c < VDIM; c++)
464 {
465 const real_t u = D[c+VDIM*0];
466 const real_t v = D[c+VDIM*1];
467 const real_t w = D[c+VDIM*2];
468 const real_t JiU = Jinv[0]*u + Jinv[1]*v + Jinv[2]*w;
469 const real_t JiV = Jinv[3]*u + Jinv[4]*v + Jinv[5]*w;
470 const real_t JiW = Jinv[6]*u + Jinv[7]*v + Jinv[8]*w;
471 if (q_layout == QVectorLayout::byVDIM)
472 {
473 der(c,0,q,e) = JiU;
474 der(c,1,q,e) = JiV;
475 der(c,2,q,e) = JiW;
476 }
477 if (q_layout == QVectorLayout::byNODES)
478 {
479 der(q,c,0,e) = JiU;
480 der(q,c,1,e) = JiV;
481 der(q,c,2,e) = JiW;
482 }
483 }
484 }
485 if (VDIM == 3 && (eval_flags & QI::DETERMINANTS))
486 {
487 // The check (VDIM == 3) should eliminate this block when VDIM is
488 // known at compile time and (VDIM != 3).
489 det(q,e) = kernels::Det<3>(D);
490 }
491 }
492 }
493 });
494}
495
496} // namespace quadrature_interpolator
497
498} // namespace internal
499
500void QuadratureInterpolator::Mult(const Vector &e_vec,
501 unsigned eval_flags,
502 Vector &q_val,
503 Vector &q_der,
504 Vector &q_det) const
505{
506 using namespace internal::quadrature_interpolator;
507
508 const int ne = fespace->GetNE();
509 if (ne == 0) { return; }
510 const FiniteElement *fe = fespace->GetFE(0);
511
512 if (fe->GetMapType() == FiniteElement::MapType::H_DIV)
513 {
514 // q_der == q_div
515 return MultHDiv(e_vec, eval_flags, q_val, q_der);
516 }
517
518 const int vdim = fespace->GetVDim();
519 const int sdim = fespace->GetMesh()->SpaceDimension();
520 const bool use_tensor_eval =
521 use_tensor_products &&
522 dynamic_cast<const TensorBasisElement*>(fe) != nullptr;
523 const IntegrationRule *ir =
524 IntRule ? IntRule : &qspace->GetElementIntRule(0);
525 const DofToQuad::Mode mode =
526 use_tensor_eval ? DofToQuad::TENSOR : DofToQuad::FULL;
527 const DofToQuad &maps = fe->GetDofToQuad(*ir, mode);
528 const int dim = maps.FE->GetDim();
529 const int nd = maps.ndof;
530 const int nq = maps.nqpt;
531 const GeometricFactors *geom = nullptr;
532 if (eval_flags & PHYSICAL_DERIVATIVES)
533 {
534 const int jacobians = GeometricFactors::JACOBIANS;
535 geom = fespace->GetMesh()->GetGeometricFactors(*ir, jacobians);
536 }
537
538 MFEM_ASSERT(!(eval_flags & DETERMINANTS) || dim == vdim ||
539 (dim == 2 && vdim == 3), "Invalid dimensions for determinants.");
540 MFEM_ASSERT(fespace->GetMesh()->GetNumGeometries(
541 fespace->GetMesh()->Dimension()) == 1,
542 "mixed meshes are not supported");
543
544 if (use_tensor_eval)
545 {
546 if (eval_flags & VALUES)
547 {
548 TensorEvalKernels::Run(dim, q_layout, vdim, nd, nq, ne, maps.B.Read(),
549 e_vec.Read(), q_val.Write(), vdim, nd, nq);
550 }
551 if (eval_flags & (DERIVATIVES | PHYSICAL_DERIVATIVES))
552 {
553 const bool phys = (eval_flags & PHYSICAL_DERIVATIVES);
554 const real_t *J = phys ? geom->J.Read() : nullptr;
555 const int s_dim = phys ? sdim : dim;
556 GradKernels::Run(dim, q_layout, phys, vdim, nd, nq, ne,
557 maps.B.Read(), maps.G.Read(), J, e_vec.Read(),
558 q_der.Write(), s_dim, vdim, nd, nq);
559 }
560 if (eval_flags & DETERMINANTS)
561 {
562 DetKernels::Run(dim, vdim, nd, nq, ne, maps.B.Read(),
563 maps.G.Read(), e_vec.Read(), q_det.Write(), nd,
564 nq, &d_buffer);
565 }
566 }
567 else // use_tensor_eval == false
568 {
569 EvalKernels::Run(dim, vdim, maps.ndof, maps.nqpt, ne,vdim,q_layout,
570 geom, maps,e_vec, q_val,q_der,q_det,eval_flags);
571 }
572}
573
574void QuadratureInterpolator::MultHDiv(const Vector &e_vec,
575 unsigned eval_flags,
576 Vector &q_val,
577 Vector &q_div) const
578{
579 const int ne = fespace->GetNE();
580 if (ne == 0) { return; }
581 MFEM_VERIFY(fespace->IsVariableOrder() == false,
582 "variable order spaces are not supported yet!");
583 const FiniteElement *fe = fespace->GetFE(0);
584 MFEM_VERIFY(fe->GetMapType() == FiniteElement::MapType::H_DIV,
585 "this method can be used only for H(div) spaces");
586 MFEM_VERIFY((eval_flags &
587 ~(VALUES | PHYSICAL_VALUES | PHYSICAL_MAGNITUDES)) == 0,
588 "only VALUES, PHYSICAL_VALUES, and PHYSICAL_MAGNITUDES"
589 " evaluations are implemented!");
590 const int dim = fespace->GetMesh()->Dimension();
591 const int sdim = fespace->GetMesh()->SpaceDimension();
592 MFEM_VERIFY((dim == 2 || dim == 3) && dim == sdim,
593 "dim = " << dim << ", sdim = " << sdim
594 << " is not supported yet!");
595 MFEM_VERIFY(fespace->GetMesh()->GetNumGeometries(dim) <= 1,
596 "mixed meshes are not supported yet!");
597 const int vdim = fespace->GetVDim();
598 MFEM_VERIFY(vdim == 1, "vdim != 1 is not supported yet!");
599 auto tfe = dynamic_cast<const VectorTensorFiniteElement *>(fe);
600 MFEM_VERIFY(tfe != nullptr, "only quad and hex elements are supported!");
601 MFEM_VERIFY(use_tensor_products,
602 "non-tensor-product evaluation are not supported yet!");
603 const bool use_tensor_eval = use_tensor_products && (tfe != nullptr);
604 const IntegrationRule *ir =
605 IntRule ? IntRule : &qspace->GetElementIntRule(0);
606 const DofToQuad::Mode mode =
607 use_tensor_eval ? DofToQuad::TENSOR : DofToQuad::FULL;
608 const DofToQuad &maps_c = tfe->GetDofToQuad(*ir, mode);
609 const DofToQuad &maps_o = tfe->GetDofToQuadOpen(*ir, mode);
610 const int nd = maps_c.ndof;
611 const int nq = maps_c.nqpt;
612 const GeometricFactors *geom = nullptr;
613 if (eval_flags & (PHYSICAL_VALUES | PHYSICAL_MAGNITUDES))
614 {
615 const int jacobians = GeometricFactors::JACOBIANS;
616 geom = fespace->GetMesh()->GetGeometricFactors(*ir, jacobians);
617 }
618 // Check that at most one of VALUES, PHYSICAL_VALUES, and PHYSICAL_MAGNITUDES
619 // is specified:
620 MFEM_VERIFY(bool(eval_flags & VALUES) + bool(eval_flags & PHYSICAL_VALUES) +
621 bool(eval_flags & PHYSICAL_MAGNITUDES) <= 1,
622 "only one of VALUES, PHYSICAL_VALUES, and PHYSICAL_MAGNITUDES"
623 " can be requested at a time!");
624 const unsigned value_eval_mode =
625 eval_flags & (VALUES | PHYSICAL_VALUES | PHYSICAL_MAGNITUDES);
626 if (value_eval_mode)
627 {
628 // For PHYSICAL_MAGNITUDES the QVectorLayouts are the same and we
629 // instantiate only QVectorLayout::byNODES:
630 const auto q_l = (eval_flags & PHYSICAL_MAGNITUDES) ?
631 QVectorLayout::byNODES : q_layout;
632 TensorEvalHDivKernels::Run(
633 // dispatch params: dim + the template params of EvalHDiv2D/3D:
634 dim, q_l, value_eval_mode, nd, nq,
635 // runtime params, see the arguments of EvalHDiv2D/3D:
636 ne, maps_o.B.Read(), maps_c.B.Read(), geom ? geom->J.Read() : nullptr,
637 e_vec.Read(), q_val.Write(), nd, nq);
638 }
639 MFEM_CONTRACT_VAR(q_div);
640}
641
642void QuadratureInterpolator::MultTranspose(unsigned eval_flags,
643 const Vector &q_val,
644 const Vector &q_der,
645 Vector &e_vec) const
646{
647 MFEM_CONTRACT_VAR(eval_flags);
648 MFEM_CONTRACT_VAR(q_val);
649 MFEM_CONTRACT_VAR(q_der);
650 MFEM_CONTRACT_VAR(e_vec);
651 MFEM_ABORT("this method is not implemented yet");
652}
653
654void QuadratureInterpolator::Values(const Vector &e_vec,
655 Vector &q_val) const
656{
657 Vector empty;
658 Mult(e_vec, VALUES, q_val, empty, empty);
659}
660
661void QuadratureInterpolator::Derivatives(const Vector &e_vec,
662 Vector &q_der) const
663{
664 Vector empty;
665 Mult(e_vec, DERIVATIVES, empty, q_der, empty);
666}
667
668void QuadratureInterpolator::PhysDerivatives(const Vector &e_vec,
669 Vector &q_der) const
670{
671 Vector empty;
672 Mult(e_vec, PHYSICAL_DERIVATIVES, empty, q_der, empty);
673}
674
675void QuadratureInterpolator::Determinants(const Vector &e_vec,
676 Vector &q_det) const
677{
678 Vector empty;
679 Mult(e_vec, DETERMINANTS, empty, empty, q_det);
680}
681
682/// @cond Suppress_Doxygen_warnings
683
684namespace
685{
686
687using namespace internal::quadrature_interpolator;
688
689using EvalKernel = QuadratureInterpolator::EvalKernelType;
690using TensorEvalKernel = QuadratureInterpolator::TensorEvalKernelType;
691using GradKernel = QuadratureInterpolator::GradKernelType;
692using CollocatedGradKernel = QuadratureInterpolator::CollocatedGradKernelType;
693
694template <QVectorLayout Q_LAYOUT>
695TensorEvalKernel FallbackTensorEvalKernel(int DIM)
696{
697 if (DIM == 1) { return Values1D<Q_LAYOUT>; }
698 else if (DIM == 2) { return Values2D<Q_LAYOUT>; }
699 else if (DIM == 3) { return Values3D<Q_LAYOUT>; }
700 else { MFEM_ABORT(""); }
701}
702
703template<QVectorLayout Q_LAYOUT, bool GRAD_PHYS>
704GradKernel GetGradKernel(int DIM)
705{
706 if (DIM == 1) { return Derivatives1D<Q_LAYOUT, GRAD_PHYS>; }
707 else if (DIM == 2) { return Derivatives2D<Q_LAYOUT, GRAD_PHYS>; }
708 else if (DIM == 3) { return Derivatives3D<Q_LAYOUT, GRAD_PHYS>; }
709 else { MFEM_ABORT(""); }
710}
711
712
713template<QVectorLayout Q_LAYOUT>
714GradKernel GetGradKernel(int DIM, bool GRAD_PHYS)
715{
716 if (GRAD_PHYS) { return GetGradKernel<Q_LAYOUT, true>(DIM); }
717 else { return GetGradKernel<Q_LAYOUT, false>(DIM); }
718}
719
720template<QVectorLayout Q_LAYOUT, bool GRAD_PHYS>
721CollocatedGradKernel GetCollocatedGradKernel(int DIM)
722{
723 if (DIM == 1) { return CollocatedDerivatives1D<Q_LAYOUT, GRAD_PHYS>; }
724 else if (DIM == 2) { return CollocatedDerivatives2D<Q_LAYOUT, GRAD_PHYS>; }
725 else if (DIM == 3) { return CollocatedDerivatives3D<Q_LAYOUT, GRAD_PHYS>; }
726 else { MFEM_ABORT(""); }
727}
728
729template<QVectorLayout Q_LAYOUT>
730CollocatedGradKernel GetCollocatedGradKernel(int DIM, bool GRAD_PHYS)
731{
732 if (GRAD_PHYS) { return GetCollocatedGradKernel<Q_LAYOUT, true>(DIM); }
733 else { return GetCollocatedGradKernel<Q_LAYOUT, false>(DIM); }
734}
735} // namespace
736
737template <int DIM, int VDIM, int ND, int NQ>
738EvalKernel QuadratureInterpolator::EvalKernels::Kernel()
739{
740 using namespace internal::quadrature_interpolator;
741 if (DIM == 1) { return Eval1D; }
742 else if (DIM == 2) { return Eval2D<VDIM,ND,NQ>; }
743 else if (DIM == 3) { return Eval3D<VDIM,ND,NQ>; }
744 else { MFEM_ABORT(""); }
745}
746
747template <int DIM>
748EvalKernel GetEvalKernelVDimFallback(int VDIM)
749{
750 using EvalKernels = QuadratureInterpolator::EvalKernels;
751 if (VDIM == 1) { return EvalKernels::Kernel<DIM,1,0,0>(); }
752 else if (VDIM == 2) { return EvalKernels::Kernel<DIM,2,0,0>(); }
753 else if (VDIM == 3) { return EvalKernels::Kernel<DIM,3,0,0>(); }
754 else { MFEM_ABORT(""); }
755}
756
757EvalKernel QuadratureInterpolator::EvalKernels::Fallback(
758 int DIM, int VDIM, int ND, int NQ)
759{
760 if (DIM == 1) { return GetEvalKernelVDimFallback<1>(VDIM); }
761 else if (DIM == 2) { return GetEvalKernelVDimFallback<2>(VDIM); }
762 else if (DIM == 3) { return GetEvalKernelVDimFallback<3>(VDIM); }
763 else { MFEM_ABORT(""); }
764}
765
766TensorEvalKernel QuadratureInterpolator::TensorEvalKernels::Fallback(
767 int DIM, QVectorLayout Q_LAYOUT, int, int, int)
768{
769 if (Q_LAYOUT == QVectorLayout::byNODES) { return FallbackTensorEvalKernel<QVectorLayout::byNODES>(DIM); }
770 else { return FallbackTensorEvalKernel<QVectorLayout::byVDIM>(DIM); }
771}
772
773GradKernel QuadratureInterpolator::GradKernels::Fallback(
774 int DIM, QVectorLayout Q_LAYOUT, bool GRAD_PHYS, int, int, int)
775{
776 if (Q_LAYOUT == QVectorLayout::byNODES) { return GetGradKernel<QVectorLayout::byNODES>(DIM, GRAD_PHYS); }
777 else { return GetGradKernel<QVectorLayout::byVDIM>(DIM, GRAD_PHYS); }
778}
779
780CollocatedGradKernel QuadratureInterpolator::CollocatedGradKernels::Fallback(
781 int DIM, QVectorLayout Q_LAYOUT, bool GRAD_PHYS, int, int)
782{
783 if (Q_LAYOUT == QVectorLayout::byNODES) { return GetCollocatedGradKernel<QVectorLayout::byNODES>(DIM, GRAD_PHYS); }
784 else { return GetCollocatedGradKernel<QVectorLayout::byVDIM>(DIM, GRAD_PHYS); }
785}
786
787/// @endcond
788
789namespace internal
790{
791namespace quadrature_interpolator
792{
793void InitEvalKernels()
794{
795 using k = QuadratureInterpolator::EvalKernels;
796 // 2D, VDIM = 1
797 k::Specialization<2,1,1,1>::Add();
798 k::Specialization<2,1,1,4>::Add();
799 // Q1
800 k::Specialization<2,1,4,4>::Add();
801 k::Specialization<2,1,4,9>::Add();
802 // Q2
803 k::Specialization<2,1,9,9>::Add();
804 k::Specialization<2,1,9,16>::Add();
805 // Q3
806 k::Specialization<2,1,16,16>::Add();
807 k::Specialization<2,1,16,25>::Add();
808 k::Specialization<2,1,16,36>::Add();
809 // Q4
810 k::Specialization<2,1,25,25>::Add();
811 k::Specialization<2,1,25,36>::Add();
812 k::Specialization<2,1,25,49>::Add();
813 k::Specialization<2,1,25,64>::Add();
814
815 // 3D, VDIM = 1
816 // Q0
817 k::Specialization<3,1,1,1>::Add();
818 k::Specialization<3,1,1,8>::Add();
819 // Q1
820 k::Specialization<3,1,8,8>::Add();
821 k::Specialization<3,1,8,27>::Add();
822 // Q2
823 k::Specialization<3,1,27,27>::Add();
824 k::Specialization<3,1,27,64>::Add();
825 // Q3
826 k::Specialization<3,1,64,64>::Add();
827 k::Specialization<3,1,64,125>::Add();
828 k::Specialization<3,1,64,216>::Add();
829 // Q4
830 k::Specialization<3,1,125,125>::Add();
831 k::Specialization<3,1,125,216>::Add();
832
833 // 2D, VDIM = 3
834 // Q0
835 k::Specialization<2,3,1,1>::Add();
836 k::Specialization<2,3,1,4>::Add();
837 // Q1
838 k::Specialization<2,3,4,4>::Add();
839 k::Specialization<2,3,4,9>::Add();
840 // Q2
841 k::Specialization<2,3,9,4>::Add();
842 k::Specialization<2,3,9,9>::Add();
843 k::Specialization<2,3,9,16>::Add();
844 k::Specialization<2,3,9,25>::Add();
845 // Q3
846 k::Specialization<2,3,16,16>::Add();
847 k::Specialization<2,3,16,25>::Add();
848 k::Specialization<2,3,16,36>::Add();
849 // Q4
850 k::Specialization<2,3,25,25>::Add();
851 k::Specialization<2,3,25,36>::Add();
852 k::Specialization<2,3,25,49>::Add();
853 k::Specialization<2,3,25,64>::Add();
854
855 // 2D, VDIM = 2
856 // Q1
857 k::Specialization<2,2,4,4>::Add();
858 k::Specialization<2,2,4,9>::Add();
859 // Q2
860 k::Specialization<2,2,9,9>::Add();
861 k::Specialization<2,2,9,16>::Add();
862 // Q3
863 k::Specialization<2,2,16,16>::Add();
864 k::Specialization<2,2,16,25>::Add();
865 k::Specialization<2,2,16,36>::Add();
866 // Q4
867 k::Specialization<2,2,25,25>::Add();
868 k::Specialization<2,2,25,36>::Add();
869 k::Specialization<2,2,25,49>::Add();
870 k::Specialization<2,2,25,64>::Add();
871
872 // 3D, VDIM = 3
873 // Q1
874 k::Specialization<3,3,8,8>::Add();
875 k::Specialization<3,3,8,27>::Add();
876 // Q2
877 k::Specialization<3,3,27,27>::Add();
878 k::Specialization<3,3,27,64>::Add();
879 k::Specialization<3,3,27,125>::Add();
880 // Q3
881 k::Specialization<3,3,64,64>::Add();
882 k::Specialization<3,3,64,125>::Add();
883 k::Specialization<3,3,64,216>::Add();
884 // Q4
885 k::Specialization<3,3,125,125>::Add();
886 k::Specialization<3,3,125,216>::Add();
887}
888
889} // namespace quadrature_Interpolator
890} // namespace internal
891
892} // namespace mfem
mfem::Array::Read
const T * Read(bool on_dev=true) const
Shortcut for mfem::Read(a.GetMemory(), a.Size(), on_dev).
Definition array.hpp:337
mfem::DofToQuad
Structure representing the matrices/tensors needed to evaluate (in reference space) the values,...
Definition fe_base.hpp:141
mfem::DofToQuad::mode
Mode mode
Describes the contents of the B, Bt, G, and Gt arrays, see Mode.
Definition fe_base.hpp:174
mfem::DofToQuad::G
Array< real_t > G
Gradients/divergences/curls of basis functions evaluated at quadrature points.
Definition fe_base.hpp:214
mfem::DofToQuad::Mode
Mode
Type of data stored in the arrays B, Bt, G, and Gt.
Definition fe_base.hpp:154
mfem::DofToQuad::FULL
@ FULL
Full multidimensional representation which does not use tensor product structure. The ordering of the...
Definition fe_base.hpp:158
mfem::DofToQuad::TENSOR
@ TENSOR
Tensor product representation using 1D matrices/tensors with dimensions using 1D number of quadrature...
Definition fe_base.hpp:165
mfem::DofToQuad::B
Array< real_t > B
Basis functions evaluated at quadrature points.
Definition fe_base.hpp:193
mfem::DofToQuad::ndof
int ndof
Number of degrees of freedom = number of basis functions. When mode is TENSOR, this is the 1D number.
Definition fe_base.hpp:178
mfem::DofToQuad::FE
const class FiniteElement * FE
The FiniteElement that created and owns this object.
Definition fe_base.hpp:145
mfem::DofToQuad::nqpt
int nqpt
Number of quadrature points. When mode is TENSOR, this is the 1D number.
Definition fe_base.hpp:182
mfem::FiniteElementSpace
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition fespace.hpp:244
mfem::FiniteElementSpace::IsVariableOrder
bool IsVariableOrder() const
Returns true if the space contains elements of varying polynomial orders.
Definition fespace.hpp:709
mfem::FiniteElementSpace::GetFE
virtual const FiniteElement * GetFE(int i) const
Returns pointer to the FiniteElement in the FiniteElementCollection associated with i'th element in t...
Definition fespace.cpp:3835
mfem::FiniteElementSpace::GetNE
int GetNE() const
Returns number of elements in the mesh.
Definition fespace.hpp:891
mfem::FiniteElementSpace::GetMesh
Mesh * GetMesh() const
Returns the mesh.
Definition fespace.hpp:679
mfem::FiniteElementSpace::GetVDim
int GetVDim() const
Returns the vector dimension of the finite element space.
Definition fespace.hpp:841
mfem::FiniteElementSpace::GetTypicalFE
const FiniteElement * GetTypicalFE() const
Return GetFE(0) if the local mesh is not empty; otherwise return a typical FE based on the Geometry t...
Definition fespace.cpp:3871
mfem::FiniteElement
Abstract class for all finite elements.
Definition fe_base.hpp:244
mfem::FiniteElement::GetDofToQuad
virtual const DofToQuad & GetDofToQuad(const IntegrationRule &ir, DofToQuad::Mode mode) const
Return a DofToQuad structure corresponding to the given IntegrationRule using the given DofToQuad::Mo...
Definition fe_base.cpp:365
mfem::FiniteElement::GetDim
int GetDim() const
Returns the reference space dimension for the finite element.
Definition fe_base.hpp:321
mfem::FiniteElement::GetMapType
int GetMapType() const
Returns the FiniteElement::MapType of the element describing how reference functions are mapped to ph...
Definition fe_base.hpp:360
mfem::GeometricFactors
Structure for storing mesh geometric factors: coordinates, Jacobians, and determinants of the Jacobia...
Definition mesh.hpp:2937
mfem::GeometricFactors::mesh
const Mesh * mesh
Definition mesh.hpp:2943
mfem::GeometricFactors::J
Vector J
Jacobians of the element transformations at all quadrature points.
Definition mesh.hpp:2976
mfem::IntegrationRule
Class for an integration rule - an Array of IntegrationPoint.
Definition intrules.hpp:100
mfem::Mesh::Dimension
int Dimension() const
Dimension of the reference space used within the elements.
Definition mesh.hpp:1216
mfem::Mesh::SpaceDimension
int SpaceDimension() const
Dimension of the physical space containing the mesh.
Definition mesh.hpp:1219
mfem::Mesh::GetGeometricFactors
const GeometricFactors * GetGeometricFactors(const IntegrationRule &ir, const int flags, MemoryType d_mt=MemoryType::DEFAULT)
Return the mesh geometric factors corresponding to the given integration rule.
Definition mesh.cpp:880
mfem::Mesh::GetNumGeometries
int GetNumGeometries(int dim) const
Return the number of geometries of the given dimension present in the mesh.
Definition mesh.cpp:7243
mfem::QuadratureInterpolator
A class that performs interpolation from an E-vector to quadrature point values and/or derivatives (Q...
Definition quadinterpolator.hpp:31
mfem::QuadratureInterpolator::CollocatedGradKernelType
void(*)(const int, const real_t *, const real_t *, const real_t *, real_t *, const int, const int, const int) CollocatedGradKernelType
Definition quadinterpolator.hpp:160
mfem::QuadratureInterpolator::use_tensor_products
bool use_tensor_products
Tensor product evaluation mode.
Definition quadinterpolator.hpp:40
mfem::QuadratureInterpolator::VALUES
@ VALUES
Evaluate the values at quadrature points.
Definition quadinterpolator.hpp:58
mfem::QuadratureInterpolator::DERIVATIVES
@ DERIVATIVES
Evaluate the derivatives at quadrature points.
Definition quadinterpolator.hpp:59
mfem::QuadratureInterpolator::PHYSICAL_DERIVATIVES
@ PHYSICAL_DERIVATIVES
Evaluate the physical derivatives.
Definition quadinterpolator.hpp:64
mfem::QuadratureInterpolator::DETERMINANTS
@ DETERMINANTS
Assuming the derivative at quadrature points form a matrix, this flag can be used to compute and stor...
Definition quadinterpolator.hpp:63
mfem::QuadratureInterpolator::EvalKernelType
void(*)(const int, const int, const QVectorLayout, const GeometricFactors *, const DofToQuad &, const Vector &, Vector &, Vector &, Vector &, const int) EvalKernelType
Definition quadinterpolator.hpp:167
mfem::QuadratureInterpolator::QuadratureInterpolator
QuadratureInterpolator(const FiniteElementSpace &fes, const IntegrationRule &ir)
Definition quadinterpolator.cpp:59
mfem::QuadratureInterpolator::GradKernelType
void(*)(const int, const real_t *, const real_t *, const real_t *, const real_t *, real_t *, const int, const int, const int, const int) GradKernelType
Definition quadinterpolator.hpp:157
mfem::QuadratureInterpolator::Mult
void Mult(const Vector &e_vec, unsigned eval_flags, Vector &q_val, Vector &q_der, Vector &q_det) const
Interpolate the E-vector e_vec to quadrature points.
Definition quadinterpolator.cpp:500
mfem::QuadratureInterpolator::Determinants
void Determinants(const Vector &e_vec, Vector &q_det) const
Compute the determinants of the derivatives (with respect to reference coordinates) of the E-vector e...
Definition quadinterpolator.cpp:675
mfem::QuadratureInterpolator::MultTranspose
void MultTranspose(unsigned eval_flags, const Vector &q_val, const Vector &q_der, Vector &e_vec) const
Perform the transpose operation of Mult(). (TODO)
Definition quadinterpolator.cpp:642
mfem::QuadratureInterpolator::TensorEvalKernelType
void(*)(const int, const real_t *, const real_t *, real_t *, const int, const int, const int) TensorEvalKernelType
Definition quadinterpolator.hpp:155
mfem::QuadratureInterpolator::Values
void Values(const Vector &e_vec, Vector &q_val) const
Interpolate the values of the E-vector e_vec at quadrature points.
Definition quadinterpolator.cpp:654
mfem::QuadratureInterpolator::Derivatives
void Derivatives(const Vector &e_vec, Vector &q_der) const
Interpolate the derivatives (with respect to reference coordinates) of the E-vector e_vec at quadratu...
Definition quadinterpolator.cpp:661
mfem::QuadratureInterpolator::q_layout
QVectorLayout q_layout
Output Q-vector layout.
Definition quadinterpolator.hpp:38
mfem::QuadratureInterpolator::PhysDerivatives
void PhysDerivatives(const Vector &e_vec, Vector &q_der) const
Interpolate the derivatives in physical space of the E-vector e_vec at quadrature points.
Definition quadinterpolator.cpp:668
mfem::QuadratureInterpolator::IntRule
const IntegrationRule * IntRule
Not owned.
Definition quadinterpolator.hpp:37
mfem::QuadratureInterpolator::d_buffer
Vector d_buffer
Auxiliary device buffer.
Definition quadinterpolator.hpp:41
mfem::QuadratureInterpolator::fespace
const FiniteElementSpace * fespace
Not owned.
Definition quadinterpolator.hpp:35
mfem::QuadratureInterpolator::qspace
const QuadratureSpace * qspace
Not owned.
Definition quadinterpolator.hpp:36
mfem::QuadratureInterpolator::MultHDiv
void MultHDiv(const Vector &e_vec, unsigned eval_flags, Vector &q_val, Vector &q_div) const
Auxiliary method called by Mult() when using H(div)-conforming space.
Definition quadinterpolator.cpp:574
mfem::QuadratureSpace
Class representing the storage layout of a QuadratureFunction.
Definition qspace.hpp:120
mfem::QuadratureSpace::GetElementIntRule
const IntegrationRule & GetElementIntRule(int idx) const
Get the IntegrationRule associated with mesh element idx.
Definition qspace.hpp:149
mfem::Vector
Vector data type.
Definition vector.hpp:82
mfem::Vector::Read
virtual const real_t * Read(bool on_dev=true) const
Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), on_dev).
Definition vector.hpp:494
mfem::Vector::UseDevice
virtual void UseDevice(bool use_dev) const
Enable execution of Vector operations using the mfem::Device.
Definition vector.hpp:144
mfem::Vector::Write
virtual real_t * Write(bool on_dev=true)
Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), on_dev).
Definition vector.hpp:502
dim
int dim
Definition ex24.cpp:53
b
real_t b
Definition lissajous.cpp:42
DIM
constexpr int DIM
Definition minimal-surface.cpp:72
mfem::kernels::CalcInverse
MFEM_HOST_DEVICE void CalcInverse(const T *data, T *inv_data)
Return the inverse of a matrix with given size and data into the matrix with data inv_data.
Definition kernels.hpp:246
mfem::kernels::Det
MFEM_HOST_DEVICE T Det(const T *data)
Compute the determinant of a square matrix of size dim with given data.
Definition kernels.hpp:237
mfem::u
real_t u(const Vector &xvec)
Definition lor_mms.hpp:22
mfem::Write
T * Write(Memory< T > &mem, int size, bool on_dev=true)
Get a pointer for write access to mem with the mfem::Device's DeviceMemoryClass, if on_dev = true,...
Definition device.hpp:358
mfem::Reshape
MFEM_HOST_DEVICE DeviceTensor< sizeof...(Dims), T > Reshape(T *ptr, Dims... dims)
Wrap a pointer as a DeviceTensor with automatically deduced template parameters.
Definition dtensor.hpp:131
mfem::forall_2D
void forall_2D(int N, int X, int Y, lambda &&body)
Definition forall.hpp:762
mfem::UsesTensorBasis
bool UsesTensorBasis(const FiniteElementSpace &fes)
Return true if the mesh contains only one topology and the elements are tensor elements.
Definition fespace.hpp:1555
mfem::QVectorLayout
QVectorLayout
Type describing possible layouts for Q-vectors.
Definition fespace.hpp:53
mfem::real_t
float real_t
Definition config.hpp:43
mfem::forall
void forall(int N, lambda &&body)
Definition forall.hpp:753
Generated on Wed Apr 9 2025 16:53:04 for MFEM by doxygen 1.11.0