Add support for using free boundary conditions in Flow

The OPM spesific keywords FREEBC[XYZ[-]] can be used to specify
boundary cells that are open. Default is a closed boundary.
This commit is contained in:
Tor Harald Sandve 2019-02-06 11:18:51 +01:00
parent 4ecda15b11
commit a5463ed1a0
3 changed files with 95 additions and 21 deletions

View File

@ -146,11 +146,18 @@ public:
if (enableTemperature || enableEnergy)
fluidState.setTemperature(initialState.temperature()[elemIdx]);
// set the phase pressures.
// set the phase pressures, invB factor and density
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
if (!FluidSystem::phaseIsActive(phaseIdx))
continue;
fluidState.setPressure(phaseIdx, initialState.press()[phaseIdx][elemIdx]);
const auto& b = FluidSystem::inverseFormationVolumeFactor(fluidState, phaseIdx, regionIdx);
fluidState.setInvB(phaseIdx, b);
const auto& rho = FluidSystem::density(fluidState, phaseIdx, regionIdx);
fluidState.setDensity(phaseIdx, rho);
}
}
}

View File

@ -361,11 +361,12 @@ protected:
unsigned timeIdx,
const FluidState& exFluidState)
{
bool enableBoundaryMassFlux = false;
const auto& problem = elemCtx.problem();
bool enableBoundaryMassFlux = problem.hasFreeBoundaryConditions();
if (!enableBoundaryMassFlux)
return;
const auto& problem = elemCtx.problem();
const auto& stencil = elemCtx.stencil(timeIdx);
const auto& scvf = stencil.boundaryFace(scvfIdx);

View File

@ -623,6 +623,8 @@ public:
}
tracerModel_.init();
readBoundaryConditions_();
}
void prefetch(const Element& elem) const
@ -823,24 +825,6 @@ public:
aquiferModel_.endTimeStep();
tracerModel_.endTimeStep();
// we no longer need the initial soluiton
if (this->simulator().episodeIndex() == 0 && !initialFluidStates_.empty()) {
// we always need to provide a temperature and if energy is not conserved, we
// use the initial one. This means we have to "salvage" the temperature from
// the initial fluid states before deleting the array.
if (!enableEnergy) {
initialTemperature_.resize(initialFluidStates_.size());
for (unsigned i = 0; i < initialFluidStates_.size(); ++i) {
const auto& fs = initialFluidStates_[i];
initialTemperature_[i] = fs.temperature(/*phaseIdx=*/0);
}
}
initialFluidStates_.clear();
}
updateCompositionChangeLimits_();
}
@ -1308,6 +1292,48 @@ public:
values.setThermalFlow(context, spaceIdx, timeIdx, initialFluidStates_[globalDofIdx]);
}
if (hasFreeBoundaryConditions()) {
unsigned indexInInside = context.intersection(spaceIdx).indexInInside();
unsigned interiorDofIdx = context.interiorScvIndex(spaceIdx, timeIdx);
unsigned globalDofIdx = context.globalSpaceIndex(interiorDofIdx, timeIdx);
switch (indexInInside) {
case 0:
{
if (freebcXMinus_[globalDofIdx])
values.setFreeFlow(context, spaceIdx, timeIdx, initialFluidStates_[globalDofIdx]);
break;
}
case 1:
if (freebcX_[globalDofIdx])
values.setFreeFlow(context, spaceIdx, timeIdx, initialFluidStates_[globalDofIdx]);
break;
case 2:
if (freebcYMinus_[globalDofIdx])
values.setFreeFlow(context, spaceIdx, timeIdx, initialFluidStates_[globalDofIdx]);
break;
case 3:
if (freebcY_[globalDofIdx])
values.setFreeFlow(context, spaceIdx, timeIdx, initialFluidStates_[globalDofIdx]);
break;
case 4:
if (freebcZMinus_[globalDofIdx])
values.setFreeFlow(context, spaceIdx, timeIdx, initialFluidStates_[globalDofIdx]);
break;
case 5:
if (freebcZ_[globalDofIdx])
values.setFreeFlow(context, spaceIdx, timeIdx, initialFluidStates_[globalDofIdx]);
break;
default:
throw std::logic_error("invalid face index for boundary condition");
}
}
}
/*!
@ -1488,6 +1514,9 @@ public:
return (oilVaporizationControl.getType() == Opm::OilVaporizationEnum::VAPPARS);
}
bool hasFreeBoundaryConditions() const
{ return hasFreeBoundaryConditions_; }
private:
bool drsdtActive_() const
{
@ -2264,6 +2293,35 @@ private:
pffDofData_.update(distFn);
}
void readBoundaryConditions_()
{
hasFreeBoundaryConditions_ = false;
readBoundaryConditionKeyword_("FREEBCX", freebcX_ );
readBoundaryConditionKeyword_("FREEBCX-", freebcXMinus_ );
readBoundaryConditionKeyword_("FREEBCY", freebcY_ );
readBoundaryConditionKeyword_("FREEBCY-", freebcYMinus_ );
readBoundaryConditionKeyword_("FREEBCZ", freebcZ_ );
readBoundaryConditionKeyword_("FREEBCZ-", freebcZMinus_ );
}
void readBoundaryConditionKeyword_(const std::string& name, std::vector<bool>& compressedData )
{
const auto& eclProps = this->simulator().vanguard().eclState().get3DProperties();
const auto& vanguard = this->simulator().vanguard();
unsigned numElems = vanguard.gridView().size(/*codim=*/0);
compressedData.resize(numElems, false);
if (eclProps.hasDeckDoubleGridProperty(name)) {
const std::vector<double>& data = eclProps.getDoubleGridProperty(name).getData();
for (unsigned elemIdx = 0; elemIdx < numElems; ++elemIdx) {
unsigned cartElemIdx = vanguard.cartesianIndex(elemIdx);
compressedData[elemIdx] = (data[cartElemIdx] > 0);
}
hasFreeBoundaryConditions_ = true;
}
}
static std::string briefDescription_;
std::vector<Scalar> porosity_;
@ -2309,6 +2367,14 @@ private:
PffGridVector<GridView, Stencil, PffDofData_, DofMapper> pffDofData_;
TracerModel tracerModel_;
bool hasFreeBoundaryConditions_;
std::vector<bool> freebcX_;
std::vector<bool> freebcXMinus_;
std::vector<bool> freebcY_;
std::vector<bool> freebcYMinus_;
std::vector<bool> freebcZ_;
std::vector<bool> freebcZMinus_;
};