xref: /libCEED/tests/t539-operator.h (revision 437c7c9068f3a84dadcfd91d5007d58f42fe2307)
1*437c7c90SJeremy L Thompson // Copyright (c) 2017-2022, Lawrence Livermore National Security, LLC and other CEED contributors.
2*437c7c90SJeremy L Thompson // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
3*437c7c90SJeremy L Thompson //
4*437c7c90SJeremy L Thompson // SPDX-License-Identifier: BSD-2-Clause
5*437c7c90SJeremy L Thompson //
6*437c7c90SJeremy L Thompson // This file is part of CEED:  http://github.com/ceed
7*437c7c90SJeremy L Thompson 
8*437c7c90SJeremy L Thompson #include <ceed.h>
9*437c7c90SJeremy L Thompson 
10*437c7c90SJeremy L Thompson CEED_QFUNCTION(apply)(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
11*437c7c90SJeremy L Thompson   // in[0] is gradient u_0, shape [2, num_comp=2, Q]
12*437c7c90SJeremy L Thompson   // in[1] is mass quadrature data, size (Q)
13*437c7c90SJeremy L Thompson   // in[2] is Poisson quadrature data, size (Q)
14*437c7c90SJeremy L Thompson   // in[3] is u_0, size (Q)
15*437c7c90SJeremy L Thompson   // in[4] is u_1, size (Q)
16*437c7c90SJeremy L Thompson   const CeedScalar *du_0 = in[0], *qd_mass = in[1], *qd_diff = in[2], *u_0 = in[3], *u_1 = in[4];
17*437c7c90SJeremy L Thompson 
18*437c7c90SJeremy L Thompson   // out[0] is output to multiply against v_0, size (Q)
19*437c7c90SJeremy L Thompson   // out[1] is output to multiply against v_1, size (Q)
20*437c7c90SJeremy L Thompson   // out[2] is output to multiply against gradient v_0, shape [2, num_comp=2, Q]
21*437c7c90SJeremy L Thompson   CeedScalar *v_0 = out[0], *v_1 = out[1], *dv_0 = out[2];
22*437c7c90SJeremy L Thompson 
23*437c7c90SJeremy L Thompson   const CeedInt num_comp_0 = 2;
24*437c7c90SJeremy L Thompson 
25*437c7c90SJeremy L Thompson   // Quadrature point loop
26*437c7c90SJeremy L Thompson   for (CeedInt i = 0; i < Q; i++) {
27*437c7c90SJeremy L Thompson     for (CeedInt c = 0; c < num_comp_0; c++) {
28*437c7c90SJeremy L Thompson       // Mass
29*437c7c90SJeremy L Thompson       v_0[i + Q * c] = qd_mass[i] * (c + 1) * u_0[i + Q * c];
30*437c7c90SJeremy L Thompson       // Diff
31*437c7c90SJeremy L Thompson       dv_0[i + Q * (0 * num_comp_0 + c)] =
32*437c7c90SJeremy L Thompson           qd_diff[i + Q * 0] * (c + 1) * du_0[i + Q * (0 * num_comp_0 + c)] + qd_diff[i + Q * 2] * du_0[i + Q * (1 * num_comp_0 + c)];
33*437c7c90SJeremy L Thompson       dv_0[i + Q * (1 * num_comp_0 + c)] =
34*437c7c90SJeremy L Thompson           qd_diff[i + Q * 2] * (c + 1) * du_0[i + Q * (0 * num_comp_0 + c)] + qd_diff[i + Q * 1] * du_0[i + Q * (1 * num_comp_0 + c)];
35*437c7c90SJeremy L Thompson     }
36*437c7c90SJeremy L Thompson     // Mass
37*437c7c90SJeremy L Thompson     v_1[i] = qd_mass[i] * u_1[i];
38*437c7c90SJeremy L Thompson   }
39*437c7c90SJeremy L Thompson 
40*437c7c90SJeremy L Thompson   return 0;
41*437c7c90SJeremy L Thompson }
42