1 // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 4 #include <ceed.h> 5 6 #include "utils.h" 7 8 // @brief QFunction to calculate the divergence of the diffusive flux 9 CEED_QFUNCTION_HELPER int ComputeDivDiffusiveFluxGeneric(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out, const CeedInt dim, 10 const CeedInt num_comps) { 11 const CeedScalar *grad_q = in[0]; 12 const CeedScalar(*q_data) = in[1]; 13 CeedScalar(*v)[CEED_Q_VLA] = (CeedScalar(*)[CEED_Q_VLA])out[0]; 14 15 CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) { 16 CeedScalar dXdx[9]; 17 18 QdataUnpack_ND(dim, Q, i, q_data, NULL, dXdx); 19 CeedPragmaSIMD for (CeedInt n = 0; n < num_comps; n++) { 20 CeedScalar grad_qn[9]; 21 22 // Get gradient into dim x dim matrix form, with orientation [flux_direction][gradient_direction] 23 // Equivalent of GradUnpackN 24 const CeedInt offset = Q * n * dim; // offset to reach nth component flux gradients 25 for (CeedInt g = 0; g < dim; g++) { 26 for (CeedInt f = 0; f < dim; f++) { 27 grad_qn[f * dim + g] = grad_q[offset + (Q * num_comps * dim) * g + Q * f + i]; 28 } 29 } 30 v[n][i] = 0; 31 DivergenceND(grad_qn, dXdx, dim, &v[n][i]); 32 } 33 } 34 return 0; 35 } 36 37 CEED_QFUNCTION(ComputeDivDiffusiveFlux3D_4)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) { 38 return ComputeDivDiffusiveFluxGeneric(ctx, Q, in, out, 3, 4); 39 } 40