MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
bilininteg_hdiv_kernels.hpp
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#ifndef MFEM_BILININTEG_HDIV_KERNELS_HPP
13#define MFEM_BILININTEG_HDIV_KERNELS_HPP
14
15#include "../../config/config.hpp"
16#include "../../general/array.hpp"
17#include "../../general/forall.hpp"
18#include "../../linalg/dtensor.hpp"
19#include "../../linalg/vector.hpp"
20#include "../bilininteg.hpp"
21
22// Piola transformation in H(div): w = (1 / det (dF)) dF \hat{w}
23// div w = (1 / det (dF)) \hat{div} \hat{w}
24
25namespace mfem
26{
27
28namespace internal
29{
30
31// PA H(div) Mass Assemble 2D kernel
32void PAHdivMassSetup2D(const int Q1D,
33 const int coeffDim,
34 const int NE,
35 const Array<real_t> &w,
36 const Vector &j,
37 Vector &coeff_,
38 Vector &op);
39
40// PA H(div) Mass Assemble 3D kernel
41void PAHdivMassSetup3D(const int Q1D,
42 const int coeffDim,
43 const int NE,
44 const Array<real_t> &w,
45 const Vector &j,
46 Vector &coeff_,
47 Vector &op);
48
49// PA H(div) Mass Diagonal 2D kernel
50void PAHdivMassAssembleDiagonal2D(const int D1D,
51 const int Q1D,
52 const int NE,
53 const bool symmetric,
54 const Array<real_t> &Bo_,
55 const Array<real_t> &Bc_,
56 const Vector &op_,
57 Vector &diag_);
58
59// PA H(div) Mass Diagonal 3D kernel
60void PAHdivMassAssembleDiagonal3D(const int D1D,
61 const int Q1D,
62 const int NE,
63 const bool symmetric,
64 const Array<real_t> &Bo_,
65 const Array<real_t> &Bc_,
66 const Vector &op_,
67 Vector &diag_);
68
69void PAHdivMassApply(const int dim,
70 const int D1D,
71 const int Q1D,
72 const int NE,
73 const bool symmetric,
74 const Array<real_t> &Bo,
75 const Array<real_t> &Bc,
76 const Array<real_t> &Bot,
77 const Array<real_t> &Bct,
78 const Vector &op,
79 const Vector &x,
80 Vector &y);
81
82// PA H(div) Mass Apply 2D kernel
83void PAHdivMassApply2D(const int D1D,
84 const int Q1D,
85 const int NE,
86 const bool symmetric,
87 const Array<real_t> &Bo_,
88 const Array<real_t> &Bc_,
89 const Array<real_t> &Bot_,
90 const Array<real_t> &Bct_,
91 const Vector &op_,
92 const Vector &x_,
93 Vector &y_);
94
95// PA H(div) Mass Apply 3D kernel
96void PAHdivMassApply3D(const int D1D,
97 const int Q1D,
98 const int NE,
99 const bool symmetric,
100 const Array<real_t> &Bo_,
101 const Array<real_t> &Bc_,
102 const Array<real_t> &Bot_,
103 const Array<real_t> &Bct_,
104 const Vector &op_,
105 const Vector &x_,
106 Vector &y_);
107
108// Shared memory PA H(div) Mass Apply 2D kernel
109template<int T_D1D = 0, int T_Q1D = 0>
110inline void SmemPAHdivMassApply2D(const int NE,
111 const bool symmetric,
112 const Array<real_t> &Bo_,
113 const Array<real_t> &Bc_,
114 const Array<real_t> &Bot_,
115 const Array<real_t> &Bct_,
116 const Vector &op_,
117 const Vector &x_,
118 Vector &y_,
119 const int d1d = 0,
120 const int q1d = 0)
121{
122 MFEM_CONTRACT_VAR(Bot_);
123 MFEM_CONTRACT_VAR(Bct_);
124
125 static constexpr int VDIM = 2;
126
127 const int D1D = T_D1D ? T_D1D : d1d;
128 const int Q1D = T_Q1D ? T_Q1D : q1d;
129
130 const auto bo = Reshape(Bo_.Read(), Q1D, D1D-1);
131 const auto bc = Reshape(Bc_.Read(), Q1D, D1D);
132 const auto D = Reshape(op_.Read(), Q1D, Q1D, symmetric ? 3 : 4, NE);
133 const auto x = Reshape(x_.Read(), D1D*(D1D-1), VDIM, NE);
134 auto y = y_.ReadWrite();
135
136 mfem::forall_3D(NE, Q1D, Q1D, VDIM, [=] MFEM_HOST_DEVICE (int e)
137 {
138 const int tidz = MFEM_THREAD_ID(z);
139
140 const int D1D = T_D1D ? T_D1D : d1d;
141 const int Q1D = T_Q1D ? T_Q1D : q1d;
142
143 constexpr int MQ1 = T_Q1D ? T_Q1D : DofQuadLimits::HDIV_MAX_Q1D;
144 constexpr int MD1 = T_D1D ? T_D1D : DofQuadLimits::HDIV_MAX_D1D;
145 constexpr int MDQ = (MQ1 > MD1) ? MQ1 : MD1;
146
147 MFEM_SHARED real_t smo[MQ1*(MD1-1)];
148 DeviceMatrix Bo(smo, D1D-1, Q1D);
149
150 MFEM_SHARED real_t smc[MQ1*MD1];
151 DeviceMatrix Bc(smc, D1D, Q1D);
152
153 MFEM_SHARED real_t sm0[VDIM*MDQ*MDQ];
154 MFEM_SHARED real_t sm1[VDIM*MDQ*MDQ];
155 DeviceMatrix X(sm0, D1D*(D1D-1), VDIM);
156 DeviceCube QD(sm1, Q1D, D1D, VDIM);
157 DeviceCube QQ(sm0, Q1D, Q1D, VDIM);
158 DeviceCube DQ(sm1, D1D, Q1D, VDIM);
159
160 // Load X, Bo and Bc into shared memory
161 MFEM_FOREACH_THREAD(vd,z,VDIM)
162 {
163 MFEM_FOREACH_THREAD(dy,y,D1D)
164 {
165 MFEM_FOREACH_THREAD(qx,x,Q1D)
166 {
167 if (qx < D1D && dy < (D1D-1))
168 {
169 X(qx + dy*D1D,vd) = x(qx+dy*D1D,vd,e);
170 }
171 if (tidz == 0)
172 {
173 if (dy < (D1D-1)) { Bo(dy,qx) = bo(qx,dy); }
174 Bc(dy,qx) = bc(qx,dy);
175 }
176 }
177 }
178 }
179 MFEM_SYNC_THREAD;
180 // Apply B operator
181 MFEM_FOREACH_THREAD(vd,z,VDIM)
182 {
183 const int nx = (vd == 0) ? D1D : D1D-1;
184 const int ny = (vd == 1) ? D1D : D1D-1;
185 DeviceCube Xxy(X, nx, ny, VDIM);
186 DeviceMatrix Bx = (vd == 0) ? Bc : Bo;
187 MFEM_FOREACH_THREAD(dy,y,ny)
188 {
189 MFEM_FOREACH_THREAD(qx,x,Q1D)
190 {
191 real_t dq = 0.0;
192 for (int dx = 0; dx < nx; ++dx)
193 {
194 dq += Xxy(dx,dy,vd) * Bx(dx,qx);
195 }
196 QD(qx,dy,vd) = dq;
197 }
198 }
199 }
200 MFEM_SYNC_THREAD;
201 MFEM_FOREACH_THREAD(vd,z,VDIM)
202 {
203 const int ny = (vd == 1) ? D1D : D1D-1;
204 DeviceMatrix By = (vd == 1) ? Bc : Bo;
205 MFEM_FOREACH_THREAD(qy,y,Q1D)
206 {
207 MFEM_FOREACH_THREAD(qx,x,Q1D)
208 {
209 real_t qq = 0.0;
210 for (int dy = 0; dy < ny; ++dy)
211 {
212 qq += QD(qx,dy,vd) * By(dy,qy);
213 }
214 QQ(qx,qy,vd) = qq;
215 }
216 }
217 }
218 MFEM_SYNC_THREAD;
219 // Apply D operator
220 if (tidz == 0)
221 {
222 MFEM_FOREACH_THREAD(qy,y,Q1D)
223 {
224 MFEM_FOREACH_THREAD(qx,x,Q1D)
225 {
226 const real_t Qx = QQ(qx,qy,0);
227 const real_t Qy = QQ(qx,qy,1);
228
229 const real_t D11 = D(qx,qy,0,e);
230 const real_t D12 = D(qx,qy,1,e);
231 const real_t D21 = symmetric ? D12 : D(qx,qy,2,e);
232 const real_t D22 = symmetric ? D(qx,qy,2,e) : D(qx,qy,3,e);
233
234 QQ(qx,qy,0) = D11*Qx + D12*Qy;
235 QQ(qx,qy,1) = D21*Qx + D22*Qy;
236 }
237 }
238 }
239 MFEM_SYNC_THREAD;
240 // Apply Bt operator
241 MFEM_FOREACH_THREAD(vd,z,VDIM)
242 {
243 const int nx = (vd == 0) ? D1D : D1D-1;
244 DeviceMatrix Btx = (vd == 0) ? Bc : Bo;
245 MFEM_FOREACH_THREAD(qy,y,Q1D)
246 {
247 MFEM_FOREACH_THREAD(dx,x,nx)
248 {
249 real_t qd = 0.0;
250 for (int qx = 0; qx < Q1D; ++qx)
251 {
252 qd += QQ(qx,qy,vd) * Btx(dx,qx);
253 }
254 DQ(dx,qy,vd) = qd;
255 }
256 }
257 }
258 MFEM_SYNC_THREAD;
259 MFEM_FOREACH_THREAD(vd,z,VDIM)
260 {
261 const int nx = (vd == 0) ? D1D : D1D-1;
262 const int ny = (vd == 1) ? D1D : D1D-1;
263 DeviceMatrix Bty = (vd == 1) ? Bc : Bo;
264 DeviceTensor<4> Yxy(y, nx, ny, VDIM, NE);
265 MFEM_FOREACH_THREAD(dy,y,ny)
266 {
267 MFEM_FOREACH_THREAD(dx,x,nx)
268 {
269 real_t dd = 0.0;
270 for (int qy = 0; qy < Q1D; ++qy)
271 {
272 dd += DQ(dx,qy,vd) * Bty(dy,qy);
273 }
274 Yxy(dx,dy,vd,e) += dd;
275 }
276 }
277 }
278 MFEM_SYNC_THREAD;
279 });
280}
281
282// Shared memory PA H(div) Mass Apply 3D kernel
283template<int T_D1D = 0, int T_Q1D = 0>
284inline void SmemPAHdivMassApply3D(const int NE,
285 const bool symmetric,
286 const Array<real_t> &Bo_,
287 const Array<real_t> &Bc_,
288 const Array<real_t> &Bot_,
289 const Array<real_t> &Bct_,
290 const Vector &op_,
291 const Vector &x_,
292 Vector &y_,
293 const int d1d = 0,
294 const int q1d = 0)
295{
296 MFEM_CONTRACT_VAR(Bot_);
297 MFEM_CONTRACT_VAR(Bct_);
298
299 static constexpr int VDIM = 3;
300
301 const int D1D = T_D1D ? T_D1D : d1d;
302 const int Q1D = T_Q1D ? T_Q1D : q1d;
303
304 const auto bo = Reshape(Bo_.Read(), Q1D, D1D-1);
305 const auto bc = Reshape(Bc_.Read(), Q1D, D1D);
306 const auto D = Reshape(op_.Read(), Q1D, Q1D, Q1D, symmetric ? 6 : 9, NE);
307 const auto x = Reshape(x_.Read(), D1D*(D1D-1)*(D1D-1), VDIM, NE);
308 auto y = y_.ReadWrite();
309
310 mfem::forall_3D(NE, Q1D, Q1D, VDIM, [=] MFEM_HOST_DEVICE (int e)
311 {
312 const int tidz = MFEM_THREAD_ID(z);
313
314 const int D1D = T_D1D ? T_D1D : d1d;
315 const int Q1D = T_Q1D ? T_Q1D : q1d;
316
317 constexpr int MQ1 = T_Q1D ? T_Q1D : DofQuadLimits::HDIV_MAX_Q1D;
318 constexpr int MD1 = T_D1D ? T_D1D : DofQuadLimits::HDIV_MAX_D1D;
319 constexpr int MDQ = (MQ1 > MD1) ? MQ1 : MD1;
320
321 MFEM_SHARED real_t smo[MQ1*(MD1-1)];
322 DeviceMatrix Bo(smo, D1D-1, Q1D);
323
324 MFEM_SHARED real_t smc[MQ1*MD1];
325 DeviceMatrix Bc(smc, D1D, Q1D);
326
327 MFEM_SHARED real_t sm0[VDIM*MDQ*MDQ*MDQ];
328 MFEM_SHARED real_t sm1[VDIM*MDQ*MDQ*MDQ];
329 DeviceMatrix X(sm0, D1D*(D1D-1)*(D1D-1), VDIM);
330 DeviceTensor<4> QDD(sm1, Q1D, D1D, D1D, VDIM);
331 DeviceTensor<4> QQD(sm0, Q1D, Q1D, D1D, VDIM);
332 DeviceTensor<4> QQQ(sm1, Q1D, Q1D, Q1D, VDIM);
333 DeviceTensor<4> DQQ(sm0, D1D, Q1D, Q1D, VDIM);
334 DeviceTensor<4> DDQ(sm1, D1D, D1D, Q1D, VDIM);
335
336 // Load X into shared memory
337 MFEM_FOREACH_THREAD(vd,z,VDIM)
338 {
339 MFEM_FOREACH_THREAD(dz,y,D1D-1)
340 {
341 MFEM_FOREACH_THREAD(dy,x,D1D-1)
342 {
343 MFEM_UNROLL(MD1)
344 for (int dx = 0; dx < D1D; ++dx)
345 {
346 X(dx+(dy+dz*(D1D-1))*D1D,vd) = x(dx+(dy+dz*(D1D-1))*D1D,vd,e);
347 }
348 }
349 }
350 }
351 // Load Bo and Bc into shared memory
352 if (tidz == 0)
353 {
354 MFEM_FOREACH_THREAD(d,y,D1D-1)
355 {
356 MFEM_FOREACH_THREAD(q,x,Q1D)
357 {
358 Bo(d,q) = bo(q,d);
359 }
360 }
361 MFEM_FOREACH_THREAD(d,y,D1D)
362 {
363 MFEM_FOREACH_THREAD(q,x,Q1D)
364 {
365 Bc(d,q) = bc(q,d);
366 }
367 }
368 }
369 MFEM_SYNC_THREAD;
370 // Apply B operator
371 MFEM_FOREACH_THREAD(vd,z,VDIM)
372 {
373 const int nx = (vd == 0) ? D1D : D1D-1;
374 const int ny = (vd == 1) ? D1D : D1D-1;
375 const int nz = (vd == 2) ? D1D : D1D-1;
376 DeviceTensor<4> Xxyz(X, nx, ny, nz, VDIM);
377 DeviceMatrix Bx = (vd == 0) ? Bc : Bo;
378 MFEM_FOREACH_THREAD(dy,y,ny)
379 {
380 MFEM_FOREACH_THREAD(qx,x,Q1D)
381 {
382 real_t u[D1D];
383 MFEM_UNROLL(MD1)
384 for (int dz = 0; dz < nz; ++dz) { u[dz] = 0.0; }
385 MFEM_UNROLL(MD1)
386 for (int dx = 0; dx < nx; ++dx)
387 {
388 MFEM_UNROLL(MD1)
389 for (int dz = 0; dz < nz; ++dz)
390 {
391 u[dz] += Xxyz(dx,dy,dz,vd) * Bx(dx,qx);
392 }
393 }
394 MFEM_UNROLL(MD1)
395 for (int dz = 0; dz < nz; ++dz) { QDD(qx,dy,dz,vd) = u[dz]; }
396 }
397 }
398 }
399 MFEM_SYNC_THREAD;
400 MFEM_FOREACH_THREAD(vd,z,VDIM)
401 {
402 const int ny = (vd == 1) ? D1D : D1D-1;
403 const int nz = (vd == 2) ? D1D : D1D-1;
404 DeviceMatrix By = (vd == 1) ? Bc : Bo;
405 MFEM_FOREACH_THREAD(qy,y,Q1D)
406 {
407 MFEM_FOREACH_THREAD(qx,x,Q1D)
408 {
409 real_t u[D1D];
410 MFEM_UNROLL(MD1)
411 for (int dz = 0; dz < nz; ++dz) { u[dz] = 0.0; }
412 MFEM_UNROLL(MD1)
413 for (int dy = 0; dy < ny; ++dy)
414 {
415 MFEM_UNROLL(MD1)
416 for (int dz = 0; dz < nz; ++dz)
417 {
418 u[dz] += QDD(qx,dy,dz,vd) * By(dy,qy);
419 }
420 }
421 MFEM_UNROLL(MD1)
422 for (int dz = 0; dz < nz; ++dz) { QQD(qx,qy,dz,vd) = u[dz]; }
423 }
424 }
425 }
426 MFEM_SYNC_THREAD;
427 MFEM_FOREACH_THREAD(vd,z,VDIM)
428 {
429 const int nz = (vd == 2) ? D1D : D1D-1;
430 DeviceMatrix Bz = (vd == 2) ? Bc : Bo;
431 MFEM_FOREACH_THREAD(qy,y,Q1D)
432 {
433 MFEM_FOREACH_THREAD(qx,x,Q1D)
434 {
435 real_t u[Q1D];
436 MFEM_UNROLL(MQ1)
437 for (int qz = 0; qz < Q1D; ++qz) { u[qz] = 0.0; }
438 MFEM_UNROLL(MD1)
439 for (int dz = 0; dz < nz; ++dz)
440 {
441 MFEM_UNROLL(MQ1)
442 for (int qz = 0; qz < Q1D; ++qz)
443 {
444 u[qz] += QQD(qx,qy,dz,vd) * Bz(dz,qz);
445 }
446 }
447 MFEM_UNROLL(MQ1)
448 for (int qz = 0; qz < Q1D; ++qz) { QQQ(qx,qy,qz,vd) = u[qz]; }
449 }
450 }
451 }
452 MFEM_SYNC_THREAD;
453 // Apply D operator
454 if (tidz == 0)
455 {
456 MFEM_FOREACH_THREAD(qy,y,Q1D)
457 {
458 MFEM_FOREACH_THREAD(qx,x,Q1D)
459 {
460 MFEM_UNROLL(MQ1)
461 for (int qz = 0; qz < Q1D; ++qz)
462 {
463 const real_t Qx = QQQ(qx,qy,qz,0);
464 const real_t Qy = QQQ(qx,qy,qz,1);
465 const real_t Qz = QQQ(qx,qy,qz,2);
466
467 const real_t D11 = D(qx,qy,qz,0,e);
468 const real_t D12 = D(qx,qy,qz,1,e);
469 const real_t D13 = D(qx,qy,qz,2,e);
470 const real_t D21 = symmetric ? D12 : D(qx,qy,qz,3,e);
471 const real_t D22 = symmetric ? D(qx,qy,qz,3,e) : D(qx,qy,qz,4,e);
472 const real_t D23 = symmetric ? D(qx,qy,qz,4,e) : D(qx,qy,qz,5,e);
473 const real_t D31 = symmetric ? D13 : D(qx,qy,qz,6,e);
474 const real_t D32 = symmetric ? D23 : D(qx,qy,qz,7,e);
475 const real_t D33 = symmetric ? D(qx,qy,qz,5,e) : D(qx,qy,qz,8,e);
476
477 QQQ(qx,qy,qz,0) = D11*Qx + D12*Qy + D13*Qz;
478 QQQ(qx,qy,qz,1) = D21*Qx + D22*Qy + D23*Qz;
479 QQQ(qx,qy,qz,2) = D31*Qx + D32*Qy + D33*Qz;
480 }
481 }
482 }
483 }
484 MFEM_SYNC_THREAD;
485 // Apply Bt operator
486 MFEM_FOREACH_THREAD(vd,z,VDIM)
487 {
488 const int nx = (vd == 0) ? D1D : D1D-1;
489 DeviceMatrix Btx = (vd == 0) ? Bc : Bo;
490 MFEM_FOREACH_THREAD(qy,y,Q1D)
491 {
492 MFEM_FOREACH_THREAD(dx,x,nx)
493 {
494 real_t u[Q1D];
495 MFEM_UNROLL(MQ1)
496 for (int qz = 0; qz < Q1D; ++qz) { u[qz] = 0.0; }
497 MFEM_UNROLL(MQ1)
498 for (int qx = 0; qx < Q1D; ++qx)
499 {
500 MFEM_UNROLL(MQ1)
501 for (int qz = 0; qz < Q1D; ++qz)
502 {
503 u[qz] += QQQ(qx,qy,qz,vd) * Btx(dx,qx);
504 }
505 }
506 MFEM_UNROLL(MQ1)
507 for (int qz = 0; qz < Q1D; ++qz) { DQQ(dx,qy,qz,vd) = u[qz]; }
508 }
509 }
510 }
511 MFEM_SYNC_THREAD;
512 MFEM_FOREACH_THREAD(vd,z,VDIM)
513 {
514 const int nx = (vd == 0) ? D1D : D1D-1;
515 const int ny = (vd == 1) ? D1D : D1D-1;
516 DeviceMatrix Bty = (vd == 1) ? Bc : Bo;
517 MFEM_FOREACH_THREAD(dy,y,ny)
518 {
519 MFEM_FOREACH_THREAD(dx,x,nx)
520 {
521 real_t u[Q1D];
522 MFEM_UNROLL(MQ1)
523 for (int qz = 0; qz < Q1D; ++qz) { u[qz] = 0.0; }
524 MFEM_UNROLL(MQ1)
525 for (int qy = 0; qy < Q1D; ++qy)
526 {
527 MFEM_UNROLL(MQ1)
528 for (int qz = 0; qz < Q1D; ++qz)
529 {
530 u[qz] += DQQ(dx,qy,qz,vd) * Bty(dy,qy);
531 }
532 }
533 MFEM_UNROLL(MQ1)
534 for (int qz = 0; qz < Q1D; ++qz) { DDQ(dx,dy,qz,vd) = u[qz]; }
535 }
536 }
537 }
538 MFEM_SYNC_THREAD;
539 MFEM_FOREACH_THREAD(vd,z,VDIM)
540 {
541 const int nx = (vd == 0) ? D1D : D1D-1;
542 const int ny = (vd == 1) ? D1D : D1D-1;
543 const int nz = (vd == 2) ? D1D : D1D-1;
544 DeviceTensor<5> Yxyz(y, nx, ny, nz, VDIM, NE);
545 DeviceMatrix Btz = (vd == 2) ? Bc : Bo;
546 MFEM_FOREACH_THREAD(dy,y,ny)
547 {
548 MFEM_FOREACH_THREAD(dx,x,nx)
549 {
550 real_t u[D1D];
551 MFEM_UNROLL(MD1)
552 for (int dz = 0; dz < nz; ++dz) { u[dz] = 0.0; }
553 MFEM_UNROLL(MQ1)
554 for (int qz = 0; qz < Q1D; ++qz)
555 {
556 MFEM_UNROLL(MD1)
557 for (int dz = 0; dz < nz; ++dz)
558 {
559 u[dz] += DDQ(dx,dy,qz,vd) * Btz(dz,qz);
560 }
561 }
562 MFEM_UNROLL(MD1)
563 for (int dz = 0; dz < nz; ++dz) { Yxyz(dx,dy,dz,vd,e) += u[dz]; }
564 }
565 }
566 }
567 MFEM_SYNC_THREAD;
568 });
569}
570
571// PA H(div) div-div Assemble 2D kernel
572void PADivDivSetup2D(const int Q1D,
573 const int NE,
574 const Array<real_t> &w,
575 const Vector &j,
576 Vector &coeff_,
577 Vector &op);
578
579// PA H(div) div-div Assemble 3D kernel
580void PADivDivSetup3D(const int Q1D,
581 const int NE,
582 const Array<real_t> &w,
583 const Vector &j,
584 Vector &coeff_,
585 Vector &op);
586
587// PA H(div) div-div Diagonal 2D kernel
588void PADivDivAssembleDiagonal2D(const int D1D,
589 const int Q1D,
590 const int NE,
591 const Array<real_t> &Bo_,
592 const Array<real_t> &Gc_,
593 const Vector &op_,
594 Vector &diag_);
595
596// PA H(div) div-div Diagonal 3D kernel
597void PADivDivAssembleDiagonal3D(const int D1D,
598 const int Q1D,
599 const int NE,
600 const Array<real_t> &Bo_,
601 const Array<real_t> &Gc_,
602 const Vector &op_,
603 Vector &diag_);
604
605// PA H(div) div-div Apply 2D kernel
606void PADivDivApply2D(const int D1D,
607 const int Q1D,
608 const int NE,
609 const Array<real_t> &Bo_,
610 const Array<real_t> &Gc_,
611 const Array<real_t> &Bot_,
612 const Array<real_t> &Gct_,
613 const Vector &op_,
614 const Vector &x_,
615 Vector &y_);
616
617// PA H(div) div-div Apply 3D kernel
618void PADivDivApply3D(const int D1D,
619 const int Q1D,
620 const int NE,
621 const Array<real_t> &Bo_,
622 const Array<real_t> &Gc_,
623 const Array<real_t> &Bot_,
624 const Array<real_t> &Gct_,
625 const Vector &op_,
626 const Vector &x_,
627 Vector &y_);
628
629// PA H(div)-L2 Assemble 2D kernel
630void PAHdivL2Setup2D(const int Q1D,
631 const int NE,
632 const Array<real_t> &w,
633 Vector &coeff_,
634 Vector &op);
635
636// PA H(div)-L2 Assemble 3D kernel
637void PAHdivL2Setup3D(const int Q1D,
638 const int NE,
639 const Array<real_t> &w,
640 Vector &coeff_,
641 Vector &op);
642
643// PA H(div)-L2 Diagonal 2D kernel
644void PAHdivL2AssembleDiagonal_ADAt_2D(const int D1D,
645 const int Q1D,
646 const int L2D1D,
647 const int NE,
648 const Array<real_t> &L2Bo_,
649 const Array<real_t> &Gct_,
650 const Array<real_t> &Bot_,
651 const Vector &op_,
652 const Vector &D_,
653 Vector &diag_);
654
655// PA H(div)-L2 Diagonal 3D kernel
656void PAHdivL2AssembleDiagonal_ADAt_3D(const int D1D,
657 const int Q1D,
658 const int L2D1D,
659 const int NE,
660 const Array<real_t> &L2Bo_,
661 const Array<real_t> &Gct_,
662 const Array<real_t> &Bot_,
663 const Vector &op_,
664 const Vector &D_,
665 Vector &diag_);
666
667// PA H(div)-L2 Apply 2D kernel
668void PAHdivL2Apply2D(const int D1D,
669 const int Q1D,
670 const int L2D1D,
671 const int NE,
672 const Array<real_t> &Bo_,
673 const Array<real_t> &Gc_,
674 const Array<real_t> &L2Bot_,
675 const Vector &op_,
676 const Vector &x_,
677 Vector &y_);
678
679// PA H(div)-L2 Apply Transpose 2D kernel
680void PAHdivL2ApplyTranspose2D(const int D1D,
681 const int Q1D,
682 const int L2D1D,
683 const int NE,
684 const Array<real_t> &L2Bo_,
685 const Array<real_t> &Gct_,
686 const Array<real_t> &Bot_,
687 const Vector &op_,
688 const Vector &x_,
689 Vector &y_);
690
691// PA H(div)-L2 Apply 3D kernel
692void PAHdivL2Apply3D(const int D1D,
693 const int Q1D,
694 const int L2D1D,
695 const int NE,
696 const Array<real_t> &Bo_,
697 const Array<real_t> &Gc_,
698 const Array<real_t> &L2Bot_,
699 const Vector &op_,
700 const Vector &x_,
701 Vector &y_);
702
703// PA H(div)-L2 Apply Transpose 3D kernel
704void PAHdivL2ApplyTranspose3D(const int D1D,
705 const int Q1D,
706 const int L2D1D,
707 const int NE,
708 const Array<real_t> &L2Bo_,
709 const Array<real_t> &Gct_,
710 const Array<real_t> &Bot_,
711 const Vector &op_,
712 const Vector &x_,
713 Vector &y_);
714
715} // namespace internal
716
717} // namespace mfem
718
719#endif
dim
int dim
Definition ex24.cpp:53
mfem::u
real_t u(const Vector &xvec)
Definition lor_mms.hpp:22
mfem::DeviceMatrix
DeviceTensor< 2, real_t > DeviceMatrix
Definition dtensor.hpp:143
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_3D
void forall_3D(int N, int X, int Y, int Z, lambda &&body)
Definition forall.hpp:774
mfem::real_t
float real_t
Definition config.hpp:43
mfem::DeviceCube
DeviceTensor< 3, real_t > DeviceCube
Definition dtensor.hpp:146
Generated on Wed Apr 9 2025 16:53:04 for MFEM by doxygen 1.11.0