support source term from deck (SOURCE)

This commit is contained in:
Tor Harald Sandve
2023-12-08 16:49:48 +01:00
parent bf55fcfd93
commit a57d04fde5

View File

@@ -1407,6 +1407,34 @@ public:
if (enableAquifers_)
aquiferModel_.addToSource(rate, globalDofIdx, timeIdx);
// Add source term from deck
const auto& sourceprop = this->simulator().vanguard().schedule()[this->episodeIndex()].sourceprop();
if (sourceprop.size() > 0) {
std::array<int,3> ijk;
this->simulator().vanguard().cartesianCoordinate(globalDofIdx, ijk);
RateVector massRate(0.0);
if ( FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
massRate[Indices::canonicalToActiveComponentIndex(oilCompIdx)] = sourceprop.rate({ijk, SourceComponent::OIL}) / this->model().dofTotalVolume(globalDofIdx);
}
if ( FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
massRate[Indices::canonicalToActiveComponentIndex(gasCompIdx)] = sourceprop.rate({ijk, SourceComponent::GAS}) / this->model().dofTotalVolume(globalDofIdx);
}
if ( FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
massRate[Indices::canonicalToActiveComponentIndex(waterCompIdx)] = sourceprop.rate({ijk, SourceComponent::WATER}) / this->model().dofTotalVolume(globalDofIdx);
}
if constexpr (enableSolvent) {
massRate[Indices::solventSaturationIdx] = sourceprop.rate({ijk, SourceComponent::SOLVENT}) / this->model().dofTotalVolume(globalDofIdx);
}
const int pvtRegionIdx = this->pvtRegionIndex(globalDofIdx);
rate.setMassRate(massRate, pvtRegionIdx);
if constexpr (enablePolymer) {
rate[Indices::polymerConcentrationIdx] = sourceprop.rate({ijk, SourceComponent::POLYMER}) / this->model().dofTotalVolume(globalDofIdx);
}
if constexpr (enableEnergy) {
rate[Indices::contiEnergyEqIdx] = sourceprop.hrate(ijk) / this->model().dofTotalVolume(globalDofIdx);
}
}
// if requested, compensate systematic mass loss for cells which were "well
// behaved" in the last time step
// Note that we don't allow for drift compensation if there are no active wells.