xref: /libCEED/examples/ceed/ex1-volume.h (revision 9ba83ac0e4b1fca39d6fa6737a318a9f0cbc172d)
1*9ba83ac0SJeremy L Thompson // Copyright (c) 2017-2026, Lawrence Livermore National Security, LLC and other CEED contributors.
23d8e8822SJeremy L Thompson // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
366087c08SValeria Barra //
43d8e8822SJeremy L Thompson // SPDX-License-Identifier: BSD-2-Clause
566087c08SValeria Barra //
63d8e8822SJeremy L Thompson // This file is part of CEED:  http://github.com/ceed
766087c08SValeria Barra 
8c0b5abf0SJeremy L Thompson #include <ceed/types.h>
9c9c2c079SJeremy L Thompson 
1066087c08SValeria Barra /// A structure used to pass additional data to f_build_mass
112b730f8bSJeremy L Thompson struct BuildContext {
122b730f8bSJeremy L Thompson   CeedInt dim, space_dim;
132b730f8bSJeremy L Thompson };
1466087c08SValeria Barra 
1566087c08SValeria Barra /// libCEED Q-function for building quadrature data for a mass operator
16d37d859eSJeremy L Thompson CEED_QFUNCTION(build_mass)(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
170a242873SJeremy L Thompson   // in[0] is Jacobians with shape [dim, dim, Q]
180a242873SJeremy L Thompson   // in[1] is quadrature weights with shape [1, Q]
190a242873SJeremy L Thompson   const CeedScalar    *w          = in[1];
20d1d35e2fSjeremylt   CeedScalar          *q_data     = out[0];
210a242873SJeremy L Thompson   struct BuildContext *build_data = (struct BuildContext *)ctx;
2266087c08SValeria Barra 
23d37d859eSJeremy L Thompson   switch (build_data->dim + 10 * build_data->space_dim) {
240a242873SJeremy L Thompson     case 11: {
250a242873SJeremy L Thompson       const CeedScalar(*J)[1][CEED_Q_VLA] = (const CeedScalar(*)[1][CEED_Q_VLA])in[0];
260a242873SJeremy L Thompson 
2766087c08SValeria Barra       // Quadrature Point Loop
280a242873SJeremy L Thompson       CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) { q_data[i] = J[0][0][i] * w[i]; }  // End of Quadrature Point Loop
290a242873SJeremy L Thompson     } break;
300a242873SJeremy L Thompson     case 22: {
310a242873SJeremy L Thompson       const CeedScalar(*J)[2][CEED_Q_VLA] = (const CeedScalar(*)[2][CEED_Q_VLA])in[0];
320a242873SJeremy L Thompson 
3366087c08SValeria Barra       // Quadrature Point Loop
342b730f8bSJeremy L Thompson       CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) {
350a242873SJeremy L Thompson         q_data[i] = (J[0][0][i] * J[1][1][i] - J[0][1][i] * J[1][0][i]) * w[i];
3666087c08SValeria Barra       }  // End of Quadrature Point Loop
370a242873SJeremy L Thompson     } break;
380a242873SJeremy L Thompson     case 33: {
390a242873SJeremy L Thompson       const CeedScalar(*J)[3][CEED_Q_VLA] = (const CeedScalar(*)[3][CEED_Q_VLA])in[0];
400a242873SJeremy L Thompson 
4166087c08SValeria Barra       // Quadrature Point Loop
422b730f8bSJeremy L Thompson       CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) {
430a242873SJeremy L Thompson         q_data[i] =
440a242873SJeremy L Thompson             (J[0][0][i] * (J[1][1][i] * J[2][2][i] - J[1][2][i] * J[2][1][i]) - J[0][1][i] * (J[1][0][i] * J[2][2][i] - J[1][2][i] * J[2][0][i]) +
450a242873SJeremy L Thompson              J[0][2][i] * (J[1][0][i] * J[2][1][i] - J[1][1][i] * J[2][0][i])) *
462b730f8bSJeremy L Thompson             w[i];
4766087c08SValeria Barra       }  // End of Quadrature Point Loop
480a242873SJeremy L Thompson     } break;
4966087c08SValeria Barra   }
500b96b02dSJeremy L Thompson   return CEED_ERROR_SUCCESS;
5166087c08SValeria Barra }
5266087c08SValeria Barra 
5366087c08SValeria Barra /// libCEED Q-function for applying a mass operator
54d37d859eSJeremy L Thompson CEED_QFUNCTION(apply_mass)(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
550a242873SJeremy L Thompson   // in[0], out[0] are solution variables with shape [1, Q]
560a242873SJeremy L Thompson   // in[1] is quadrature data with shape [1, Q]
57d1d35e2fSjeremylt   const CeedScalar *u = in[0], *q_data = in[1];
5866087c08SValeria Barra   CeedScalar       *v = out[0];
5966087c08SValeria Barra 
6066087c08SValeria Barra   // Quadrature Point Loop
612b730f8bSJeremy L Thompson   CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) { v[i] = q_data[i] * u[i]; }  // End of Quadrature Point Loop
620b96b02dSJeremy L Thompson   return CEED_ERROR_SUCCESS;
6366087c08SValeria Barra }
64