xref: /libCEED/tests/t522-operator.h (revision 4f69910b6e3819988a1446e35e0e85e74672bc23)
1 // Copyright (c) 2017-2025, Lawrence Livermore National Security, LLC and other CEED contributors.
2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
3 //
4 // SPDX-License-Identifier: BSD-2-Clause
5 //
6 // This file is part of CEED:  http://github.com/ceed
7 
8 #include <ceed/types.h>
9 
10 CEED_QFUNCTION(setup)(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
11   const CeedScalar *qw = in[0], *J = in[1];
12   CeedScalar       *qd = out[0];
13 
14   for (CeedInt i = 0; i < Q; i++) {
15     // Qdata stored in Voigt convention
16     // J: 0 2   qd: 0 2   adj(J):  J22 -J12
17     //    1 3       2 1           -J21  J11
18     const CeedScalar J11 = J[i + Q * 0];
19     const CeedScalar J21 = J[i + Q * 1];
20     const CeedScalar J12 = J[i + Q * 2];
21     const CeedScalar J22 = J[i + Q * 3];
22     const CeedScalar w   = qw[i] / (J11 * J22 - J21 * J12);
23     qd[i + Q * 0]        = w * (J12 * J12 + J22 * J22);
24     qd[i + Q * 1]        = w * (J11 * J11 + J21 * J21);
25     qd[i + Q * 2]        = -w * (J11 * J12 + J21 * J22);
26   }  // End of Quadrature Point Loop
27   return 0;
28 }
29 
30 CEED_QFUNCTION(diff)(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
31   const CeedScalar *qd = in[0], *ug = in[1];
32   CeedScalar       *vg = out[0];
33 
34   for (CeedInt i = 0; i < Q; i++) {
35     // Read spatial derivatives of u
36     const CeedScalar du[2] = {ug[i + Q * 0], ug[i + Q * 1]};
37 
38     // Read qdata (dXdxdXdxT symmetric matrix)
39     // Stored in Voigt convention
40     // 0 2
41     // 2 1
42     const CeedScalar dXdxdXdxT[2][2] = {
43         {qd[i + 0 * Q], qd[i + 2 * Q]},
44         {qd[i + 2 * Q], qd[i + 1 * Q]}
45     };
46     // j = direction of vg
47     for (int j = 0; j < 2; j++) vg[i + j * Q] = (du[0] * dXdxdXdxT[0][j] + du[1] * dXdxdXdxT[1][j]);
48   }
49   return 0;
50 }
51