xref: /libCEED/include/ceed/jit-source/gallery/ceed-poisson2dapply.h (revision 94b7b29b41ad8a17add4c577886859ef16f89dec)
1a0154adeSJed Brown // Copyright (c) 2017-2022, Lawrence Livermore National Security, LLC and other CEED contributors.
2a0154adeSJed Brown // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
3a0154adeSJed Brown //
4a0154adeSJed Brown // SPDX-License-Identifier: BSD-2-Clause
5a0154adeSJed Brown //
6a0154adeSJed Brown // This file is part of CEED:  http://github.com/ceed
7a0154adeSJed Brown 
8a0154adeSJed Brown /**
9a0154adeSJed Brown   @brief Ceed QFunction for applying the 2D Poisson operator
10a0154adeSJed Brown **/
11a0154adeSJed Brown 
12*94b7b29bSJeremy L Thompson #ifndef CEED_POISSON2DAPPLY_H
13*94b7b29bSJeremy L Thompson #define CEED_POISSON2DAPPLY_H
14a0154adeSJed Brown 
15c9c2c079SJeremy L Thompson #include <ceed.h>
16c9c2c079SJeremy L Thompson 
172b730f8bSJeremy L Thompson CEED_QFUNCTION(Poisson2DApply)(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
18a0154adeSJed Brown   // in[0] is gradient u, shape [2, nc=1, Q]
19a0154adeSJed Brown   // in[1] is quadrature data, size (3*Q)
202b730f8bSJeremy L Thompson   const CeedScalar(*ug)[CEED_Q_VLA] = (const CeedScalar(*)[CEED_Q_VLA])in[0], (*q_data)[CEED_Q_VLA] = (const CeedScalar(*)[CEED_Q_VLA])in[1];
21a0154adeSJed Brown   // out[0] is output to multiply against gradient v, shape [2, nc=1, Q]
22a0154adeSJed Brown   CeedScalar(*vg)[CEED_Q_VLA] = (CeedScalar(*)[CEED_Q_VLA])out[0];
23a0154adeSJed Brown 
24a0154adeSJed Brown   const CeedInt dim = 2;
25a0154adeSJed Brown 
26a0154adeSJed Brown   // Quadrature point loop
272b730f8bSJeremy L Thompson   CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) {
28a0154adeSJed Brown     // Read qdata (dXdxdXdxT symmetric matrix)
29a0154adeSJed Brown     // Stored in Voigt convention
30a0154adeSJed Brown     // 0 2
31a0154adeSJed Brown     // 2 1
322b730f8bSJeremy L Thompson     const CeedScalar dXdxdXdxT[2][2] = {
332b730f8bSJeremy L Thompson         {q_data[0][i], q_data[2][i]},
342b730f8bSJeremy L Thompson         {q_data[2][i], q_data[1][i]}
35a0154adeSJed Brown     };
36a0154adeSJed Brown 
37a0154adeSJed Brown     // Apply Poisson operator
38a0154adeSJed Brown     // j = direction of vg
392b730f8bSJeremy L Thompson     for (CeedInt j = 0; j < dim; j++) vg[j][i] = (ug[0][i] * dXdxdXdxT[0][j] + ug[1][i] * dXdxdXdxT[1][j]);
40a0154adeSJed Brown   }  // End of Quadrature Point Loop
41a0154adeSJed Brown 
42a0154adeSJed Brown   return CEED_ERROR_SUCCESS;
43a0154adeSJed Brown }
44a0154adeSJed Brown 
45*94b7b29bSJeremy L Thompson #endif  // CEED_POISSON2DAPPLY_H
46