OilfieldToolsNet/opm-common
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge pull request #408 from akva2/use_compressed_props

Browse Source 
Use compressed field properties
This commit is contained in:
Markus Blatt 2020-03-02 09:46:16 +01:00 committed by GitHub
commit fb803f1fec
4 changed files with 67 additions and 98 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

60
opm/material/fluidmatrixinteractions/EclEpsGridProperties.hpp
View File

@ -59,23 +59,10 @@ namespace Opm {
*/
namespace {
template <typename T>
std::vector<T> compressed_copy(const std::vector<T>& global_vector, const std::vector<int>& compressedToCartesianElemIdx) {
std::vector<T> compressed = std::vector<T>(compressedToCartesianElemIdx.size());
for (std::size_t active_index = 0; active_index < compressedToCartesianElemIdx.size(); active_index++) {
auto global_index = compressedToCartesianElemIdx[active_index];
compressed[active_index] = global_vector[global_index];
}
return compressed;
}
std::vector<double> try_get(const FieldPropsManager& fp, const std::string& keyword, const std::vector<int>& compressedToCartesianElemIdx) {
std::vector<double> try_get(const FieldPropsManager& fp, const std::string& keyword) {
if (fp.has_double(keyword))
return compressed_copy(fp.get_global_double(keyword), compressedToCartesianElemIdx);
return fp.get_double(keyword);
return {};
}
@ -91,48 +78,47 @@ public:
#if HAVE_ECL_INPUT
EclEpsGridProperties(const Opm::EclipseState& eclState,
bool useImbibition,
const std::vector<int>& compressedToCartesianElemIdx)
bool useImbibition)
{
std::string kwPrefix = useImbibition?"I":"";
const auto& fp = eclState.fieldProps();
if (useImbibition)
compressed_satnum = compressed_copy(fp.get_global_int("IMBNUM"), compressedToCartesianElemIdx);
compressed_satnum = fp.get_int("IMBNUM");
else
compressed_satnum = compressed_copy(fp.get_global_int("SATNUM"), compressedToCartesianElemIdx);
compressed_satnum = fp.get_int("SATNUM");
this->compressed_swl = try_get( fp, kwPrefix+"SWL", compressedToCartesianElemIdx);
this->compressed_sgl = try_get( fp, kwPrefix+"SGL", compressedToCartesianElemIdx);
this->compressed_swcr = try_get( fp, kwPrefix+"SWCR", compressedToCartesianElemIdx);
this->compressed_sgcr = try_get( fp, kwPrefix+"SGCR", compressedToCartesianElemIdx);
this->compressed_sowcr = try_get( fp, kwPrefix+"SOWCR", compressedToCartesianElemIdx);
this->compressed_sogcr = try_get( fp, kwPrefix+"SOGCR", compressedToCartesianElemIdx);
this->compressed_swu = try_get( fp, kwPrefix+"SWU", compressedToCartesianElemIdx);
this->compressed_sgu = try_get( fp, kwPrefix+"SGU", compressedToCartesianElemIdx);
this->compressed_pcw = try_get( fp, kwPrefix+"PCW", compressedToCartesianElemIdx);
this->compressed_pcg = try_get( fp, kwPrefix+"PCG", compressedToCartesianElemIdx);
this->compressed_krw = try_get( fp, kwPrefix+"KRW", compressedToCartesianElemIdx);
this->compressed_kro = try_get( fp, kwPrefix+"KRO", compressedToCartesianElemIdx);
this->compressed_krg = try_get( fp, kwPrefix+"KRG", compressedToCartesianElemIdx);
this->compressed_swl = try_get( fp, kwPrefix+"SWL");
this->compressed_sgl = try_get( fp, kwPrefix+"SGL");
this->compressed_swcr = try_get( fp, kwPrefix+"SWCR");
this->compressed_sgcr = try_get( fp, kwPrefix+"SGCR");
this->compressed_sowcr = try_get( fp, kwPrefix+"SOWCR");
this->compressed_sogcr = try_get( fp, kwPrefix+"SOGCR");
this->compressed_swu = try_get( fp, kwPrefix+"SWU");
this->compressed_sgu = try_get( fp, kwPrefix+"SGU");
this->compressed_pcw = try_get( fp, kwPrefix+"PCW");
this->compressed_pcg = try_get( fp, kwPrefix+"PCG");
this->compressed_krw = try_get( fp, kwPrefix+"KRW");
this->compressed_kro = try_get( fp, kwPrefix+"KRO");
this->compressed_krg = try_get( fp, kwPrefix+"KRG");
// _may_ be needed to calculate the Leverett capillary pressure scaling factor
if (fp.has_double("PORO"))
this->compressed_poro = compressed_copy(fp.get_global_double("PORO"), compressedToCartesianElemIdx);
this->compressed_poro = fp.get_double("PORO");
if (fp.has_double("PERMX"))
this->compressed_permx = compressed_copy(fp.get_global_double("PERMX"), compressedToCartesianElemIdx);
this->compressed_permx = fp.get_double("PERMX");
else
this->compressed_permx = std::vector<double>(this->compressed_satnum.size());
if (fp.has_double("PERMY"))
this->compressed_permy= compressed_copy(fp.get_global_double("PERMY"), compressedToCartesianElemIdx);
this->compressed_permy = fp.get_double("PERMY");
else
this->compressed_permy= this->compressed_permx;
this->compressed_permy = this->compressed_permx;
if (fp.has_double("PERMZ"))
this->compressed_permz= compressed_copy(fp.get_global_double("PERMZ"), compressedToCartesianElemIdx);
this->compressed_permz = fp.get_double("PERMZ");
else
this->compressed_permz= this->compressed_permx;
}

18
opm/material/fluidmatrixinteractions/EclMaterialLawManager.hpp
View File

@ -148,8 +148,7 @@ public:
}
}
void initParamsForElements(const EclipseState& eclState,
const std::vector<int>& compressedToCartesianElemIdx)
void initParamsForElements(const EclipseState& eclState, size_t numCompressedElems)
{
// get the number of saturation regions
const size_t numSatRegions = eclState.runspec().tabdims().getNumSatTables();
@ -169,15 +168,13 @@ public:
readOilWaterEffectiveParameters_(oilWaterEffectiveParamVector_, eclState, satRegionIdx);
}
unsigned numCompressedElems = static_cast<unsigned>(compressedToCartesianElemIdx.size());
// copy the SATNUM grid property. in some cases this is not necessary, but it
// should not require much memory anyway...
satnumRegionArray_.resize(numCompressedElems);
if (eclState.fieldProps().has_int("SATNUM")) {
const auto& satnumRawData = eclState.fieldProps().get_global_int("SATNUM");
const auto& satnumRawData = eclState.fieldProps().get_int("SATNUM");
for (unsigned elemIdx = 0; elemIdx < numCompressedElems; ++elemIdx) {
unsigned cartesianElemIdx = static_cast<unsigned>(compressedToCartesianElemIdx[elemIdx]);
satnumRegionArray_[elemIdx] = satnumRawData[cartesianElemIdx] - 1;
satnumRegionArray_[elemIdx] = satnumRawData[elemIdx] - 1;
}
}
else
@ -187,10 +184,9 @@ public:
// the same as the saturation region (SATNUM)
imbnumRegionArray_ = satnumRegionArray_;
if (eclState.fieldProps().has_int("IMBNUM")) {
const auto& imbnumRawData = eclState.fieldProps().get_global_int("IMBNUM");
const auto& imbnumRawData = eclState.fieldProps().get_int("IMBNUM");
for (unsigned elemIdx = 0; elemIdx < numCompressedElems; ++elemIdx) {
int cartesianElemIdx = compressedToCartesianElemIdx[elemIdx];
imbnumRegionArray_[elemIdx] = imbnumRawData[cartesianElemIdx] - 1;
imbnumRegionArray_[elemIdx] = imbnumRawData[elemIdx] - 1;
}
}
@ -213,7 +209,7 @@ public:
oilWaterScaledImbPointsVector.resize(numCompressedElems);
}
EclEpsGridProperties epsGridProperties(eclState, false, compressedToCartesianElemIdx);
EclEpsGridProperties epsGridProperties(eclState, false);
for (unsigned elemIdx = 0; elemIdx < numCompressedElems; ++elemIdx) {
readGasOilScaledPoints_(gasOilScaledInfoVector,
@ -233,7 +229,7 @@ public:
}
if (enableHysteresis()) {
EclEpsGridProperties epsImbGridProperties(eclState, true, compressedToCartesianElemIdx);
EclEpsGridProperties epsImbGridProperties(eclState, true);
for (unsigned elemIdx = 0; elemIdx < numCompressedElems; ++elemIdx) {
readGasOilScaledPoints_(gasOilScaledImbInfoVector,
gasOilScaledImbPointsVector,

73
opm/material/thermal/EclThermalLawManager.hpp
View File

@ -69,8 +69,7 @@ public:
thermalConductivityApproach_ = ThermalConductionLawParams::undefinedApproach;
}
void initParamsForElements(const Opm::EclipseState& eclState,
const std::vector<int>& compressedToCartesianElemIdx)
void initParamsForElements(const Opm::EclipseState& eclState, size_t numElems)
{
const auto& fp = eclState.fieldProps();
const auto& tableManager = eclState.getTableManager();
@ -79,16 +78,16 @@ public:
bool has_thc = fp.has_double("THCROCK") || fp.has_double("THCOIL") || fp.has_double("THCGAS") || fp.has_double("THCWATER");
if (has_heatcr)
initHeatcr_(eclState, compressedToCartesianElemIdx);
initHeatcr_(eclState, numElems);
else if (tableManager.hasTables("SPECROCK"))
initSpecrock_(eclState, compressedToCartesianElemIdx);
initSpecrock_(eclState, numElems);
else
initNullRockEnergy_();
if (has_thconr)
initThconr_(eclState, compressedToCartesianElemIdx);
initThconr_(eclState, numElems);
else if (has_thc)
initThc_(eclState, compressedToCartesianElemIdx);
initThc_(eclState, numElems);
else
initNullCond_();
}
@ -139,7 +138,7 @@ private:
* \brief Initialize the parameters for the solid energy law using using HEATCR and friends.
*/
void initHeatcr_(const Opm::EclipseState& eclState,
const std::vector<int>& compressedToCartesianElemIdx)
size_t numElems)
{
solidEnergyApproach_ = SolidEnergyLawParams::heatcrApproach;
// actually the value of the reference temperature does not matter for energy
@ -147,18 +146,16 @@ private:
HeatcrLawParams::setReferenceTemperature(FluidSystem::surfaceTemperature);
const auto& fp = eclState.fieldProps();
const std::vector<double>& heatcrData = fp.get_global_double("HEATCR");
const std::vector<double>& heatcrtData = fp.get_global_double("HEATCRT");
unsigned numElems = compressedToCartesianElemIdx.size();
const std::vector<double>& heatcrData = fp.get_double("HEATCR");
const std::vector<double>& heatcrtData = fp.get_double("HEATCRT");
solidEnergyLawParams_.resize(numElems);
for (unsigned elemIdx = 0; elemIdx < numElems; ++elemIdx) {
auto& elemParam = solidEnergyLawParams_[elemIdx];
elemParam.setSolidEnergyApproach(SolidEnergyLawParams::heatcrApproach);
auto& heatcrElemParams = elemParam.template getRealParams<SolidEnergyLawParams::heatcrApproach>();
unsigned cartesianElemIdx = compressedToCartesianElemIdx[elemIdx];
heatcrElemParams.setReferenceRockHeatCapacity(heatcrData[cartesianElemIdx]);
heatcrElemParams.setDRockHeatCapacity_dT(heatcrtData[cartesianElemIdx]);
heatcrElemParams.setReferenceRockHeatCapacity(heatcrData[elemIdx]);
heatcrElemParams.setDRockHeatCapacity_dT(heatcrtData[elemIdx]);
heatcrElemParams.finalize();
elemParam.finalize();
}
@ -168,20 +165,18 @@ private:
* \brief Initialize the parameters for the solid energy law using using SPECROCK and friends.
*/
void initSpecrock_(const Opm::EclipseState& eclState,
const std::vector<int>& compressedToCartesianElemIdx)
size_t numElems)
{
solidEnergyApproach_ = SolidEnergyLawParams::specrockApproach;
// initialize the element index -> SATNUM index mapping
const auto& fp = eclState.fieldProps();
const std::vector<int>& satnumData = fp.get_global_int("SATNUM");
elemToSatnumIdx_.resize(compressedToCartesianElemIdx.size());
for (unsigned elemIdx = 0; elemIdx < compressedToCartesianElemIdx.size(); ++ elemIdx) {
unsigned cartesianElemIdx = compressedToCartesianElemIdx[elemIdx];
const std::vector<int>& satnumData = fp.get_int("SATNUM");
elemToSatnumIdx_.resize(numElems);
for (unsigned elemIdx = 0; elemIdx < numElems; ++ elemIdx) {
// satnumData contains Fortran-style indices, i.e., they start with 1 instead
// of 0!
elemToSatnumIdx_[elemIdx] = satnumData[cartesianElemIdx] - 1;
elemToSatnumIdx_[elemIdx] = satnumData[elemIdx] - 1;
}
// internalize the SPECROCK table
unsigned numSatRegions = eclState.runspec().tabdims().getNumSatTables();
@ -219,7 +214,7 @@ private:
* \brief Initialize the parameters for the thermal conduction law using THCONR and friends.
*/
void initThconr_(const Opm::EclipseState& eclState,
const std::vector<int>& compressedToCartesianElemIdx)
size_t numElems)
{
thermalConductivityApproach_ = ThermalConductionLawParams::thconrApproach;
@ -227,21 +222,19 @@ private:
std::vector<double> thconrData;
std::vector<double> thconsfData;
if (fp.has_double("THCONR"))
thconrData = fp.get_global_double("THCONR");
thconrData = fp.get_double("THCONR");
if (fp.has_double("THCONSF"))
thconsfData = fp.get_global_double("THCONSF");
thconsfData = fp.get_double("THCONSF");
unsigned numElems = compressedToCartesianElemIdx.size();
thermalConductionLawParams_.resize(numElems);
for (unsigned elemIdx = 0; elemIdx < numElems; ++elemIdx) {
auto& elemParams = thermalConductionLawParams_[elemIdx];
elemParams.setThermalConductionApproach(ThermalConductionLawParams::thconrApproach);
auto& thconrElemParams = elemParams.template getRealParams<ThermalConductionLawParams::thconrApproach>();
int cartElemIdx = compressedToCartesianElemIdx[elemIdx];
double thconr = thconrData.empty() ? 0.0 : thconrData[cartElemIdx];
double thconsf = thconsfData.empty() ? 0.0 : thconsfData[cartElemIdx];
double thconr = thconrData.empty() ? 0.0 : thconrData[elemIdx];
double thconsf = thconsfData.empty() ? 0.0 : thconsfData[elemIdx];
thconrElemParams.setReferenceTotalThermalConductivity(thconr);
thconrElemParams.setDTotalThermalConductivity_dSg(thconsf);
@ -254,7 +247,7 @@ private:
* \brief Initialize the parameters for the thermal conduction law using THCROCK and friends.
*/
void initThc_(const Opm::EclipseState& eclState,
const std::vector<int>& compressedToCartesianElemIdx)
size_t numElems)
{
thermalConductivityApproach_ = ThermalConductionLawParams::thcApproach;
@ -262,35 +255,33 @@ private:
std::vector<double> thcrockData;
std::vector<double> thcoilData;
std::vector<double> thcgasData;
std::vector<double> thcwaterData = fp.get_global_double("THCWATER");
std::vector<double> thcwaterData = fp.get_double("THCWATER");
if (fp.has_double("THCROCK"))
thcrockData = fp.get_global_double("THCROCK");
thcrockData = fp.get_double("THCROCK");
if (fp.has_double("THCOIL"))
thcoilData = fp.get_global_double("THCOIL");
thcoilData = fp.get_double("THCOIL");
if (fp.has_double("THCGAS"))
thcgasData = fp.get_global_double("THCGAS");
thcgasData = fp.get_double("THCGAS");
if (fp.has_double("THCWATER"))
thcwaterData = fp.get_global_double("THCWATER");
thcwaterData = fp.get_double("THCWATER");
const std::vector<double>& poroData = fp.get_global_double("PORO");
const std::vector<double>& poroData = fp.get_double("PORO");
unsigned numElems = compressedToCartesianElemIdx.size();
thermalConductionLawParams_.resize(numElems);
for (unsigned elemIdx = 0; elemIdx < numElems; ++elemIdx) {
auto& elemParams = thermalConductionLawParams_[elemIdx];
elemParams.setThermalConductionApproach(ThermalConductionLawParams::thcApproach);
auto& thcElemParams = elemParams.template getRealParams<ThermalConductionLawParams::thcApproach>();
int cartElemIdx = compressedToCartesianElemIdx[elemIdx];
thcElemParams.setPorosity(poroData[cartElemIdx]);
double thcrock = thcrockData.empty() ? 0.0 : thcrockData[cartElemIdx];
double thcoil = thcoilData.empty() ? 0.0 : thcoilData[cartElemIdx];
double thcgas = thcgasData.empty() ? 0.0 : thcgasData[cartElemIdx];
double thcwater = thcwaterData.empty() ? 0.0 : thcwaterData[cartElemIdx];
thcElemParams.setPorosity(poroData[elemIdx]);
double thcrock = thcrockData.empty() ? 0.0 : thcrockData[elemIdx];
double thcoil = thcoilData.empty() ? 0.0 : thcoilData[elemIdx];
double thcgas = thcgasData.empty() ? 0.0 : thcgasData[elemIdx];
double thcwater = thcwaterData.empty() ? 0.0 : thcwaterData[elemIdx];
thcElemParams.setThcrock(thcrock);
thcElemParams.setThcoil(thcoil);
thcElemParams.setThcgas(thcgas);

14
tests/test_eclmateriallawmanager.cpp
View File

@ -447,14 +447,10 @@ inline void testAll()
const auto& eclGrid = eclState.getInputGrid();
size_t n = eclGrid.getCartesianSize();
std::vector<int> compressedToCartesianIdx(n);
for (size_t i = 0; i < n; ++ i)
compressedToCartesianIdx[i] = static_cast<int>(i);
MaterialLawManager materialLawManager;
materialLawManager.initFromDeck(deck, eclState);
materialLawManager.initParamsForElements(eclState, compressedToCartesianIdx);
materialLawManager.initParamsForElements(eclState, n);
if (materialLawManager.enableEndPointScaling())
throw std::logic_error("Discrepancy between the deck and the EclMaterialLawManager");
@ -468,7 +464,7 @@ inline void testAll()
Opm::EclMaterialLawManager<MaterialTraits> fam2MaterialLawManager;
fam2MaterialLawManager.initFromDeck(fam2Deck, fam2EclState);
fam2MaterialLawManager.initParamsForElements(fam2EclState, compressedToCartesianIdx);
fam2MaterialLawManager.initParamsForElements(fam2EclState, n);
if (fam2MaterialLawManager.enableEndPointScaling())
throw std::logic_error("Discrepancy between the deck and the EclMaterialLawManager");
@ -481,7 +477,7 @@ inline void testAll()
Opm::EclMaterialLawManager<MaterialTraits> hysterMaterialLawManager;
hysterMaterialLawManager.initFromDeck(hysterDeck, hysterEclState);
hysterMaterialLawManager.initParamsForElements(hysterEclState, compressedToCartesianIdx);
hysterMaterialLawManager.initParamsForElements(hysterEclState, n);
if (hysterMaterialLawManager.enableEndPointScaling())
throw std::logic_error("Discrepancy between the deck and the EclMaterialLawManager");
@ -574,14 +570,14 @@ inline void testAll()
MaterialLawManager fam1materialLawManager;
fam1materialLawManager.initFromDeck(fam1Deck, fam1EclState);
fam1materialLawManager.initParamsForElements(fam1EclState, compressedToCartesianIdx);
fam1materialLawManager.initParamsForElements(fam1EclState, n);
const auto fam2Deck = parser.parseString(fam2DeckStringGasOil);
const Opm::EclipseState fam2EclState(fam2Deck);
Opm::EclMaterialLawManager<MaterialTraits> fam2MaterialLawManager;
fam2MaterialLawManager.initFromDeck(fam2Deck, fam2EclState);
fam2MaterialLawManager.initParamsForElements(fam2EclState, compressedToCartesianIdx);
fam2MaterialLawManager.initParamsForElements(fam2EclState, n);
for (unsigned elemIdx = 0; elemIdx < n; ++ elemIdx) {
for (int i = 0; i < 100; ++ i) {