Merge pull request #5978 from akva2/compositional_container

Add CompositionalContainer, a container for compositional data output
This commit is contained in:
Bård Skaflestad 2025-02-13 12:35:20 +01:00 committed by GitHub
commit 2d9aac3043
No known key found for this signature in database
GPG Key ID: B5690EEEBB952194
7 changed files with 323 additions and 105 deletions

View File

@ -90,6 +90,7 @@ list (APPEND MAIN_SOURCE_FILES
opm/simulators/flow/BlackoilModelParameters.cpp
opm/simulators/flow/BlackoilModelConvergenceMonitor.cpp
opm/simulators/flow/CollectDataOnIORank.cpp
opm/simulators/flow/CompositionalContainer.cpp
opm/simulators/flow/ConvergenceOutputConfiguration.cpp
opm/simulators/flow/EclGenericWriter.cpp
opm/simulators/flow/ExtboContainer.cpp
@ -821,6 +822,7 @@ list (APPEND PUBLIC_HEADER_FILES
opm/simulators/flow/BlackoilModelProperties.hpp
opm/simulators/flow/CollectDataOnIORank.hpp
opm/simulators/flow/CollectDataOnIORank_impl.hpp
opm/simulators/flow/CompositionalContainer.hpp
opm/simulators/flow/ConvergenceOutputConfiguration.hpp
opm/simulators/flow/countGlobalCells.hpp
opm/simulators/flow/CpGridVanguard.hpp

View File

@ -0,0 +1,183 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 2 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
#include <config.h>
#include <opm/simulators/flow/CompositionalContainer.hpp>
#include <opm/material/fluidsystems/GenericOilGasFluidSystem.hpp>
#include <opm/output/data/Solution.hpp>
#include <algorithm>
#include <tuple>
#include <fmt/format.h>
namespace Opm {
template<class FluidSystem>
void CompositionalContainer<FluidSystem>::
allocate(const unsigned bufferSize,
std::map<std::string, int>& rstKeywords)
{
if (auto& zmf = rstKeywords["ZMF"]; zmf > 0) {
this->allocated_ = true;
zmf = 0;
for (int i = 0; i < numComponents; ++i) {
moleFractions_[i].resize(bufferSize, 0.0);
}
}
if (auto& xmf = rstKeywords["XMF"]; xmf > 0 && FluidSystem::phaseIsActive(oilPhaseIdx)) {
this->allocated_ = true;
xmf = 0;
for (int i = 0; i < numComponents; ++i) {
phaseMoleFractions_[oilPhaseIdx][i].resize(bufferSize, 0.0);
}
}
if (auto& ymf = rstKeywords["YMF"]; ymf > 0 && FluidSystem::phaseIsActive(gasPhaseIdx)) {
this->allocated_ = true;
ymf = 0;
for (int i = 0; i < numComponents; ++i) {
phaseMoleFractions_[gasPhaseIdx][i].resize(bufferSize, 0.0);
}
}
}
template<class FluidSystem>
void CompositionalContainer<FluidSystem>::
assignGasFractions(const unsigned globalDofIdx,
const AssignFunction& fractions)
{
if (phaseMoleFractions_[gasPhaseIdx][0].empty()) {
return;
}
std::for_each(phaseMoleFractions_[gasPhaseIdx].begin(),
phaseMoleFractions_[gasPhaseIdx].end(),
[globalDofIdx, &fractions, c = 0](auto& comp) mutable
{ comp[globalDofIdx] = fractions(c++); });
}
template<class FluidSystem>
void CompositionalContainer<FluidSystem>::
assignMoleFractions(const unsigned globalDofIdx,
const AssignFunction& fractions)
{
if (moleFractions_.empty()) {
return;
}
std::for_each(moleFractions_.begin(), moleFractions_.end(),
[&fractions, globalDofIdx, c = 0](auto& comp) mutable
{ comp[globalDofIdx] = fractions(c++); });
}
template<class FluidSystem>
void CompositionalContainer<FluidSystem>::
assignOilFractions(const unsigned globalDofIdx,
const AssignFunction& fractions)
{
if (phaseMoleFractions_[oilPhaseIdx][0].empty()) {
return;
}
std::for_each(phaseMoleFractions_[oilPhaseIdx].begin(),
phaseMoleFractions_[oilPhaseIdx].end(),
[globalDofIdx, &fractions, c = 0](auto& comp) mutable
{ comp[globalDofIdx] = fractions(c++); });
}
template<class FluidSystem>
void CompositionalContainer<FluidSystem>::
outputRestart(data::Solution& sol,
ScalarBuffer& oil_saturation)
{
using DataEntry =
std::tuple<std::string, UnitSystem::measure, std::vector<Scalar>&>;
auto doInsert = [&sol](DataEntry& entry,
const data::TargetType target)
{
if (std::get<2>(entry).empty()) {
return;
}
sol.insert(std::get<std::string>(entry),
std::get<UnitSystem::measure>(entry),
std::move(std::get<2>(entry)),
target);
};
auto entries = std::vector<DataEntry>{};
// ZMF
if (!moleFractions_[0].empty()) {
for (int i = 0; i < numComponents; ++i) {
const auto name = fmt::format("ZMF{}", i + 1); // Generate ZMF1, ZMF2, ...
entries.emplace_back(name, UnitSystem::measure::identity, moleFractions_[i]);
}
}
// XMF
if (!phaseMoleFractions_[oilPhaseIdx][0].empty()) {
for (int i = 0; i < numComponents; ++i) {
const auto name = fmt::format("XMF{}", i + 1); // Generate XMF1, XMF2, ...
entries.emplace_back(name, UnitSystem::measure::identity,
phaseMoleFractions_[oilPhaseIdx][i]);
}
}
// YMF
if (!phaseMoleFractions_[gasPhaseIdx][0].empty()) {
for (int i = 0; i < numComponents; ++i) {
const auto name = fmt::format("YMF{}", i + 1); // Generate YMF1, YMF2, ...
entries.emplace_back(name, UnitSystem::measure::identity,
phaseMoleFractions_[gasPhaseIdx][i]);
}
}
if (!oil_saturation.empty()) {
entries.emplace_back("SOIL", UnitSystem::measure::identity, oil_saturation);
}
std::for_each(entries.begin(), entries.end(),
[&doInsert](auto& array)
{ doInsert(array, data::TargetType::RESTART_SOLUTION); });
this->allocated_ = false;
}
#define INSTANTIATE_COMP(NUM) \
template<class T> using FS##NUM = GenericOilGasFluidSystem<T, NUM>; \
template class CompositionalContainer<FS##NUM<double>>;
INSTANTIATE_COMP(0)
INSTANTIATE_COMP(2)
INSTANTIATE_COMP(3)
INSTANTIATE_COMP(4)
INSTANTIATE_COMP(5)
INSTANTIATE_COMP(6)
INSTANTIATE_COMP(7)
} // namespace Opm

View File

@ -0,0 +1,83 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 2 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
/*!
* \file
* \copydoc Opm::OutputBlackOilModule
*/
#ifndef OPM_COMPOSITIONAL_CONTAINER_HPP
#define OPM_COMPOSITIONAL_CONTAINER_HPP
#include <array>
#include <functional>
#include <map>
#include <string>
#include <vector>
namespace Opm {
namespace data { class Solution; }
template<class FluidSystem>
class CompositionalContainer
{
using Scalar = typename FluidSystem::Scalar;
using ScalarBuffer = std::vector<Scalar>;
static constexpr int numComponents = FluidSystem::numComponents;
static constexpr int numPhases = FluidSystem::numPhases;
static constexpr int gasPhaseIdx = FluidSystem::gasPhaseIdx;
static constexpr int oilPhaseIdx = FluidSystem::oilPhaseIdx;
static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
public:
void allocate(const unsigned bufferSize,
std::map<std::string, int>& rstKeywords);
using AssignFunction = std::function<Scalar(const unsigned)>;
void assignGasFractions(const unsigned globalDofIdx,
const AssignFunction& fractions);
void assignMoleFractions(const unsigned globalDofIdx,
const AssignFunction& fractions);
void assignOilFractions(const unsigned globalDofIdx,
const AssignFunction& fractions);
void outputRestart(data::Solution& sol,
ScalarBuffer& oil_saturation);
bool allocated() const
{ return allocated_; }
private:
bool allocated_ = false;
// total mole fractions for each component
std::array<ScalarBuffer, numComponents> moleFractions_;
// mole fractions for each component in each phase
std::array<std::array<ScalarBuffer, numComponents>, numPhases> phaseMoleFractions_;
};
} // namespace Opm
#endif // OPM_COMPOSITIONAL_CONTAINER_HPP

View File

@ -104,8 +104,7 @@ GenericOutputBlackoilModule(const EclipseState& eclState,
bool enableBrine,
bool enableSaltPrecipitation,
bool enableExtbo,
bool enableMICP,
bool isCompositional)
bool enableMICP)
: eclState_(eclState)
, schedule_(schedule)
, summaryState_(summaryState)
@ -124,7 +123,6 @@ GenericOutputBlackoilModule(const EclipseState& eclState,
, enableSaltPrecipitation_(enableSaltPrecipitation)
, enableExtbo_(enableExtbo)
, enableMICP_(enableMICP)
, isCompositional_(isCompositional)
, local_data_valid_(false)
{
const auto& fp = eclState_.fieldProps();
@ -527,38 +525,6 @@ assignToSolution(data::Solution& sol)
DataEntry{"TMULT_RC", UnitSystem::measure::identity, rockCompTransMultiplier_},
};
// basically, for compositional, we can not use std::array for this. We need to generate the ZMF1, ZMF2, and so on
// and also, we need to map these values.
// TODO: the following should go to a function
if (this->isCompositional_) {
auto compositionalEntries = std::vector<DataEntry>{};
{
// ZMF
for (int i = 0; i < numComponents; ++i) {
const auto name = fmt::format("ZMF{}", i + 1); // Generate ZMF1, ZMF2, ...
compositionalEntries.emplace_back(name, UnitSystem::measure::identity, moleFractions_[i]);
}
// XMF
for (int i = 0; i < numComponents; ++i) {
const auto name = fmt::format("XMF{}", i + 1); // Generate XMF1, XMF2, ...
compositionalEntries.emplace_back(name, UnitSystem::measure::identity,
phaseMoleFractions_[oilPhaseIdx][i]);
}
// YMF
for (int i = 0; i < numComponents; ++i) {
const auto name = fmt::format("YMF{}", i + 1); // Generate YMF1, YMF2, ...
compositionalEntries.emplace_back(name, UnitSystem::measure::identity,
phaseMoleFractions_[gasPhaseIdx][i]);
}
}
for (auto& array: compositionalEntries) {
doInsert(array, data::TargetType::RESTART_SOLUTION);
}
}
for (auto& array : baseSolutionVector) {
doInsert(array, data::TargetType::RESTART_SOLUTION);
}
@ -602,15 +568,6 @@ assignToSolution(data::Solution& sol)
data::TargetType::RESTART_SOLUTION);
}
if (this->isCompositional_ && FluidSystem::phaseIsActive(oilPhaseIdx) &&
! this->saturation_[oilPhaseIdx].empty())
{
sol.insert("SOIL", UnitSystem::measure::identity,
std::move(this->saturation_[oilPhaseIdx]),
data::TargetType::RESTART_SOLUTION);
}
if ((eclState_.runspec().co2Storage() || eclState_.runspec().h2Storage()) && !rsw_.empty()) {
auto mfrac = std::vector<double>(this->rsw_.size(), 0.0);
@ -834,10 +791,14 @@ doAllocBuffers(const unsigned bufferSize,
const bool enableWettingHysteresis,
const unsigned numTracers,
const std::vector<bool>& enableSolTracers,
const unsigned numOutputNnc)
const unsigned numOutputNnc,
std::map<std::string, int> rstKeywords)
{
if (rstKeywords.empty()) {
rstKeywords = schedule_.rst_keywords(reportStepNum);
}
// Output RESTART_OPM_EXTENDED only when explicitly requested by user.
std::map<std::string, int> rstKeywords = schedule_.rst_keywords(reportStepNum);
for (auto& [keyword, should_write] : rstKeywords) {
if (this->isOutputCreationDirective_(keyword)) {
// 'BASIC', 'FREQ' and similar. Don't attempt to create
@ -1293,30 +1254,6 @@ doAllocBuffers(const unsigned bufferSize,
overburdenPressure_.resize(bufferSize, 0.0);
}
if (this->isCompositional_) {
if (rstKeywords["ZMF"] > 0) {
rstKeywords["ZMF"] = 0;
for (int i = 0; i < numComponents; ++i) {
moleFractions_[i].resize(bufferSize, 0.0);
}
}
if (rstKeywords["XMF"] > 0 && FluidSystem::phaseIsActive(oilPhaseIdx)) {
rstKeywords["XMF"] = 0;
for (int i = 0; i < numComponents; ++i) {
phaseMoleFractions_[oilPhaseIdx][i].resize(bufferSize, 0.0);
}
}
if (rstKeywords["YMF"] > 0 && FluidSystem::phaseIsActive(gasPhaseIdx)) {
rstKeywords["YMF"] = 0;
for (int i = 0; i < numComponents; ++i) {
phaseMoleFractions_[gasPhaseIdx][i].resize(bufferSize, 0.0);
}
}
}
//Warn for any unhandled keyword
if (log) {
for (auto& keyValue: rstKeywords) {

View File

@ -319,8 +319,7 @@ protected:
bool enableBrine,
bool enableSaltPrecipitation,
bool enableExtbo,
bool enableMICP,
bool isCompositional = false);
bool enableMICP);
void doAllocBuffers(unsigned bufferSize,
unsigned reportStepNum,
@ -329,11 +328,12 @@ protected:
const bool isRestart,
const bool vapparsActive = false,
const bool enablePCHysteresis = false,
const bool enableNonWettingHysteresis =false,
const bool enableNonWettingHysteresis = false,
const bool enableWettingHysteresis = false,
unsigned numTracers = 0,
const std::vector<bool>& enableSolTracers = {},
unsigned numOutputNnc = 0);
unsigned numOutputNnc = 0,
std::map<std::string, int> rstKeywords = {});
void makeRegionSum(Inplace& inplace,
const std::string& region_name,
@ -390,7 +390,6 @@ protected:
bool enableSaltPrecipitation_{false};
bool enableExtbo_{false};
bool enableMICP_{false};
bool isCompositional_{false};
bool forceDisableFipOutput_{false};
bool forceDisableFipresvOutput_{false};
@ -468,10 +467,6 @@ protected:
std::array<ScalarBuffer, numPhases> viscosity_;
std::array<ScalarBuffer, numPhases> relativePermeability_;
// total mole fractions for each component
std::array<ScalarBuffer, numComponents> moleFractions_;
// mole fractions for each component in each phase
std::array<std::array<ScalarBuffer, numComponents>, numPhases> phaseMoleFractions_;
std::vector<ScalarBuffer> freeTracerConcentrations_;
std::vector<ScalarBuffer> solTracerConcentrations_;

View File

@ -32,15 +32,20 @@
#include <opm/simulators/utils/moduleVersion.hpp>
#include <opm/common/Exceptions.hpp>
#include <opm/common/ErrorMacros.hpp>
#include <opm/common/TimingMacros.hpp>
#include <opm/common/OpmLog/OpmLog.hpp>
#include <opm/input/eclipse/EclipseState/SummaryConfig/SummaryConfig.hpp>
#include <opm/material/common/Valgrind.hpp>
#include <opm/models/blackoil/blackoilproperties.hh>
#include <opm/models/common/multiphasebaseproperties.hh>
#include <opm/models/utils/parametersystem.hpp>
#include <opm/models/utils/propertysystem.hh>
#include <opm/simulators/flow/CompositionalContainer.hpp>
#include <opm/simulators/flow/FlowBaseVanguard.hpp>
#include <opm/simulators/flow/GenericOutputBlackoilModule.hpp>
@ -53,8 +58,7 @@
#include <vector>
namespace Opm
{
namespace Opm {
// forward declaration
template <class TypeTag>
@ -102,8 +106,7 @@ public:
getPropValue<TypeTag, Properties::EnableBrine>(),
getPropValue<TypeTag, Properties::EnableSaltPrecipitation>(),
getPropValue<TypeTag, Properties::EnableExtbo>(),
getPropValue<TypeTag, Properties::EnableMICP>(),
true)
getPropValue<TypeTag, Properties::EnableMICP>())
, simulator_(simulator)
{
for (auto& region_pair : this->regions_) {
@ -155,11 +158,24 @@ public:
return;
}
this->doAllocBuffers(bufferSize,
reportStepNum,
substep,
log,
isRestart);
auto rstKeywords = this->schedule_.rst_keywords(reportStepNum);
this->compC_.allocate(bufferSize, rstKeywords);
this->doAllocBuffers(bufferSize, reportStepNum, substep, log, isRestart,
/* vapparsActive =*/ false,
/* enablePCHysteresis = */ false,
/* enableNonWettingHysteresis =*/ false,
/* enableWettingHysteresis =*/ false,
/* numTracers = */ 0,
/* enableSoltracers =*/ {},
/* numOutputNnc =*/ 0,
std::move(rstKeywords));
}
void assignToSolution(data::Solution& sol)
{
this->compC_.outputRestart(sol, this->saturation_[oilPhaseIdx]);
BaseType::assignToSolution(sol);
}
/*!
@ -187,20 +203,21 @@ public:
Valgrind::CheckDefined(this->saturation_[phaseIdx][globalDofIdx]);
}
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
if (this->moleFractions_[compIdx].empty()) continue;
if (this->compC_.allocated()) {
this->compC_.assignMoleFractions(globalDofIdx,
[&fs](const unsigned compIdx)
{ return getValue(fs.moleFraction(compIdx)); });
this->moleFractions_[compIdx][globalDofIdx] = getValue(fs.moleFraction(compIdx));
}
// XMF and YMF
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
if (this->phaseMoleFractions_[oilPhaseIdx][compIdx].empty()) continue;
this->phaseMoleFractions_[oilPhaseIdx][compIdx][globalDofIdx] = getValue(fs.moleFraction(oilPhaseIdx, compIdx));
}
if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
if (this->phaseMoleFractions_[gasPhaseIdx][compIdx].empty()) continue;
this->phaseMoleFractions_[gasPhaseIdx][compIdx][globalDofIdx] = getValue(fs.moleFraction(gasPhaseIdx, compIdx));
this->compC_.assignGasFractions(globalDofIdx,
[&fs](const unsigned compIdx)
{ return getValue(fs.moleFraction(gasPhaseIdx, compIdx)); });
}
if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
this->compC_.assignOilFractions(globalDofIdx,
[&fs](const unsigned compIdx)
{ return getValue(fs.moleFraction(oilPhaseIdx, compIdx)); });
}
}
@ -332,6 +349,7 @@ private:
}
const Simulator& simulator_;
CompositionalContainer<FluidSystem> compC_;
};
} // namespace Opm

View File

@ -118,13 +118,13 @@ template<class Scalar> class WellContributions;
// TODO: where we should put these types, WellInterface or Well Model?
// or there is some other strategy, like TypeTag
typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
typedef Dune::BlockVector<VectorBlockType> BVector;
using VectorBlockType = Dune::FieldVector<Scalar, numEq>;
using BVector = Dune::BlockVector<VectorBlockType>;
typedef BlackOilPolymerModule<TypeTag> PolymerModule;
typedef BlackOilMICPModule<TypeTag> MICPModule;
using PolymerModule = BlackOilPolymerModule<TypeTag>;
using MICPModule = BlackOilMICPModule<TypeTag>;
// For the conversion between the surface volume rate and resrevoir voidage rate
// For the conversion between the surface volume rate and reservoir voidage rate
using RateConverterType = RateConverter::
SurfaceToReservoirVoidage<FluidSystem, std::vector<int> >;