1 // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 4 /// @file 5 /// QFunctions for the `bc_freestream` and `bc_outflow` boundary conditions 6 #include "bc_freestream_type.h" 7 #include "newtonian_state.h" 8 #include "newtonian_types.h" 9 #include "riemann_solver.h" 10 11 // ***************************************************************************** 12 // Freestream Boundary Condition 13 // ***************************************************************************** 14 CEED_QFUNCTION_HELPER int Freestream(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out, StateVariable state_var, 15 RiemannFluxType flux_type) { 16 const FreestreamContext context = (FreestreamContext)ctx; 17 const CeedScalar(*q)[CEED_Q_VLA] = (const CeedScalar(*)[CEED_Q_VLA])in[0]; 18 const CeedScalar(*q_data_sur) = in[2]; 19 CeedScalar(*v)[CEED_Q_VLA] = (CeedScalar(*)[CEED_Q_VLA])out[0]; 20 CeedScalar(*jac_data_sur) = context->newtonian_ctx.is_implicit ? out[1] : NULL; 21 22 const NewtonianIdealGasContext newt_ctx = &context->newtonian_ctx; 23 const bool is_implicit = newt_ctx->is_implicit; 24 25 CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) { 26 const CeedScalar qi[5] = {q[0][i], q[1][i], q[2][i], q[3][i], q[4][i]}; 27 const State s = StateFromQ(newt_ctx, qi, state_var); 28 29 CeedScalar wdetJb, normal[3]; 30 QdataBoundaryUnpack_3D(Q, i, q_data_sur, &wdetJb, NULL, normal); 31 wdetJb *= is_implicit ? -1. : 1.; 32 33 StateConservative flux; 34 switch (flux_type) { 35 case RIEMANN_HLL: 36 flux = RiemannFlux_HLL(newt_ctx, s, context->S_infty, normal); 37 break; 38 case RIEMANN_HLLC: 39 flux = RiemannFlux_HLLC(newt_ctx, s, context->S_infty, normal); 40 break; 41 } 42 CeedScalar Flux[5]; 43 UnpackState_U(flux, Flux); 44 for (CeedInt j = 0; j < 5; j++) v[j][i] = -wdetJb * Flux[j]; 45 46 if (is_implicit) StoredValuesPack(Q, i, 0, 5, qi, jac_data_sur); 47 } 48 return 0; 49 } 50 51 CEED_QFUNCTION(Freestream_Conserv_HLL)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 52 return Freestream(ctx, Q, in, out, STATEVAR_CONSERVATIVE, RIEMANN_HLL); 53 } 54 55 CEED_QFUNCTION(Freestream_Prim_HLL)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 56 return Freestream(ctx, Q, in, out, STATEVAR_PRIMITIVE, RIEMANN_HLL); 57 } 58 59 CEED_QFUNCTION(Freestream_Entropy_HLL)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 60 return Freestream(ctx, Q, in, out, STATEVAR_ENTROPY, RIEMANN_HLL); 61 } 62 63 CEED_QFUNCTION(Freestream_Conserv_HLLC)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 64 return Freestream(ctx, Q, in, out, STATEVAR_CONSERVATIVE, RIEMANN_HLLC); 65 } 66 67 CEED_QFUNCTION(Freestream_Prim_HLLC)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 68 return Freestream(ctx, Q, in, out, STATEVAR_PRIMITIVE, RIEMANN_HLLC); 69 } 70 71 CEED_QFUNCTION(Freestream_Entropy_HLLC)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 72 return Freestream(ctx, Q, in, out, STATEVAR_ENTROPY, RIEMANN_HLLC); 73 } 74 75 CEED_QFUNCTION_HELPER int Freestream_Jacobian(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out, StateVariable state_var, 76 RiemannFluxType flux_type) { 77 const CeedScalar(*dq)[CEED_Q_VLA] = (const CeedScalar(*)[CEED_Q_VLA])in[0]; 78 const CeedScalar(*q_data_sur) = in[2]; 79 const CeedScalar(*jac_data_sur) = in[4]; 80 81 CeedScalar(*v)[CEED_Q_VLA] = (CeedScalar(*)[CEED_Q_VLA])out[0]; 82 83 const FreestreamContext context = (FreestreamContext)ctx; 84 const NewtonianIdealGasContext newt_ctx = &context->newtonian_ctx; 85 const bool is_implicit = newt_ctx->is_implicit; 86 const State dS_infty = {0}; 87 88 CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) { 89 CeedScalar wdetJb, normal[3]; 90 QdataBoundaryUnpack_3D(Q, i, q_data_sur, &wdetJb, NULL, normal); 91 wdetJb *= is_implicit ? -1. : 1.; 92 93 CeedScalar qi[5], dqi[5]; 94 StoredValuesUnpack(Q, i, 0, 5, jac_data_sur, qi); 95 for (int j = 0; j < 5; j++) dqi[j] = dq[j][i]; 96 State s = StateFromQ(newt_ctx, qi, state_var); 97 State ds = StateFromQ_fwd(newt_ctx, s, dqi, state_var); 98 99 StateConservative dflux; 100 switch (flux_type) { 101 case RIEMANN_HLL: 102 dflux = RiemannFlux_HLL_fwd(newt_ctx, s, ds, context->S_infty, dS_infty, normal); 103 break; 104 case RIEMANN_HLLC: 105 dflux = RiemannFlux_HLLC_fwd(newt_ctx, s, ds, context->S_infty, dS_infty, normal); 106 break; 107 } 108 CeedScalar dFlux[5]; 109 UnpackState_U(dflux, dFlux); 110 for (CeedInt j = 0; j < 5; j++) v[j][i] = -wdetJb * dFlux[j]; 111 } 112 return 0; 113 } 114 115 CEED_QFUNCTION(Freestream_Jacobian_Conserv_HLL)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 116 return Freestream_Jacobian(ctx, Q, in, out, STATEVAR_CONSERVATIVE, RIEMANN_HLL); 117 } 118 119 CEED_QFUNCTION(Freestream_Jacobian_Prim_HLL)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 120 return Freestream_Jacobian(ctx, Q, in, out, STATEVAR_PRIMITIVE, RIEMANN_HLL); 121 } 122 123 CEED_QFUNCTION(Freestream_Jacobian_Entropy_HLL)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 124 return Freestream_Jacobian(ctx, Q, in, out, STATEVAR_ENTROPY, RIEMANN_HLL); 125 } 126 127 CEED_QFUNCTION(Freestream_Jacobian_Conserv_HLLC)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 128 return Freestream_Jacobian(ctx, Q, in, out, STATEVAR_CONSERVATIVE, RIEMANN_HLLC); 129 } 130 131 CEED_QFUNCTION(Freestream_Jacobian_Prim_HLLC)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 132 return Freestream_Jacobian(ctx, Q, in, out, STATEVAR_PRIMITIVE, RIEMANN_HLLC); 133 } 134 135 CEED_QFUNCTION(Freestream_Jacobian_Entropy_HLLC)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 136 return Freestream_Jacobian(ctx, Q, in, out, STATEVAR_ENTROPY, RIEMANN_HLLC); 137 } 138