Skip to content

Commit

Permalink
Add option specifying whether to include BC terms in gradient calc (#198
Browse files Browse the repository at this point in the history
)

Eventually useBCinGrad = true will be the default (and only) option.
But, during development, it is useful to be able to turn it off so we
can ensure we don't introduce any unintended effects by running the
regression tests (which will have to have their 'truth' solutions
regenerated once this change is ready.
  • Loading branch information
trevilo committed May 9, 2023
1 parent 10fa769 commit 1cfdbd0
Show file tree
Hide file tree
Showing 4 changed files with 28 additions and 16 deletions.
4 changes: 3 additions & 1 deletion src/M2ulPhyS.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -600,7 +600,7 @@ void M2ulPhyS::initVariables() {

gradUp_A = new GradNonLinearForm(gradUpfes, intRules, dim, num_equation);
gradUp_A->AddInteriorFaceIntegrator(new GradFaceIntegrator(intRules, dim, num_equation));
gradUp_A->AddBdrFaceIntegrator(new GradFaceIntegrator(intRules, dim, num_equation, bcIntegrator));
gradUp_A->AddBdrFaceIntegrator(new GradFaceIntegrator(intRules, dim, num_equation, bcIntegrator, config.useBCinGrad));

rhsOperator =
new RHSoperator(iter, dim, num_equation, order, eqSystem, max_char_speed, intRules, intRuleType, d_fluxClass,
Expand Down Expand Up @@ -2874,6 +2874,8 @@ void M2ulPhyS::parseReactionInputs() {
}

void M2ulPhyS::parseBCInputs() {
tpsP->getInput("boundaryConditions/useBCinGrad", config.useBCinGrad, false);

// number of BC regions defined
int numWalls, numInlets, numOutlets;
tpsP->getInput("boundaryConditions/numWalls", numWalls, 0);
Expand Down
32 changes: 18 additions & 14 deletions src/faceGradientIntegration.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -34,8 +34,8 @@

// Implementation of class FaceIntegrator
GradFaceIntegrator::GradFaceIntegrator(IntegrationRules *_intRules, const int _dim, const int _num_equation,
BCintegrator *bc)
: dim(_dim), num_equation(_num_equation), intRules(_intRules), bc_(bc) {}
BCintegrator *bc, bool useBCinGrad)
: dim(_dim), num_equation(_num_equation), intRules(_intRules), bc_(bc), useBCinGrad_(useBCinGrad) {}

void GradFaceIntegrator::AssembleFaceVector(const FiniteElement &el1, const FiniteElement &el2,
FaceElementTransformations &Tr, const Vector &elfun, Vector &elvect) {
Expand Down Expand Up @@ -96,18 +96,22 @@ void GradFaceIntegrator::AssembleFaceVector(const FiniteElement &el1, const Fini
if (Tr.Elem2No < 0) {
assert(bc_ != NULL);

const int attr = Tr.Attribute;
std::unordered_map<int, BoundaryCondition *>::const_iterator ibc = bc_->inletBCmap.find(attr);
std::unordered_map<int, BoundaryCondition *>::const_iterator obc = bc_->outletBCmap.find(attr);
std::unordered_map<int, BoundaryCondition *>::const_iterator wbc = bc_->wallBCmap.find(attr);
if (ibc != bc_->inletBCmap.end()) {
ibc->second->computeBdrPrimitiveStateForGradient(iUp1, iUp2);
}
if (obc != bc_->outletBCmap.end()) {
obc->second->computeBdrPrimitiveStateForGradient(iUp1, iUp2);
}
if (wbc != bc_->wallBCmap.end()) {
wbc->second->computeBdrPrimitiveStateForGradient(iUp1, iUp2);
if (useBCinGrad_) {
const int attr = Tr.Attribute;
std::unordered_map<int, BoundaryCondition *>::const_iterator ibc = bc_->inletBCmap.find(attr);
std::unordered_map<int, BoundaryCondition *>::const_iterator obc = bc_->outletBCmap.find(attr);
std::unordered_map<int, BoundaryCondition *>::const_iterator wbc = bc_->wallBCmap.find(attr);
if (ibc != bc_->inletBCmap.end()) {
ibc->second->computeBdrPrimitiveStateForGradient(iUp1, iUp2);
}
if (obc != bc_->outletBCmap.end()) {
obc->second->computeBdrPrimitiveStateForGradient(iUp1, iUp2);
}
if (wbc != bc_->wallBCmap.end()) {
wbc->second->computeBdrPrimitiveStateForGradient(iUp1, iUp2);
}
} else {
iUp2 = iUp1;
}
} else {
elfun2_mat.MultTranspose(shape2, iUp2);
Expand Down
5 changes: 4 additions & 1 deletion src/faceGradientIntegration.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -45,10 +45,13 @@ class GradFaceIntegrator : public NonlinearFormIntegrator {
const int dim;
const int num_equation;
IntegrationRules *intRules;

BCintegrator *bc_; // NB: GradFaceIntegrator is a friend of BCintegrator
const bool useBCinGrad_;

public:
GradFaceIntegrator(IntegrationRules *_intRules, const int _dim, const int _num_equation, BCintegrator *bc = NULL);
GradFaceIntegrator(IntegrationRules *_intRules, const int _dim, const int _num_equation, BCintegrator *bc = NULL,
bool useBCinGrad = false);

virtual void AssembleFaceVector(const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Tr,
const Vector &elfun, Vector &elvect);
Expand Down
3 changes: 3 additions & 0 deletions src/run_configuration.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -183,6 +183,9 @@ class RunConfiguration {
std::vector<WallData> wallBC;
std::vector<pair<int, WallType>> wallPatchType;

// Apply BCs in gradient calculation
bool useBCinGrad;

// Passive scalar data
Array<passiveScalarData*> arrayPassiveScalar;

Expand Down

0 comments on commit 1cfdbd0

Please sign in to comment.