xref: /libCEED/examples/fluids/qfunctions/velocity_gradient_projection.h (revision b6972d7456611f84b0e462eb1490bcb662442e6a)
1 // Copyright (c) 2017-2023, 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 #ifndef velocity_gradient_projection_h
9 #define velocity_gradient_projection_h
10 
11 #include <ceed.h>
12 
13 #include "newtonian_state.h"
14 #include "newtonian_types.h"
15 #include "utils.h"
16 
17 CEED_QFUNCTION_HELPER int VelocityGradientProjectionRHS(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out,
18                                                         StateVariable state_var) {
19   const CeedScalar(*q)[CEED_Q_VLA]         = (const CeedScalar(*)[CEED_Q_VLA])in[0];
20   const CeedScalar(*Grad_q)[5][CEED_Q_VLA] = (const CeedScalar(*)[5][CEED_Q_VLA])in[1];
21   const CeedScalar(*q_data)                = in[2];
22   CeedScalar(*v)[CEED_Q_VLA]               = (CeedScalar(*)[CEED_Q_VLA])out[0];
23 
24   NewtonianIdealGasContext context = (NewtonianIdealGasContext)ctx;
25 
26   CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) {
27     const CeedScalar qi[5] = {q[0][i], q[1][i], q[2][i], q[3][i], q[4][i]};
28     CeedScalar       wdetJ, dXdx[3][3];
29     QdataUnpack_3D(Q, i, q_data, &wdetJ, dXdx);
30 
31     const State s = StateFromQ(context, qi, state_var);
32     State       grad_s[3];
33     for (CeedInt k = 0; k < 3; k++) {
34       CeedScalar dqi[5];
35       for (CeedInt j = 0; j < 5; j++) {
36         dqi[j] = Grad_q[0][j][i] * dXdx[0][k] + Grad_q[1][j][i] * dXdx[1][k] + Grad_q[2][j][i] * dXdx[2][k];
37       }
38       grad_s[k] = StateFromQ_fwd(context, s, dqi, state_var);
39     }
40 
41     CeedScalar grad_velocity[3][3];
42     VelocityGradient(grad_s, grad_velocity);
43 
44     for (CeedInt j = 0; j < 3; j++) {
45       for (CeedInt k = 0; k < 3; k++) {
46         v[j * 3 + k][i] = wdetJ * grad_velocity[j][k];
47       }
48     }
49   }
50   return 0;
51 }
52 
53 CEED_QFUNCTION(VelocityGradientProjectionRHS_Conserv)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
54   return VelocityGradientProjectionRHS(ctx, Q, in, out, STATEVAR_CONSERVATIVE);
55 }
56 
57 CEED_QFUNCTION(VelocityGradientProjectionRHS_Prim)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
58   return VelocityGradientProjectionRHS(ctx, Q, in, out, STATEVAR_PRIMITIVE);
59 }
60 #endif  // velocity_gradient_projection_h
61