1 // Copyright (c) 2017-2026, 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 /// @file
9 /// Mass operator for Navier-Stokes example using PETSc
10 #include <ceed/types.h>
11 #ifndef CEED_RUNNING_JIT_PASS
12 #include <math.h>
13 #endif
14
15 // *****************************************************************************
16 // This QFunction applies the mass matrix to five interlaced fields.
17 //
18 // Inputs:
19 // u - Input vector at quadrature points
20 // q_data - Quadrature weights
21 //
22 // Output:
23 // v - Output vector at quadrature points
24 //
25 // *****************************************************************************
Mass_N(void * ctx,CeedInt Q,const CeedScalar * const * in,CeedScalar * const * out,const CeedInt N)26 CEED_QFUNCTION_HELPER int Mass_N(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out, const CeedInt N) {
27 const CeedScalar(*u)[CEED_Q_VLA] = (const CeedScalar(*)[CEED_Q_VLA])in[0];
28 const CeedScalar(*q_data) = in[1];
29 CeedScalar(*v)[CEED_Q_VLA] = (CeedScalar(*)[CEED_Q_VLA])out[0];
30
31 CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) {
32 CeedPragmaSIMD for (CeedInt j = 0; j < N; j++) { v[j][i] = q_data[i] * u[j][i]; }
33 }
34 return 0;
35 }
36
Mass_1(void * ctx,CeedInt Q,const CeedScalar * const * in,CeedScalar * const * out)37 CEED_QFUNCTION(Mass_1)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { return Mass_N(ctx, Q, in, out, 1); }
38
Mass_5(void * ctx,CeedInt Q,const CeedScalar * const * in,CeedScalar * const * out)39 CEED_QFUNCTION(Mass_5)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { return Mass_N(ctx, Q, in, out, 5); }
40
Mass_7(void * ctx,CeedInt Q,const CeedScalar * const * in,CeedScalar * const * out)41 CEED_QFUNCTION(Mass_7)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { return Mass_N(ctx, Q, in, out, 7); }
42
Mass_9(void * ctx,CeedInt Q,const CeedScalar * const * in,CeedScalar * const * out)43 CEED_QFUNCTION(Mass_9)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { return Mass_N(ctx, Q, in, out, 9); }
44
Mass_22(void * ctx,CeedInt Q,const CeedScalar * const * in,CeedScalar * const * out)45 CEED_QFUNCTION(Mass_22)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { return Mass_N(ctx, Q, in, out, 22); }
46