1 // SPDX-FileCopyrightText: Copyright (c) 2017-2025, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 #include <ceed/types.h> 4 5 #include "../newtonian_state.h" 6 #include "../numerics.h" 7 #include "../utils.h" 8 9 typedef struct CflPe_SpanStatsContext_ *CflPe_SpanStatsContext; 10 struct CflPe_SpanStatsContext_ { 11 CeedScalar solution_time; 12 CeedScalar previous_time; 13 CeedScalar diffusion_coeff; 14 CeedScalar timestep; 15 struct NewtonianIdealGasContext_ gas; 16 }; 17 18 CEED_QFUNCTION_HELPER int ChildStatsCollection_CflPe(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out, 19 StateVariable state_var, CeedInt dim) { 20 const CeedScalar(*q)[CEED_Q_VLA] = (const CeedScalar(*)[CEED_Q_VLA])in[0]; 21 const CeedScalar(*q_data) = in[1]; 22 CeedScalar(*v)[CEED_Q_VLA] = (CeedScalar(*)[CEED_Q_VLA])out[0]; 23 24 CflPe_SpanStatsContext context = (CflPe_SpanStatsContext)ctx; 25 NewtonianIdealGasContext gas = &context->gas; 26 CeedScalar delta_t = context->solution_time - context->previous_time; 27 28 CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) { 29 const CeedScalar qi[5] = {q[0][i], q[1][i], q[2][i], q[3][i], q[4][i]}; 30 const State s = StateFromQ(gas, qi, state_var); 31 CeedScalar Pe, cfl, wdetJ; 32 33 switch (dim) { 34 case 2: { 35 CeedScalar dXdx[2][2], gijd_mat[2][2] = {{0.}}; 36 37 QdataUnpack_2D(Q, i, q_data, &wdetJ, dXdx); 38 wdetJ = wdetJ * delta_t; 39 40 MatMat2(dXdx, dXdx, CEED_TRANSPOSE, CEED_NOTRANSPOSE, gijd_mat); 41 // (1/2)^2 to account for reference element size; for length 1 square/cube element, gij should be identity matrix 42 ScaleN((CeedScalar *)gijd_mat, 0.25, Square(dim)); 43 44 cfl = CalculateCFL_2D(s.Y.velocity, context->timestep, gijd_mat); 45 Pe = CalculatePe_2D(s.Y.velocity, context->diffusion_coeff, gijd_mat); 46 } break; 47 case 3: { 48 CeedScalar dXdx[3][3], gijd_mat[3][3] = {{0.}}; 49 50 QdataUnpack_3D(Q, i, q_data, &wdetJ, dXdx); 51 wdetJ = wdetJ * delta_t; 52 53 MatMat3(dXdx, dXdx, CEED_TRANSPOSE, CEED_NOTRANSPOSE, gijd_mat); 54 // (1/2)^2 to account for reference element size; for length 1 square/cube element, gij should be identity matrix 55 ScaleN((CeedScalar *)gijd_mat, 0.25, Square(dim)); 56 57 cfl = CalculateCFL_3D(s.Y.velocity, context->timestep, gijd_mat); 58 Pe = CalculatePe_3D(s.Y.velocity, context->diffusion_coeff, gijd_mat); 59 } break; 60 } 61 62 v[0][i] = wdetJ * cfl; 63 v[1][i] = wdetJ * Square(cfl); 64 v[2][i] = wdetJ * Cube(cfl); 65 v[3][i] = wdetJ * Pe; 66 v[4][i] = wdetJ * Square(Pe); 67 v[5][i] = wdetJ * Cube(Pe); 68 } 69 return 0; 70 } 71 72 CEED_QFUNCTION(ChildStatsCollection_3D_Conserv)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 73 return ChildStatsCollection_CflPe(ctx, Q, in, out, STATEVAR_CONSERVATIVE, 3); 74 } 75 76 CEED_QFUNCTION(ChildStatsCollection_3D_Prim)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 77 return ChildStatsCollection_CflPe(ctx, Q, in, out, STATEVAR_PRIMITIVE, 3); 78 } 79 80 CEED_QFUNCTION(ChildStatsCollection_3D_Entropy)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 81 return ChildStatsCollection_CflPe(ctx, Q, in, out, STATEVAR_ENTROPY, 3); 82 } 83 84 CEED_QFUNCTION(ChildStatsCollection_2D_Conserv)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 85 return ChildStatsCollection_CflPe(ctx, Q, in, out, STATEVAR_CONSERVATIVE, 2); 86 } 87 88 CEED_QFUNCTION(ChildStatsCollection_2D_Prim)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 89 return ChildStatsCollection_CflPe(ctx, Q, in, out, STATEVAR_PRIMITIVE, 2); 90 } 91 92 CEED_QFUNCTION(ChildStatsCollection_2D_Entropy)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 93 return ChildStatsCollection_CflPe(ctx, Q, in, out, STATEVAR_ENTROPY, 2); 94 } 95