OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

output XMF ,YMF and ZMF in UNRST file.

Browse Source 
currently, it looks the value for ZMF is based on XMF and YMF.
ideally, the ZMF is part of the primary varialbes, it does not need to
be calculated based on XMF and YMF.

there is some refactoring that we can do here.
This commit is contained in:
Kai Bao 2024-10-14 14:57:46 +02:00
parent 4ba1cd53d7
commit ad587e77ad
7 changed files with 116 additions and 36 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
opm/models/ptflash/flashintensivequantities.hh
View File

@ -113,6 +113,8 @@ public:
const Scalar flashTolerance = Parameters::Get<Parameters::FlashTolerance<Scalar>>();
const int flashVerbosity = Parameters::Get<Parameters::FlashVerbosity>();
const std::string flashTwoPhaseMethod = Parameters::Get<Parameters::FlashTwoPhaseMethod>();
// TODO: the formulation here is still to begin with XMF and YMF values to derive ZMF value
// TODO: we should check how we update ZMF in the newton update, since it is the primary variables.
// extract the total molar densities of the components
ComponentVector z(0.);
@ -132,6 +134,10 @@ public:
z /= sumz;
}
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
fluidState_.setMoleFraction(compIdx, z[compIdx]);
}
Evaluation p = priVars.makeEvaluation(pressure0Idx, timeIdx);
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
fluidState_.setPressure(phaseIdx, p);

16
opm/models/ptflash/flashprimaryvariables.hh
View File

@ -120,22 +120,10 @@ public:
// the energy module
EnergyModule::setPriVarTemperatures(*this, fluidState);
// determine the component fractions
Dune::FieldVector<Scalar, numComponents> z(0.0);
Scalar sumMoles = 0.0;
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
Scalar tmp = Opm::getValue(fluidState.molarity(phaseIdx, compIdx) * fluidState.saturation(phaseIdx));
z[compIdx] += Opm::max(tmp, 1e-8);
sumMoles += tmp;
}
}
z /= sumMoles;
for (int i = 0; i < numComponents - 1; ++i)
(*this)[z0Idx + i] = z[i];
(*this)[z0Idx + i] = getValue(fluidState.moleFraction(i));
(*this)[pressure0Idx] = Opm::getValue(fluidState.pressure(0));
(*this)[pressure0Idx] = getValue(fluidState.pressure(0));
}
/*!

43
opm/simulators/flow/FlowProblemComp.hpp
View File

@ -326,25 +326,26 @@ public:
{
const unsigned globalDofIdx = context.globalSpaceIndex(spaceIdx, timeIdx);
const auto& initial_fs = initialFluidStates_[globalDofIdx];
Opm::CompositionalFluidState<Evaluation, FluidSystem> fs;
using ComponentVector = Dune::FieldVector<Evaluation, numComponents>;
Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
// TODO: the current approach is assuming we begin with XMF and YMF.
// TODO: maybe we should make it begin with ZMF
using ComponentVector = Dune::FieldVector<Scalar, numComponents>;
for (unsigned p = 0; p < numPhases; ++p) { // TODO: assuming the phaseidx continuous
ComponentVector evals;
auto& last_eval = evals[numComponents - 1];
ComponentVector vals;
auto& last_eval = vals[numComponents - 1];
last_eval = 1.;
for (unsigned c = 0; c < numComponents - 1; ++c) {
const auto val = initial_fs.moleFraction(p, c);
const Evaluation eval = Evaluation::createVariable(val, c + 1);
evals[c] = eval;
last_eval -= eval;
vals[c] = val;
last_eval -= val;
}
for (unsigned c = 0; c < numComponents; ++c) {
fs.setMoleFraction(p, c, evals[c]);
fs.setMoleFraction(p, c, vals[c]);
}
// pressure
const auto p_val = initial_fs.pressure(p);
fs.setPressure(p, Evaluation::createVariable(p_val, 0));
fs.setPressure(p, p_val);
const auto sat_val = initial_fs.saturation(p);
fs.setSaturation(p, sat_val);
@ -354,7 +355,7 @@ public:
}
{
typename FluidSystem::template ParameterCache<Evaluation> paramCache;
typename FluidSystem::template ParameterCache<Scalar> paramCache;
paramCache.updatePhase(fs, FluidSystem::oilPhaseIdx);
paramCache.updatePhase(fs, FluidSystem::gasPhaseIdx);
fs.setDensity(FluidSystem::oilPhaseIdx, FluidSystem::density(fs, paramCache, FluidSystem::oilPhaseIdx));
@ -363,13 +364,29 @@ public:
fs.setViscosity(FluidSystem::gasPhaseIdx, FluidSystem::viscosity(fs, paramCache, FluidSystem::gasPhaseIdx));
}
// determine the component fractions
Dune::FieldVector<Scalar, numComponents> z(0.0);
Scalar sumMoles = 0.0;
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
Scalar tmp = Opm::getValue(fs.molarity(phaseIdx, compIdx) * fs.saturation(phaseIdx));
z[compIdx] += Opm::max(tmp, 1e-8);
sumMoles += tmp;
}
}
z /= sumMoles;
// Set initial K and L
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
const Evaluation Ktmp = fs.wilsonK_(compIdx);
const auto& Ktmp = fs.wilsonK_(compIdx);
fs.setKvalue(compIdx, Ktmp);
}
const Evaluation& Ltmp = -1.0;
for (unsigned compIdx = 0; compIdx < numComponents - 1; ++compIdx) {
fs.setMoleFraction(compIdx, z[compIdx]);
}
const Scalar& Ltmp = -1.0;
fs.setLvalue(Ltmp);
values.assignNaive(fs);
@ -540,6 +557,8 @@ protected:
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
const std::size_t data_idx = compIdx * numDof + dofIdx;
const Scalar zmf = zmfData[data_idx];
dofFluidState.setMoleFraction(compIdx, zmf);
if (gas_active) {
const auto ymf = (dofFluidState.saturation(FluidSystem::gasPhaseIdx) > 0.) ? zmf : Scalar{0};
dofFluidState.setMoleFraction(FluidSystem::gasPhaseIdx, compIdx, ymf);

49
opm/simulators/flow/GenericOutputCompositionalModule.cpp
View File

@ -25,8 +25,8 @@
#include <opm/common/OpmLog/OpmLog.hpp>
#include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
#include <opm/material/fluidsystems/BlackOilDefaultIndexTraits.hpp>
// #include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
// #include <opm/material/fluidsystems/BlackOilDefaultIndexTraits.hpp>
#include <opm/grid/common/CommunicationUtils.hpp>
@ -590,6 +590,32 @@ assignToSolution(data::Solution& sol)
DataEntry{"UREA", UnitSystem::measure::density, cUrea_},
};
auto compositionalEntries = std::vector<DataEntry>{};
{
// ZMF
for (int i = 0; i < numComponents; ++i) {
const std::string name = "ZMF" + std::to_string(i + 1); // Generate ZMF1, ZMF2, ...
compositionalEntries.emplace_back(name, UnitSystem::measure::identity, moleFractions_[i]);
}
// XMF
for (int i = 0; i < numComponents; ++i) {
const std::string name = "XMF" + std::to_string(i + 1); // Generate XMF1, XMF2, ...
compositionalEntries.emplace_back(name, UnitSystem::measure::identity, phaseMoleFractions_[oilPhaseIdx][i]);
}
// YMF
for (int i = 0; i < numComponents; ++i) {
const std::string name = "YMF" + std::to_string(i + 1); // Generate YMF1, YMF2, ...
compositionalEntries.emplace_back(name, UnitSystem::measure::identity, phaseMoleFractions_[gasPhaseIdx][i]);
}
}
for (const auto& array : compositionalEntries) {
doInsert(array, data::TargetType::RESTART_SOLUTION);
}
// 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.
@ -1459,8 +1485,27 @@ doAllocBuffers(const unsigned bufferSize,
}
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) {

1
opm/simulators/flow/GenericOutputCompositionalModule.hpp
View File

@ -536,6 +536,7 @@ protected:
std::array<ScalarBuffer, numPhases> relativePermeability_;
std::array<ScalarBuffer, numComponents> moleFractions_;
std::array<std::array<ScalarBuffer, numComponents>, numPhases> phaseMoleFractions_;
std::vector<ScalarBuffer> freeTracerConcentrations_;
std::vector<ScalarBuffer> solTracerConcentrations_;

13
opm/simulators/flow/OutputBlackoilModule.hpp
View File

@ -65,13 +65,16 @@
#include <utility>
#include <vector>
//namespace Opm::Parameters {
// TODO: EclWriter how includes both OutputBlackoilModule and OutputCompositionalModule, so the following are defined twice
// Will see whether some refactoring are needed here.
// namespace Opm::Parameters {
//
//struct ForceDisableFluidInPlaceOutput { static constexpr bool value = false; };
//struct ForceDisableResvFluidInPlaceOutput { static constexpr bool value = false; };
//
//} // namespace Opm::Parameters
// struct ForceDisableFluidInPlaceOutput { static constexpr bool value = false; };
// struct ForceDisableResvFluidInPlaceOutput { static constexpr bool value = false; };
//
// } // namespace Opm::Parameters
namespace Opm
{

24
opm/simulators/flow/OutputCompositionalModule.hpp
View File

@ -39,8 +39,8 @@
#include <opm/material/common/Valgrind.hpp>
#include <opm/material/fluidmatrixinteractions/EclEpsScalingPoints.hpp>
#include <opm/material/fluidstates/BlackOilFluidState.hpp>
#include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
// #include <opm/material/fluidstates/BlackOilFluidState.hpp>
// #include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
#include <opm/models/blackoil/blackoilproperties.hh>
#include <opm/models/discretization/common/fvbaseproperties.hh>
@ -106,6 +106,7 @@ class OutputCompositionalModule : public GenericOutputCompositionalModule<GetPro
enum { conti0EqIdx = Indices::conti0EqIdx };
enum { numPhases = FluidSystem::numPhases };
enum { numComponents = FluidSystem::numComponents };
enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
@ -144,7 +145,7 @@ public:
return collectToIORank.isCartIdxOnThisRank(idx);
};
this->setupBlockData(isCartIdxOnThisRank);
this->setupBlockData(isCartIdxOnThisRank);
this->forceDisableFipOutput_ =
Parameters::Get<Parameters::ForceDisableFluidInPlaceOutput>();
@ -296,6 +297,23 @@ public:
Valgrind::CheckDefined(this->saturation_[phaseIdx][globalDofIdx]);
}
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
if (this->moleFractions_[compIdx].empty()) continue;
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));
}
}
// if (this->regionAvgDensity_.has_value()) {
// // Note: We intentionally exclude effects of rock
// // compressibility by using referencePorosity() here.