Skip to content

Commit

Permalink
Merge pull request #1571 from CEED/jrwrigh/createqdata_function
Browse files Browse the repository at this point in the history 
fluids: Add CreateQData and CreateQDataBoundary functions
  • Loading branch information
@jrwrigh
jrwrigh committed Apr 30, 2024
2 parents 7dc3a04 + cfb075a commit 7467d97
Show file tree
Hide file tree
Showing 10 changed files with 318 additions and 145 deletions.
5 changes: 0 additions & 5 deletions examples/fluids/navierstokes.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -285,17 +285,12 @@ int main(int argc, char **argv) {
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_freestream_jacobian.qfunction_context));
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_slip.qfunction_context));
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_slip_jacobian.qfunction_context));
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->setup_sur.qfunction_context));
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->setup_vol.qfunction_context));
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->ics.qfunction_context));
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_vol_rhs.qfunction_context));
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_vol_ifunction.qfunction_context));
PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_vol_ijacobian.qfunction_context));
}

// -- QFunctions
PetscCallCeed(ceed, CeedQFunctionDestroy(&ceed_data->qf_setup_sur));

// -- Operators
PetscCall(OperatorApplyContextDestroy(ceed_data->op_ics_ctx));
PetscCall(OperatorApplyContextDestroy(user->op_rhs_ctx));
Expand Down
11 changes: 8 additions & 3 deletions examples/fluids/navierstokes.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -155,7 +155,6 @@ struct CeedData_private {
CeedBasis basis_x, basis_q;
CeedElemRestriction elem_restr_x, elem_restr_q, elem_restr_qd_i;
OperatorApplyContext op_ics_ctx;
CeedQFunction qf_setup_sur;
};

typedef struct {
Expand Down Expand Up @@ -291,8 +290,8 @@ typedef struct ProblemData_private *ProblemData;
struct ProblemData_private {
CeedInt dim, q_data_size_vol, q_data_size_sur, jac_data_size_sur;
CeedScalar dm_scale;
ProblemQFunctionSpec setup_vol, setup_sur, ics, apply_vol_rhs, apply_vol_ifunction, apply_vol_ijacobian, apply_inflow, apply_outflow,
apply_freestream, apply_slip, apply_inflow_jacobian, apply_outflow_jacobian, apply_freestream_jacobian, apply_slip_jacobian;
ProblemQFunctionSpec ics, apply_vol_rhs, apply_vol_ifunction, apply_vol_ijacobian, apply_inflow, apply_outflow, apply_freestream, apply_slip,
apply_inflow_jacobian, apply_outflow_jacobian, apply_freestream_jacobian, apply_slip_jacobian;
bool non_zero_time;
PetscBool bc_from_ics, use_strong_bc_ceed, uses_newtonian;
PetscErrorCode (*print_info)(User, ProblemData, AppCtx);
Expand Down Expand Up @@ -349,6 +348,12 @@ PetscErrorCode CreateBasisFromPlex(Ceed ceed, DM dm, DMLabel domain_label, CeedI

PetscErrorCode SetupLibceed(Ceed ceed, CeedData ceed_data, DM dm, User user, AppCtx app_ctx, ProblemData problem, SimpleBC bc);

PetscErrorCode QDataGet(Ceed ceed, DM dm, DMLabel domain_label, PetscInt label_value, CeedElemRestriction elem_restr_x, CeedBasis basis_x,
CeedVector x_coord, CeedElemRestriction *elem_restr_qd, CeedVector *q_data, CeedInt *q_data_size);
PetscErrorCode QDataGetNumComponents(DM dm, CeedInt *q_data_size);
PetscErrorCode QDataBoundaryGet(Ceed ceed, DM dm, DMLabel domain_label, PetscInt label_value, CeedElemRestriction elem_restr_x, CeedBasis basis_x,
CeedVector x_coord, CeedElemRestriction *elem_restr_qd, CeedVector *q_data, CeedInt *q_data_size);
PetscErrorCode QDataBoundaryGetNumComponents(DM dm, CeedInt *q_data_size);
// -----------------------------------------------------------------------------
// Time-stepping functions
// -----------------------------------------------------------------------------
Expand Down
14 changes: 0 additions & 14 deletions examples/fluids/problems/advection.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,8 +14,6 @@
#include <petscdm.h>

#include "../navierstokes.h"
#include "../qfunctions/setupgeo.h"
#include "../qfunctions/setupgeo2d.h"

// @brief Create CeedOperator for stabilized mass KSP for explicit timestepping
//
Expand Down Expand Up @@ -106,12 +104,6 @@ PetscErrorCode NS_ADVECTION(ProblemData problem, DM dm, void *ctx, SimpleBC bc)
switch (dim) {
case 2:
problem->dim = 2;
problem->q_data_size_vol = 5;
problem->q_data_size_sur = 3;
problem->setup_vol.qfunction = Setup2d;
problem->setup_vol.qfunction_loc = Setup2d_loc;
problem->setup_sur.qfunction = SetupBoundary2d;
problem->setup_sur.qfunction_loc = SetupBoundary2d_loc;
problem->ics.qfunction = ICsAdvection2d;
problem->ics.qfunction_loc = ICsAdvection2d_loc;
problem->apply_vol_rhs.qfunction = RHS_Advection2d;
Expand All @@ -125,12 +117,6 @@ PetscErrorCode NS_ADVECTION(ProblemData problem, DM dm, void *ctx, SimpleBC bc)
break;
case 3:
problem->dim = 3;
problem->q_data_size_vol = 10;
problem->q_data_size_sur = 10;
problem->setup_vol.qfunction = Setup;
problem->setup_vol.qfunction_loc = Setup_loc;
problem->setup_sur.qfunction = SetupBoundary;
problem->setup_sur.qfunction_loc = SetupBoundary_loc;
problem->ics.qfunction = ICsAdvection;
problem->ics.qfunction_loc = ICsAdvection_loc;
problem->apply_vol_rhs.qfunction = RHS_Advection;
Expand Down
7 changes: 0 additions & 7 deletions examples/fluids/problems/eulervortex.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,6 @@
#include <petscdm.h>

#include "../navierstokes.h"
#include "../qfunctions/setupgeo.h"

PetscErrorCode NS_EULER_VORTEX(ProblemData problem, DM dm, void *ctx, SimpleBC bc) {
EulerTestType euler_test;
Expand All @@ -33,12 +32,6 @@ PetscErrorCode NS_EULER_VORTEX(ProblemData problem, DM dm, void *ctx, SimpleBC b
// SET UP DENSITY_CURRENT
// ------------------------------------------------------
problem->dim = 3;
problem->q_data_size_vol = 10;
problem->q_data_size_sur = 10;
problem->setup_vol.qfunction = Setup;
problem->setup_vol.qfunction_loc = Setup_loc;
problem->setup_sur.qfunction = SetupBoundary;
problem->setup_sur.qfunction_loc = SetupBoundary_loc;
problem->ics.qfunction = ICsEuler;
problem->ics.qfunction_loc = ICsEuler_loc;
problem->apply_vol_rhs.qfunction = Euler;
Expand Down
17 changes: 5 additions & 12 deletions examples/fluids/problems/newtonian.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,6 @@
#include <petscdm.h>

#include "../navierstokes.h"
#include "../qfunctions/setupgeo.h"

// For use with PetscOptionsEnum
static const char *const StateVariables[] = {"CONSERVATIVE", "PRIMITIVE", "StateVariable", "STATEVAR_", NULL};
Expand Down Expand Up @@ -136,17 +135,11 @@ PetscErrorCode NS_NEWTONIAN_IG(ProblemData problem, DM dm, void *ctx, SimpleBC b
// ------------------------------------------------------
// Setup Generic Newtonian IG Problem
// ------------------------------------------------------
problem->dim = 3;
problem->q_data_size_vol = 10;
problem->q_data_size_sur = 10;
problem->jac_data_size_sur = 11;
problem->setup_vol.qfunction = Setup;
problem->setup_vol.qfunction_loc = Setup_loc;
problem->setup_sur.qfunction = SetupBoundary;
problem->setup_sur.qfunction_loc = SetupBoundary_loc;
problem->non_zero_time = PETSC_FALSE;
problem->print_info = PRINT_NEWTONIAN;
problem->uses_newtonian = PETSC_TRUE;
problem->dim = 3;
problem->jac_data_size_sur = 11;
problem->non_zero_time = PETSC_FALSE;
problem->print_info = PRINT_NEWTONIAN;
problem->uses_newtonian = PETSC_TRUE;

// ------------------------------------------------------
// Create the libCEED context
Expand Down
7 changes: 0 additions & 7 deletions examples/fluids/problems/shocktube.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,6 @@
#include <petscdm.h>

#include "../navierstokes.h"
#include "../qfunctions/setupgeo.h"

PetscErrorCode NS_SHOCKTUBE(ProblemData problem, DM dm, void *ctx, SimpleBC bc) {
SetupContextShock setup_context;
Expand All @@ -35,12 +34,6 @@ PetscErrorCode NS_SHOCKTUBE(ProblemData problem, DM dm, void *ctx, SimpleBC bc)
// SET UP SHOCKTUBE
// ------------------------------------------------------
problem->dim = 3;
problem->q_data_size_vol = 10;
problem->q_data_size_sur = 4;
problem->setup_vol.qfunction = Setup;
problem->setup_vol.qfunction_loc = Setup_loc;
problem->setup_sur.qfunction = SetupBoundary;
problem->setup_sur.qfunction_loc = SetupBoundary_loc;
problem->ics.qfunction = ICsShockTube;
problem->ics.qfunction_loc = ICsShockTube_loc;
problem->apply_vol_rhs.qfunction = EulerShockTube;
Expand Down
55 changes: 54 additions & 1 deletion examples/fluids/qfunctions/setupgeo2d.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,6 @@
/// @file
/// Geometric factors (2D) for Navier-Stokes example using PETSc
#include <ceed.h>
#include <math.h>
#include "setupgeo_helpers.h"
#include "utils.h"

Expand Down Expand Up @@ -98,3 +97,57 @@ CEED_QFUNCTION(SetupBoundary2d)(void *ctx, CeedInt Q, const CeedScalar *const *i
}
return 0;
}

// *****************************************************************************
// This QFunction sets up the geometric factor required for integration when reference coordinates are in 2D and the physical coordinates are in 3D
//
// Reference (parent) 2D coordinates: X
// Physical (current) 3D coordinates: x
// Change of coordinate matrix:
// dxdX_{i,j} = dx_i/dX_j (indicial notation) [3 * 2]
// Inverse change of coordinate matrix:
// dXdx_{i,j} = dX_i/dx_j (indicial notation) [2 * 3]
//
// (J1,J2,J3) is given by the cross product of the columns of dxdX_{i,j}
//
// detJb is the magnitude of (J1,J2,J3)
//
// dXdx is calculated via Moore–Penrose inverse:
//
// dX_i/dx_j = (dxdX^T dxdX)^(-1) dxdX
// = (dx_l/dX_i * dx_l/dX_k)^(-1) dx_j/dX_k
//
// All quadrature data is stored in 10 field vector of quadrature data.
//
// We require the determinant of the Jacobian to properly compute integrals of
// the form: int( u v )
//
// Stored: w detJb
// in q_data_sur[0]
//
// Normal vector = (J1,J2,J3) / detJb
//
// Stored: (J1,J2,J3) / detJb
//
// Stored: dXdx_{i,j}
// in q_data_sur[1:6] as
// [dXdx_11 dXdx_12 dXdx_13]
// [dXdx_21 dXdx_22 dXdx_23]
// *****************************************************************************
CEED_QFUNCTION(Setup2D_3Dcoords)(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) {
const CeedScalar(*J)[3][CEED_Q_VLA] = (const CeedScalar(*)[3][CEED_Q_VLA])in[0];
const CeedScalar(*w) = in[1];
CeedScalar(*q_data_sur) = out[0];

CeedPragmaSIMD for (CeedInt i = 0; i < Q; i++) {
CeedScalar detJb, normal[3], dXdx[2][3];

NormalVectorFromdxdX_3D(Q, i, J, normal, &detJb);
InvertBoundaryMappingJacobian_3D(Q, i, J, dXdx);
const CeedScalar wdetJ = w[i] * detJb;

StoredValuesPack(Q, i, 0, 1, &wdetJ, q_data_sur);
StoredValuesPack(Q, i, 1, 6, (const CeedScalar *)dXdx, q_data_sur);
}
return 0;
}

0 comments on commit 7467d97

Please sign in to comment.