From 82658c92d07f93721353e8052a998dd94c274469 Mon Sep 17 00:00:00 2001 From: Robert Kloefkorn Date: Thu, 9 Feb 2017 13:54:42 +0100 Subject: [PATCH] Removal of SimulatorFullyImplicitBlackoilOutputEbos.{h,c}pp. All simulators now use SimulationDataContainer to store intermediate data that is passed to the output Solution container. This is in cases not the most efficient way, but it's unified to avoid errors from code duplication. --- CMakeLists_files.cmake | 2 - opm/autodiff/BlackoilModelBase.hpp | 27 +- opm/autodiff/BlackoilModelBase_impl.hpp | 3 - opm/autodiff/BlackoilModelEbos.hpp | 179 ++++- opm/autodiff/BlackoilModelEnums.hpp | 41 + opm/autodiff/BlackoilSequentialModel.hpp | 6 + opm/autodiff/FlowMainEbos.hpp | 5 +- .../SimulatorFullyImplicitBlackoilEbos.hpp | 7 +- .../SimulatorFullyImplicitBlackoilOutput.hpp | 260 ++++-- ...mulatorFullyImplicitBlackoilOutputEbos.cpp | 307 ------- ...mulatorFullyImplicitBlackoilOutputEbos.hpp | 749 ------------------ .../WellStateFullyImplicitBlackoil.hpp | 8 +- ...FullyImplicitCompressiblePolymerSolver.hpp | 35 +- 13 files changed, 425 insertions(+), 1204 deletions(-) delete mode 100644 opm/autodiff/SimulatorFullyImplicitBlackoilOutputEbos.cpp delete mode 100644 opm/autodiff/SimulatorFullyImplicitBlackoilOutputEbos.hpp diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index c65af9dc4..2c821009f 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -35,7 +35,6 @@ list (APPEND MAIN_SOURCE_FILES opm/autodiff/ImpesTPFAAD.cpp opm/autodiff/moduleVersion.cpp opm/autodiff/SimulatorFullyImplicitBlackoilOutput.cpp - opm/autodiff/SimulatorFullyImplicitBlackoilOutputEbos.cpp opm/autodiff/SimulatorIncompTwophaseAd.cpp opm/autodiff/TransportSolverTwophaseAd.cpp opm/autodiff/BlackoilPropsAdFromDeck.cpp @@ -203,7 +202,6 @@ list (APPEND PUBLIC_HEADER_FILES opm/autodiff/SimulatorFullyImplicitBlackoilSolvent_impl.hpp opm/autodiff/SimulatorFullyImplicitBlackoilMultiSegment.hpp opm/autodiff/SimulatorFullyImplicitBlackoilMultiSegment_impl.hpp - opm/autodiff/SimulatorFullyImplicitBlackoilOutputEbos.hpp opm/autodiff/SimulatorIncompTwophaseAd.hpp opm/autodiff/SimulatorSequentialBlackoil.hpp opm/autodiff/TransportSolverTwophaseAd.hpp diff --git a/opm/autodiff/BlackoilModelBase.hpp b/opm/autodiff/BlackoilModelBase.hpp index 4162ede67..1d52fd63f 100644 --- a/opm/autodiff/BlackoilModelBase.hpp +++ b/opm/autodiff/BlackoilModelBase.hpp @@ -39,6 +39,8 @@ #include #include +#include + #include struct Wells; @@ -50,7 +52,6 @@ namespace Opm { class RockCompressibility; class NewtonIterationBlackoilInterface; class VFPProperties; - class SimulationDataContainer; /// Traits to encapsulate the types used by classes using or /// extending this model. Forward declared here, must be @@ -92,24 +93,21 @@ namespace Opm { ADB mob; // Phase mobility (per cell) }; - struct SimulatorData { + struct SimulatorData : public Opm::FIPDataEnums { SimulatorData(int num_phases); - enum FipId { - FIP_AQUA = Opm::Water, - FIP_LIQUID = Opm::Oil, - FIP_VAPOUR = Opm::Gas, - FIP_DISSOLVED_GAS = 3, - FIP_VAPORIZED_OIL = 4, - FIP_PV = 5, //< Pore volume - FIP_WEIGHTED_PRESSURE = 6 - }; + using Opm::FIPDataEnums :: FipId ; + using Opm::FIPDataEnums :: fipValues ; + std::vector rq; ADB rsSat; // Saturated gas-oil ratio ADB rvSat; // Saturated oil-gas ratio - std::array fip; + + std::array fip; }; + typedef Opm::FIPData FIPDataType; + typedef typename ModelTraits::ReservoirState ReservoirState; typedef typename ModelTraits::WellState WellState; typedef typename ModelTraits::ModelParameters ModelParameters; @@ -264,6 +262,11 @@ namespace Opm { return sd_; } + /// Return fluid-in-place data (for output functionality) + FIPDataType getFIPData() const { + return FIPDataType( sd_.fip ); + } + /// Compute fluid in place. /// \param[in] ReservoirState /// \param[in] FIPNUM for active cells not global cells. diff --git a/opm/autodiff/BlackoilModelBase_impl.hpp b/opm/autodiff/BlackoilModelBase_impl.hpp index 4dd3b0bba..56bec6762 100644 --- a/opm/autodiff/BlackoilModelBase_impl.hpp +++ b/opm/autodiff/BlackoilModelBase_impl.hpp @@ -51,7 +51,6 @@ #include #include #include -#include #include @@ -389,7 +388,6 @@ typedef Eigen::Array int BlackoilModelBase:: @@ -445,7 +443,6 @@ typedef Eigen::Array BlackoilModelBase:: ReservoirResidualQuant::ReservoirResidualQuant() diff --git a/opm/autodiff/BlackoilModelEbos.hpp b/opm/autodiff/BlackoilModelEbos.hpp index 1ad40735c..8245e1336 100644 --- a/opm/autodiff/BlackoilModelEbos.hpp +++ b/opm/autodiff/BlackoilModelEbos.hpp @@ -143,19 +143,7 @@ namespace Opm { using RateConverterType = RateConverter:: SurfaceToReservoirVoidage >; - struct FIPData { - enum FipId { - FIP_AQUA = Opm::Water, - FIP_LIQUID = Opm::Oil, - FIP_VAPOUR = Opm::Gas, - FIP_DISSOLVED_GAS = 3, - FIP_VAPORIZED_OIL = 4, - FIP_PV = 5, //< Pore volume - FIP_WEIGHTED_PRESSURE = 6 - }; - static const int fipValues = FIP_WEIGHTED_PRESSURE + 1 ; - std::array, fipValues> fip; - }; + typedef Opm::FIPData FIPDataType; // --------- Public methods --------- @@ -1022,7 +1010,7 @@ namespace Opm { const auto& pv = geo_.poreVolume(); const int maxnp = Opm::BlackoilPhases::MaxNumPhases; - for (int i = 0; i> values(dims, std::vector(FIPData::fipValues,0.0)); + std::vector> values(dims, std::vector(FIPDataType::fipValues,0.0)); std::vector hcpv(dims, 0.0); std::vector pres(dims, 0.0); @@ -1082,8 +1070,8 @@ namespace Opm { for (int c = 0; c < nc; ++c) { const int region = fipnum[c] - 1; if (region != -1) { - values[region][FIPData::FIP_DISSOLVED_GAS] += fip_.fip[FIPData::FIP_DISSOLVED_GAS][c]; - values[region][FIPData::FIP_VAPORIZED_OIL] += fip_.fip[FIPData::FIP_VAPORIZED_OIL][c]; + values[region][FIPDataType::FIP_DISSOLVED_GAS] += fip_.fip[FIPDataType::FIP_DISSOLVED_GAS][c]; + values[region][FIPDataType::FIP_VAPORIZED_OIL] += fip_.fip[FIPDataType::FIP_VAPORIZED_OIL][c]; } } } @@ -1102,7 +1090,7 @@ namespace Opm { const int region = fipnum[c] - 1; if (region != -1) { - fip_.fip[FIPData::FIP_PV][c] = pv[c]; + fip_.fip[FIPDataType::FIP_PV][c] = pv[c]; const auto& intQuants = *ebosSimulator_.model().cachedIntensiveQuantities(c, /*timeIdx=*/0); const auto& fs = intQuants.fluidState(); const double hydrocarbon = fs.saturation(FluidSystem::oilPhaseIdx).value() + fs.saturation(FluidSystem::gasPhaseIdx).value(); @@ -1110,13 +1098,13 @@ namespace Opm { //Compute hydrocarbon pore volume weighted average pressure. //If we have no hydrocarbon in region, use pore volume weighted average pressure instead if (hcpv[region] != 0) { - fip_.fip[FIPData::FIP_WEIGHTED_PRESSURE][c] = pv[c] * fs.pressure(FluidSystem::oilPhaseIdx).value() * hydrocarbon / hcpv[region]; + fip_.fip[FIPDataType::FIP_WEIGHTED_PRESSURE][c] = pv[c] * fs.pressure(FluidSystem::oilPhaseIdx).value() * hydrocarbon / hcpv[region]; } else { - fip_.fip[FIPData::FIP_WEIGHTED_PRESSURE][c] = pres[region] / pv[c]; + fip_.fip[FIPDataType::FIP_WEIGHTED_PRESSURE][c] = pres[region] / pv[c]; } - values[region][FIPData::FIP_PV] += fip_.fip[FIPData::FIP_PV][c]; - values[region][FIPData::FIP_WEIGHTED_PRESSURE] += fip_.fip[FIPData::FIP_WEIGHTED_PRESSURE][c]; + values[region][FIPDataType::FIP_PV] += fip_.fip[FIPDataType::FIP_PV][c]; + values[region][FIPDataType::FIP_WEIGHTED_PRESSURE] += fip_.fip[FIPDataType::FIP_WEIGHTED_PRESSURE][c]; } } } @@ -1131,7 +1119,7 @@ namespace Opm { auto comm = pinfo.communicator(); // Compute the global dims value and resize values accordingly. dims = comm.max(dims); - values.resize(dims, std::vector(FIPData::fipValues,0.0)); + values.resize(dims, std::vector(FIPDataType::fipValues,0.0)); //Accumulate phases for each region for (int phase = 0; phase < maxnp; ++phase) { @@ -1148,8 +1136,8 @@ namespace Opm { for (int c = 0; c < nc; ++c) { const int region = fipnum[c] - 1; if (region != -1 && mask[c]) { - values[region][FIPData::FIP_DISSOLVED_GAS] += fip_.fip[FIPData::FIP_DISSOLVED_GAS][c]; - values[region][FIPData::FIP_VAPORIZED_OIL] += fip_.fip[FIPData::FIP_VAPORIZED_OIL][c]; + values[region][FIPDataType::FIP_DISSOLVED_GAS] += fip_.fip[FIPDataType::FIP_DISSOLVED_GAS][c]; + values[region][FIPDataType::FIP_VAPORIZED_OIL] += fip_.fip[FIPDataType::FIP_VAPORIZED_OIL][c]; } } } @@ -1174,19 +1162,19 @@ namespace Opm { for (int c = 0; c < nc; ++c) { const int region = fipnum[c] - 1; if (region != -1 && mask[c]) { - fip_.fip[FIPData::FIP_PV][c] = pv[c]; + fip_.fip[FIPDataType::FIP_PV][c] = pv[c]; const auto& intQuants = *ebosSimulator_.model().cachedIntensiveQuantities(c, /*timeIdx=*/0); const auto& fs = intQuants.fluidState(); const double hydrocarbon = fs.saturation(FluidSystem::oilPhaseIdx).value() + fs.saturation(FluidSystem::gasPhaseIdx).value(); if (hcpv[region] != 0) { - fip_.fip[FIPData::FIP_WEIGHTED_PRESSURE][c] = pv[c] * fs.pressure(FluidSystem::oilPhaseIdx).value() * hydrocarbon / hcpv[region]; + fip_.fip[FIPDataType::FIP_WEIGHTED_PRESSURE][c] = pv[c] * fs.pressure(FluidSystem::oilPhaseIdx).value() * hydrocarbon / hcpv[region]; } else { - fip_.fip[FIPData::FIP_WEIGHTED_PRESSURE][c] = pres[region] / pv[c]; + fip_.fip[FIPDataType::FIP_WEIGHTED_PRESSURE][c] = pres[region] / pv[c]; } - values[region][FIPData::FIP_PV] += fip_.fip[FIPData::FIP_PV][c]; - values[region][FIPData::FIP_WEIGHTED_PRESSURE] += fip_.fip[FIPData::FIP_WEIGHTED_PRESSURE][c]; + values[region][FIPDataType::FIP_PV] += fip_.fip[FIPDataType::FIP_PV][c]; + values[region][FIPDataType::FIP_WEIGHTED_PRESSURE] += fip_.fip[FIPDataType::FIP_WEIGHTED_PRESSURE][c]; } } @@ -1206,7 +1194,130 @@ namespace Opm { return values; } - const FIPData& getFIPData() const { + SimulationDataContainer getSimulatorData () const + { + typedef std::vector VectorType; + + const auto& ebosModel = ebosSimulator().model(); + const auto& phaseUsage = fluid_.phaseUsage(); + + // extract everything which can possibly be written to disk + const int numCells = ebosModel.numGridDof(); + const int numPhazes = numPhases(); + + SimulationDataContainer simData( numCells, 0, numPhazes ); + + //Get shorthands for water, oil, gas + const int aqua_active = phaseUsage.phase_used[Opm::PhaseUsage::Aqua]; + const int liquid_active = phaseUsage.phase_used[Opm::PhaseUsage::Liquid]; + const int vapour_active = phaseUsage.phase_used[Opm::PhaseUsage::Vapour]; + + const int aqua_pos = phaseUsage.phase_pos[ Opm::PhaseUsage::Aqua ]; + const int liquid_pos = phaseUsage.phase_pos[ Opm::PhaseUsage::Liquid ]; + const int vapour_pos = phaseUsage.phase_pos[ Opm::PhaseUsage::Vapour ]; + + VectorType zero; + + VectorType& pressureOil = simData.pressure(); + VectorType& temperature = simData.temperature(); + VectorType& saturation = simData.saturation(); + + // WATER + if( aqua_active ) { + simData.registerCellData("1OVERBW", 1 ); + simData.registerCellData("WAT_DEN", 1 ); + simData.registerCellData("WAT_VISC", 1 ); + simData.registerCellData("WATKR", 1 ); + } + + VectorType& bWater = aqua_active ? simData.getCellData( "1OVERBW" ) : zero; + VectorType& rhoWater = aqua_active ? simData.getCellData( "WAT_DEN" ) : zero; + VectorType& muWater = aqua_active ? simData.getCellData( "WAT_VISC" ) : zero; + VectorType& krWater = aqua_active ? simData.getCellData( "WATKR" ) : zero; + + // OIL + if( liquid_active ) { + simData.registerCellData("1OVERBO", 1 ); + simData.registerCellData("OIL_DEN", 1 ); + simData.registerCellData("OIL_VISC", 1 ); + simData.registerCellData("OILKR", 1 ); + } + + VectorType& bOil = liquid_active ? simData.getCellData( "1OVERBO" ) : zero; + VectorType& rhoOil = liquid_active ? simData.getCellData( "OIL_DEN" ) : zero; + VectorType& muOil = liquid_active ? simData.getCellData( "OIL_VISC" ) : zero; + VectorType& krOil = liquid_active ? simData.getCellData( "OILKR" ) : zero; + + // GAS + if( vapour_active ) { + simData.registerCellData("1OVERBG", 1 ); + simData.registerCellData("GAS_DEN", 1 ); + simData.registerCellData("GAS_VISC", 1 ); + simData.registerCellData("GASKR", 1 ); + } + + VectorType& bGas = vapour_active ? simData.getCellData( "1OVERBG" ) : zero; + VectorType& rhoGas = vapour_active ? simData.getCellData( "GAS_DEN" ) : zero; + VectorType& muGas = vapour_active ? simData.getCellData( "GAS_VISC" ) : zero; + VectorType& krGas = vapour_active ? simData.getCellData( "GASKR" ) : zero; + + simData.registerCellData( BlackoilState::GASOILRATIO, 1 ); + simData.registerCellData( BlackoilState::RV, 1 ); + simData.registerCellData("RSSAT", 1 ); + simData.registerCellData("RVSAT", 1 ); + + VectorType& Rs = simData.getCellData( BlackoilState::GASOILRATIO ); + VectorType& Rv = simData.getCellData( BlackoilState::RV ); + VectorType& RsSat = simData.getCellData( "RSSAT" ); + VectorType& RvSat = simData.getCellData( "RVSAT" ); + + for (int cellIdx = 0; cellIdx < numCells; ++cellIdx) { + const auto& intQuants = *ebosModel.cachedIntensiveQuantities(cellIdx, /*timeIdx=*/0); + const auto& fs = intQuants.fluidState(); + + const int satIdx = cellIdx * numPhazes; + + pressureOil[cellIdx] = fs.pressure(FluidSystem::oilPhaseIdx).value(); + + temperature[cellIdx] = fs.temperature(FluidSystem::oilPhaseIdx).value(); + + if (aqua_active) { + saturation[ satIdx + aqua_pos ] = fs.saturation(FluidSystem::waterPhaseIdx).value(); + bWater[cellIdx] = fs.invB(FluidSystem::waterPhaseIdx).value(); + rhoWater[cellIdx] = fs.density(FluidSystem::waterPhaseIdx).value(); + muWater[cellIdx] = fs.viscosity(FluidSystem::waterPhaseIdx).value(); + krWater[cellIdx] = intQuants.relativePermeability(FluidSystem::waterPhaseIdx).value(); + } + if (vapour_active) { + saturation[ satIdx + vapour_pos ] = fs.saturation(FluidSystem::gasPhaseIdx).value(); + bGas[cellIdx] = fs.invB(FluidSystem::gasPhaseIdx).value(); + rhoGas[cellIdx] = fs.density(FluidSystem::gasPhaseIdx).value(); + muGas[cellIdx] = fs.viscosity(FluidSystem::gasPhaseIdx).value(); + krGas[cellIdx] = intQuants.relativePermeability(FluidSystem::gasPhaseIdx).value(); + Rs[cellIdx] = fs.Rs().value(); + Rv[cellIdx] = fs.Rv().value(); + RsSat[cellIdx] = FluidSystem::saturatedDissolutionFactor(fs, + FluidSystem::oilPhaseIdx, + intQuants.pvtRegionIndex(), + /*maxOilSaturation=*/1.0).value(); + RvSat[cellIdx] = FluidSystem::saturatedDissolutionFactor(fs, + FluidSystem::gasPhaseIdx, + intQuants.pvtRegionIndex(), + /*maxOilSaturation=*/1.0).value(); + } + + // oil is always active + saturation[ satIdx + liquid_pos ] = fs.saturation(FluidSystem::oilPhaseIdx).value(); + bOil[cellIdx] = fs.invB(FluidSystem::oilPhaseIdx).value(); + rhoOil[cellIdx] = fs.density(FluidSystem::oilPhaseIdx).value(); + muOil[cellIdx] = fs.viscosity(FluidSystem::oilPhaseIdx).value(); + krOil[cellIdx] = intQuants.relativePermeability(FluidSystem::oilPhaseIdx).value(); + } + + return std::move(simData); + } + + const FIPDataType& getFIPData() const { return fip_; } @@ -1254,7 +1365,7 @@ namespace Opm { std::vector> residual_norms_history_; double current_relaxation_; BVector dx_old_; - mutable FIPData fip_; + mutable FIPDataType fip_; diff --git a/opm/autodiff/BlackoilModelEnums.hpp b/opm/autodiff/BlackoilModelEnums.hpp index c5ffe42e0..3b01eff06 100644 --- a/opm/autodiff/BlackoilModelEnums.hpp +++ b/opm/autodiff/BlackoilModelEnums.hpp @@ -22,6 +22,9 @@ #include +#include +#include + namespace Opm { @@ -45,6 +48,44 @@ namespace Opm Next // For extension. }; + struct FIPDataEnums { + + enum FipId { + FIP_AQUA = Opm::Water, + FIP_LIQUID = Opm::Oil, + FIP_VAPOUR = Opm::Gas, + FIP_DISSOLVED_GAS = 3, + FIP_VAPORIZED_OIL = 4, + FIP_PV = 5, //< Pore volume + FIP_WEIGHTED_PRESSURE = 6 + }; + + static const int fipValues = FIP_WEIGHTED_PRESSURE + 1 ; + }; + + class FIPData : public FIPDataEnums + { + public: + typedef std::vector VectorType; + + using FIPDataEnums :: FipId; + using FIPDataEnums :: fipValues ; + std::array< VectorType, fipValues> fip; + + // default constructor + FIPData() {} + + // initialize from array of Eigen vectors (or std::vectors) + template + explicit FIPData( const std::array< V, fipValues>& otherFip ) + { + // copy fip vector from V to std::vector + for( int i=0; i PressureModel; typedef BlackoilTransportModel TransportModel; typedef typename TransportModel::SimulatorData SimulatorData; + typedef typename TransportModel::FIPDataType FIPDataType; /// Construct the model. It will retain references to the /// arguments of this functions, and they are expected to @@ -276,6 +277,11 @@ namespace Opm { return transport_solver_.model().getSimulatorData(); } + /// Return fluid-in-place data (for output functionality) + FIPDataType getFIPData() const { + return transport_solver_.model().getFIPData(); + } + protected: typedef NonlinearSolver PressureSolver; diff --git a/opm/autodiff/FlowMainEbos.hpp b/opm/autodiff/FlowMainEbos.hpp index 1837d2c70..b681cc573 100644 --- a/opm/autodiff/FlowMainEbos.hpp +++ b/opm/autodiff/FlowMainEbos.hpp @@ -556,12 +556,11 @@ namespace Opm // output_writer_ void setupOutputWriter() { - const PhaseUsage pu = Opm::phaseUsageFromDeck(deck()); output_writer_.reset(new OutputWriter(grid(), param_, eclState(), - *eclIO_, - pu)); + std::move( eclIO_ ), + Opm::phaseUsageFromDeck(deck())) ); } // Run the simulator. diff --git a/opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp b/opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp index 871062fa0..fea0439a5 100644 --- a/opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp +++ b/opm/autodiff/SimulatorFullyImplicitBlackoilEbos.hpp @@ -22,7 +22,8 @@ #define OPM_SIMULATORFULLYIMPLICITBLACKOILEBOS_HEADER_INCLUDED //#include -#include +//#include +#include #include #include #include @@ -59,7 +60,7 @@ public: typedef WellStateFullyImplicitBlackoilDense WellState; typedef BlackoilState ReservoirState; - typedef BlackoilOutputWriterEbos OutputWriter; + typedef BlackoilOutputWriter OutputWriter; typedef BlackoilModelEbos Model; typedef BlackoilModelParameters ModelParameters; typedef NonlinearSolver Solver; @@ -98,7 +99,7 @@ public: const bool has_disgas, const bool has_vapoil, const EclipseState& /* eclState */, - BlackoilOutputWriterEbos& output_writer, + OutputWriter& output_writer, const std::unordered_set& defunct_well_names) : ebosSimulator_(ebosSimulator), param_(param), diff --git a/opm/autodiff/SimulatorFullyImplicitBlackoilOutput.hpp b/opm/autodiff/SimulatorFullyImplicitBlackoilOutput.hpp index 049918756..b25a681d4 100644 --- a/opm/autodiff/SimulatorFullyImplicitBlackoilOutput.hpp +++ b/opm/autodiff/SimulatorFullyImplicitBlackoilOutput.hpp @@ -280,11 +280,11 @@ namespace Opm const int desiredReportStep); - template - void initFromRestartFile(const PhaseUsage& phaseusage, + template + void initFromRestartFile(const PhaseUsage& phaseUsage, const Grid& grid, SimulationDataContainer& simulatorstate, - WellStateFullyImplicitBlackoil& wellstate); + WellStateFullyImplicitBlackOel& wellstate); bool isRestart() const; @@ -349,7 +349,7 @@ namespace Opm Opm::OpmLog::warning("Parallel Output Config", "Velocity output for matlab is broken in parallel."); } - + if( parallelOutput_->isIORank() ) { if ( output_matlab ) @@ -395,13 +395,13 @@ namespace Opm } - template + template inline void BlackoilOutputWriter:: - initFromRestartFile( const PhaseUsage& phaseusage, + initFromRestartFile( const PhaseUsage& phaseUsage, const Grid& grid, SimulationDataContainer& simulatorstate, - WellStateFullyImplicitBlackoil& wellstate) + WellStateFullyImplicitBlackOel& wellstate) { std::map solution_keys {{"PRESSURE" , UnitSystem::measure::pressure}, {"SWAT" , UnitSystem::measure::identity}, @@ -430,11 +430,11 @@ namespace Opm std::unordered_set()); const Wells* wells = wellsmanager.c_wells(); - wellstate.resize(wells, simulatorstate); //Resize for restart step + wellstate.resize(wells, simulatorstate, phaseUsage ); //Resize for restart step auto state = eclIO_->loadRestart(solution_keys); - solutionToSim( state.first, phaseusage, simulatorstate ); - wellsToState( state.second, phaseusage, wellstate ); + solutionToSim( state.first, phaseUsage, simulatorstate ); + wellsToState( state.second, phaseUsage, wellstate ); } @@ -443,50 +443,130 @@ namespace Opm namespace detail { - /** - * Converts an ADB::V into a standard vector by copy - */ - inline std::vector adbVToDoubleVector(const Opm::AutoDiffBlock::V& adb_v) { - std::vector vec(adb_v.data(), adb_v.data() + adb_v.size()); - return vec; + + template + void addToSimData( SimulationDataContainer& simData, + const std::string& name, + const V& vec ) + { + typedef std::vector< double > OutputVectorType; + + // get data map + auto& dataMap = simData.cellData(); + + // insert name,vector into data map + dataMap.insert( std::make_pair( name, OutputVectorType( vec.data(), vec.data() + vec.size() ) ) ); + } + + template + void addToSimData( SimulationDataContainer& simData, + const std::string& name, + const AutoDiffBlock& adb ) + { + // forward value of ADB to output + addToSimData( simData, name, adb.value() ); } - /** - * Converts an ADB into a standard vector by copy - */ - inline std::vector adbToDoubleVector(const Opm::AutoDiffBlock& adb) { - return adbVToDoubleVector(adb.value()); + // this method basically converts all Eigen vectors to std::vectors + // stored in a SimulationDataContainer + template + SimulationDataContainer + convertToSimulationDataContainer( const SimulatorData& sd, + const SimulationDataContainer& localState, + const Opm::PhaseUsage& phaseUsage ) + { + // copy local state and then add missing data + SimulationDataContainer simData( localState ); + + //Get shorthands for water, oil, gas + const int aqua_active = phaseUsage.phase_used[Opm::PhaseUsage::Aqua]; + const int liquid_active = phaseUsage.phase_used[Opm::PhaseUsage::Liquid]; + const int vapour_active = phaseUsage.phase_used[Opm::PhaseUsage::Vapour]; + + const int aqua_idx = phaseUsage.phase_pos[Opm::PhaseUsage::Aqua]; + const int liquid_idx = phaseUsage.phase_pos[Opm::PhaseUsage::Liquid]; + const int vapour_idx = phaseUsage.phase_pos[Opm::PhaseUsage::Vapour]; + + // WATER + if( aqua_active ) { + addToSimData( simData, "1OVERBW", sd.rq[aqua_idx].b ); + addToSimData( simData, "WAT_DEN", sd.rq[aqua_idx].rho ); + addToSimData( simData, "WAT_VISC", sd.rq[aqua_idx].mu ); + addToSimData( simData, "WATKR", sd.rq[aqua_idx].kr ); + } + + // OIL + if( liquid_active ) { + addToSimData( simData, "1OVERBO", sd.rq[liquid_idx].b ); + addToSimData( simData, "OIL_DEN", sd.rq[liquid_idx].rho ); + addToSimData( simData, "OIL_VISC", sd.rq[liquid_idx].mu ); + addToSimData( simData, "OILKR", sd.rq[liquid_idx].kr ); + } + + // GAS + if( vapour_active ) { + addToSimData( simData, "1OVERBG", sd.rq[vapour_idx].b ); + addToSimData( simData, "GAS_DEN", sd.rq[vapour_idx].rho ); + addToSimData( simData, "GAS_VISC", sd.rq[vapour_idx].mu ); + addToSimData( simData, "GASKR", sd.rq[vapour_idx].kr ); + } + + // RS and RV + addToSimData( simData, "RSSAT", sd.rsSat ); + addToSimData( simData, "RVSAT", sd.rvSat ); + + return std::move( simData ); } + // in case the data is already in a SimulationDataContainer no + // conversion is needed + inline + SimulationDataContainer&& + convertToSimulationDataContainer( SimulationDataContainer&& sd, + const SimulationDataContainer& , + const Opm::PhaseUsage& ) + { + return std::move( sd ); + } /** * Returns the data requested in the restartConfig */ template void getRestartData(data::Solution& output, + SimulationDataContainer&& sd, const Opm::PhaseUsage& phaseUsage, const Model& physicalModel, const RestartConfig& restartConfig, const int reportStepNum, - const bool log) { - - const typename Model::SimulatorData& sd = physicalModel.getSimulatorData(); + const bool log) + { //Get the value of each of the keys for the restart keywords std::map rstKeywords = restartConfig.getRestartKeywords(reportStepNum); for (auto& keyValue : rstKeywords) { keyValue.second = restartConfig.getKeyword(keyValue.first, reportStepNum); } - //Get shorthands for water, oil, gas - const int aqua_active = phaseUsage.phase_used[Opm::PhaseUsage::Aqua]; - const int liquid_active = phaseUsage.phase_used[Opm::PhaseUsage::Liquid]; - const int vapour_active = phaseUsage.phase_used[Opm::PhaseUsage::Vapour]; - - const int aqua_idx = phaseUsage.phase_pos[Opm::PhaseUsage::Aqua]; - const int liquid_idx = phaseUsage.phase_pos[Opm::PhaseUsage::Liquid]; - const int vapour_idx = phaseUsage.phase_pos[Opm::PhaseUsage::Vapour]; + const bool aqua_active = sd.hasCellData("1OVERBW"); + const bool liquid_active = sd.hasCellData("1OVERBO"); + const bool vapour_active = sd.hasCellData("1OVERBG"); + assert( aqua_active == (sd.hasCellData("WAT_DEN") && + sd.hasCellData("WAT_VISC") && + sd.hasCellData("WATKR") + ) + ); + assert( liquid_active == (sd.hasCellData("OIL_DEN") && + sd.hasCellData("OIL_VISC") && + sd.hasCellData("OILKR") + ) + ); + assert( vapour_active == (sd.hasCellData("GAS_DEN") && + sd.hasCellData("GAS_VISC") && + sd.hasCellData("GASKR") + ) + ); /** * Formation volume factors for water, oil, gas @@ -495,21 +575,21 @@ namespace Opm rstKeywords["BW"] = 0; output.insert("1OVERBW", Opm::UnitSystem::measure::water_inverse_formation_volume_factor, - adbToDoubleVector(sd.rq[aqua_idx].b), + std::move( sd.getCellData("1OVERBW") ), data::TargetType::RESTART_AUXILLARY); } if (liquid_active && rstKeywords["BO"] > 0) { rstKeywords["BO"] = 0; output.insert("1OVERBO", Opm::UnitSystem::measure::oil_inverse_formation_volume_factor, - adbToDoubleVector(sd.rq[liquid_idx].b), + std::move( sd.getCellData("1OVERBO") ), data::TargetType::RESTART_AUXILLARY); } if (vapour_active && rstKeywords["BG"] > 0) { rstKeywords["BG"] = 0; output.insert("1OVERBG", Opm::UnitSystem::measure::gas_inverse_formation_volume_factor, - adbToDoubleVector(sd.rq[vapour_idx].b), + std::move( sd.getCellData("1OVERBG") ), data::TargetType::RESTART_AUXILLARY); } @@ -521,19 +601,19 @@ namespace Opm if (aqua_active) { output.insert("WAT_DEN", Opm::UnitSystem::measure::density, - adbToDoubleVector(sd.rq[aqua_idx].rho), + std::move( sd.getCellData("WAT_DEN") ), data::TargetType::RESTART_AUXILLARY); } if (liquid_active) { output.insert("OIL_DEN", Opm::UnitSystem::measure::density, - adbToDoubleVector(sd.rq[liquid_idx].rho), + std::move( sd.getCellData("OIL_DEN") ), data::TargetType::RESTART_AUXILLARY); } if (vapour_active) { output.insert("GAS_DEN", Opm::UnitSystem::measure::density, - adbToDoubleVector(sd.rq[vapour_idx].rho), + std::move( sd.getCellData("GAS_DEN") ), data::TargetType::RESTART_AUXILLARY); } } @@ -546,19 +626,19 @@ namespace Opm if (aqua_active) { output.insert("WAT_VISC", Opm::UnitSystem::measure::viscosity, - adbToDoubleVector(sd.rq[aqua_idx].mu), + std::move( sd.getCellData("WAT_VISC") ), data::TargetType::RESTART_AUXILLARY); } if (liquid_active) { output.insert("OIL_VISC", Opm::UnitSystem::measure::viscosity, - adbToDoubleVector(sd.rq[liquid_idx].mu), + std::move( sd.getCellData("OIL_VISC") ), data::TargetType::RESTART_AUXILLARY); } if (vapour_active) { output.insert("GAS_VISC", Opm::UnitSystem::measure::viscosity, - adbToDoubleVector(sd.rq[vapour_idx].mu), + std::move( sd.getCellData("GAS_VISC") ), data::TargetType::RESTART_AUXILLARY); } } @@ -567,11 +647,12 @@ namespace Opm * Relative permeabilities for water, oil, gas */ if (aqua_active && rstKeywords["KRW"] > 0) { - if (sd.rq[aqua_idx].kr.size() > 0) { + auto& krWater = sd.getCellData("WATKR"); + if (krWater.size() > 0) { rstKeywords["KRW"] = 0; - output.insert("WATKR", + output.insert("WATKR", // WAT_KR ??? Opm::UnitSystem::measure::identity, - adbToDoubleVector(sd.rq[aqua_idx].kr), + std::move( krWater ), data::TargetType::RESTART_AUXILLARY); } else { @@ -583,11 +664,12 @@ namespace Opm } } if (liquid_active && rstKeywords["KRO"] > 0) { - if (sd.rq[liquid_idx].kr.size() > 0) { + auto& krOil = sd.getCellData("OILKR"); + if (krOil.size() > 0) { rstKeywords["KRO"] = 0; output.insert("OILKR", Opm::UnitSystem::measure::identity, - adbToDoubleVector(sd.rq[liquid_idx].kr), + std::move( krOil ), data::TargetType::RESTART_AUXILLARY); } else { @@ -599,11 +681,12 @@ namespace Opm } } if (vapour_active && rstKeywords["KRG"] > 0) { - if (sd.rq[vapour_idx].kr.size() > 0) { + auto& krGas = sd.getCellData("GASKR"); + if (krGas.size() > 0) { rstKeywords["KRG"] = 0; output.insert("GASKR", Opm::UnitSystem::measure::identity, - adbToDoubleVector(sd.rq[vapour_idx].kr), + std::move( krGas ), data::TargetType::RESTART_AUXILLARY); } else { @@ -622,14 +705,14 @@ namespace Opm rstKeywords["RSSAT"] = 0; output.insert("RSSAT", Opm::UnitSystem::measure::gas_oil_ratio, - adbToDoubleVector(sd.rsSat), + std::move( sd.getCellData("RSSAT") ), data::TargetType::RESTART_AUXILLARY); } if (vapour_active && liquid_active && rstKeywords["RVSAT"] > 0) { rstKeywords["RVSAT"] = 0; output.insert("RVSAT", Opm::UnitSystem::measure::oil_gas_ratio, - adbToDoubleVector(sd.rvSat), + std::move( sd.getCellData("RVSAT") ), data::TargetType::RESTART_AUXILLARY); } @@ -637,7 +720,7 @@ namespace Opm /** * Bubble point and dew point pressures */ - if (log && vapour_active && + if (log && vapour_active && liquid_active && rstKeywords["PBPD"] > 0) { rstKeywords["PBPD"] = 0; Opm::OpmLog::warning("Bubble/dew point pressure output unsupported", @@ -683,9 +766,10 @@ namespace Opm const Model& physicalModel, const SummaryConfig& summaryConfig) { - typedef Opm::AutoDiffBlock ADB; + typedef typename Model::FIPDataType FIPDataType; + typedef typename FIPDataType::VectorType VectorType; - const typename Model::SimulatorData& sd = physicalModel.getSimulatorData(); + FIPDataType fd = physicalModel.getFIPData(); //Get shorthands for water, oil, gas const int aqua_active = phaseUsage.phase_used[Opm::PhaseUsage::Aqua]; @@ -699,60 +783,79 @@ namespace Opm if (aqua_active && hasFRBKeyword(summaryConfig, "WIP")) { output.insert("WIP", Opm::UnitSystem::measure::volume, - adbVToDoubleVector(sd.fip[Model::SimulatorData::FIP_AQUA]), + std::move( fd.fip[ FIPDataType::FIP_AQUA ] ), data::TargetType::SUMMARY ); } if (liquid_active) { - const ADB::V& oipl = sd.fip[Model::SimulatorData::FIP_LIQUID]; - const ADB::V& oipg = vapour_active ? sd.fip[Model::SimulatorData::FIP_VAPORIZED_OIL] : ADB::V(); - const ADB::V& oip = vapour_active ? oipl + oipg : oipl; + const VectorType& oipl = fd.fip[FIPDataType::FIP_LIQUID]; + VectorType oip ( oipl ); + const size_t size = oip.size(); + + const VectorType& oipg = vapour_active ? fd.fip[FIPDataType::FIP_VAPORIZED_OIL] : VectorType(size, 0.0); + if( vapour_active ) + { + // oip = oipl + oipg + for( size_t i=0; iisIORank(); - + if( output_ ) { - localCellData = simToSolution(localState, phaseUsage_); // Get "normal" data (SWAT, PRESSURE, ...); - detail::getRestartData( localCellData, phaseUsage_, physicalModel, - restartConfig, reportStepNum, logMessages ); + { + // get all data that need to be included in output from the model + // for flow_legacy and polymer this is a struct holding the data + // while for flow_ebos a SimulationDataContainer is returned + // this is addressed in the above specialized methods + SimulationDataContainer sd = + detail::convertToSimulationDataContainer( physicalModel.getSimulatorData(), localState, phaseUsage_ ); + + localCellData = simToSolution( sd, phaseUsage_); // Get "normal" data (SWAT, PRESSURE, ...); + + detail::getRestartData( localCellData, std::move(sd), phaseUsage_, physicalModel, + restartConfig, reportStepNum, logMessages ); + // sd will be invalid after getRestartData has been called + } detail::getSummaryData( localCellData, phaseUsage_, physicalModel, summaryConfig ); } - + writeTimeStepWithCellProperties(timer, localState, localCellData, localWellState, substep); } } diff --git a/opm/autodiff/SimulatorFullyImplicitBlackoilOutputEbos.cpp b/opm/autodiff/SimulatorFullyImplicitBlackoilOutputEbos.cpp deleted file mode 100644 index 0606f1b36..000000000 --- a/opm/autodiff/SimulatorFullyImplicitBlackoilOutputEbos.cpp +++ /dev/null @@ -1,307 +0,0 @@ -/* - Copyright (c) 2014 SINTEF ICT, Applied Mathematics. - Copyright (c) 2015-2016 IRIS AS - - 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 3 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 . -*/ -#include "config.h" - -#include "SimulatorFullyImplicitBlackoilOutputEbos.hpp" - -#include - -#include -#include -#include -#include -#include -#include -#include -#include - -#include -#include - -#include -#include -#include - -#include - -//For OutputWriterHelper -#include -#include - - -#ifdef HAVE_OPM_GRID -#include -#include -#include -#include -#endif -namespace Opm -{ - namespace detail { - - struct WriterCallEbos : public ThreadHandle :: ObjectInterface - { - BlackoilOutputWriterEbos& writer_; - std::unique_ptr< SimulatorTimerInterface > timer_; - const SimulationDataContainer state_; - const WellStateFullyImplicitBlackoil wellState_; - data::Solution simProps_; - const bool substep_; - - explicit WriterCallEbos( BlackoilOutputWriterEbos& writer, - const SimulatorTimerInterface& timer, - const SimulationDataContainer& state, - const WellStateFullyImplicitBlackoil& wellState, - const data::Solution& simProps, - bool substep ) - : writer_( writer ), - timer_( timer.clone() ), - state_( state ), - wellState_( wellState ), - simProps_( simProps ), - substep_( substep ) - { - } - - // callback to writer's serial writeTimeStep method - void run () - { - // write data - writer_.writeTimeStepSerial( *timer_, state_, wellState_, simProps_, substep_ ); - } - }; - } - - void - BlackoilOutputWriterEbos:: - writeTimeStepWithoutCellProperties( - const SimulatorTimerInterface& timer, - const SimulationDataContainer& localState, - const WellStateFullyImplicitBlackoil& localWellState, - bool substep) - { - data::Solution noCellProperties; - writeTimeStepWithCellProperties(timer, localState, localWellState, noCellProperties, substep); - } - - - - - - void - BlackoilOutputWriterEbos:: - writeTimeStepWithCellProperties( - const SimulatorTimerInterface& timer, - const SimulationDataContainer& localState, - const WellStateFullyImplicitBlackoil& localWellState, - const data::Solution& sol, - bool substep) - { - bool isIORank = output_ ; - if( parallelOutput_ && parallelOutput_->isParallel() ) - { - // If this is not the initial write and no substep, then the well - // state used in the computation is actually the one of the last - // step. We need that well state for the gathering. Otherwise - // It an exception with a message like "global state does not - // contain well ..." might be thrown. - int wellStateStepNumber = ( ! substep && timer.reportStepNum() > 0) ? - (timer.reportStepNum() - 1) : timer.reportStepNum(); - // collect all solutions to I/O rank - isIORank = parallelOutput_->collectToIORank( localState, localWellState, sol, wellStateStepNumber ); - } - - const SimulationDataContainer& state = (parallelOutput_ && parallelOutput_->isParallel() ) ? parallelOutput_->globalReservoirState() : localState; - const WellStateFullyImplicitBlackoil& wellState = (parallelOutput_ && parallelOutput_->isParallel() ) ? parallelOutput_->globalWellState() : localWellState; - - // serial output is only done on I/O rank - int err = 0; - std::string emsg; - if( isIORank ) - { - if( asyncOutput_ ) { - // dispatch the write call to the extra thread - asyncOutput_->dispatch( detail::WriterCallEbos( *this, timer, state, wellState, sol, substep ) ); - } - else { - // just write the data to disk - try { - writeTimeStepSerial( timer, state, wellState, sol, substep ); - } catch (std::runtime_error& msg) { - err = 1; - emsg = msg.what(); - } - } - } - - if (!asyncOutput_) { -#if HAVE_MPI - MPI_Bcast(&err, 1, MPI_INT, 0, MPI_COMM_WORLD); -#endif - if (err) { - if (isIORank) { - throw std::runtime_error(emsg); - } else { - throw std::runtime_error("I/O process encountered problems."); - } - } - } - } - - - - void - BlackoilOutputWriterEbos:: - writeTimeStepSerial(const SimulatorTimerInterface& timer, - const SimulationDataContainer& state, - const WellStateFullyImplicitBlackoil& wellState, - const data::Solution& sol, - bool substep) - { - // ECL output - if (output()) - { - const auto& initConfig = eclipseState_.getInitConfig(); - if (initConfig.restartRequested() && ((initConfig.getRestartStep()) == (timer.currentStepNum()))) { - std::cout << "Skipping restart write in start of step " << timer.currentStepNum() << std::endl; - } else { - data::Solution combined_sol = simToSolution(state, phaseUsage_); // Get "normal" data (SWAT, PRESSURE, ...) - combined_sol.insert(sol.begin(), sol.end()); // ... insert "extra" data (KR, VISC, ...) - eclIO_.writeTimeStep(timer.reportStepNum(), - substep, - timer.simulationTimeElapsed(), - combined_sol, - wellState.report(phaseUsage_)); - } - } - - // write backup file - if( backupfile_.is_open() ) - { - int reportStep = timer.reportStepNum(); - int currentTimeStep = timer.currentStepNum(); - if( (reportStep == currentTimeStep || // true for SimulatorTimer - currentTimeStep == 0 || // true for AdaptiveSimulatorTimer at reportStep - timer.done() ) // true for AdaptiveSimulatorTimer at reportStep - && lastBackupReportStep_ != reportStep ) // only backup report step once - { - // store report step - lastBackupReportStep_ = reportStep; - // write resport step number - backupfile_.write( (const char *) &reportStep, sizeof(int) ); - - try { - backupfile_ << state; - - const WellStateFullyImplicitBlackoil& boWellState = static_cast< const WellStateFullyImplicitBlackoil& > (wellState); - backupfile_ << boWellState; - } - catch ( const std::bad_cast& e ) - { - } - - backupfile_ << std::flush; - } - } // end backup - } - - void - BlackoilOutputWriterEbos:: - restore(SimulatorTimerInterface& timer, - BlackoilState& state, - WellStateFullyImplicitBlackoilDense& wellState, - const std::string& filename, - const int desiredResportStep ) - { - std::ifstream restorefile( filename.c_str() ); - if( restorefile ) - { - std::cout << "============================================================================"<. -*/ -#ifndef OPM_SIMULATORFULLYIMPLICITBLACKOILOUTPUTEBOS_HEADER_INCLUDED -#define OPM_SIMULATORFULLYIMPLICITBLACKOILOUTPUTEBOS_HEADER_INCLUDED -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#include - -#include -#include -#include - -#include -#include -#include - -#include -#include - - -#include -#include -#include -#include -#include -#include - -#include - -#ifdef HAVE_OPM_GRID -#include -#endif -namespace Opm -{ - class BlackoilState; - - /** \brief Wrapper class for VTK, Matlab, and ECL output. */ - class BlackoilOutputWriterEbos - { - - public: - // constructor creating different sub writers - template - BlackoilOutputWriterEbos(const Grid& grid, - const parameter::ParameterGroup& param, - const EclipseState& eclipseState, - EclipseIO& eclIO, - const Opm::PhaseUsage &phaseUsage); - - /*! - * \brief Write a blackoil reservoir state to disk for later inspection with - * visualization tools like ResInsight. This function will extract the - * requested output cell properties specified by the RPTRST keyword - * and write these to file. - */ - template - void writeTimeStep(const SimulatorTimerInterface& timer, - const SimulationDataContainer& reservoirState, - const Opm::WellStateFullyImplicitBlackoil& wellState, - const Model& physicalModel, - bool substep = false); - - - /*! - * \brief Write a blackoil reservoir state to disk for later inspection with - * visualization tools like ResInsight. This function will write all - * CellData in simProps to the file as well. - */ - void writeTimeStepWithCellProperties( - const SimulatorTimerInterface& timer, - const SimulationDataContainer& reservoirState, - const Opm::WellStateFullyImplicitBlackoil& wellState, - const data::Solution& sol, - bool substep = false); - - /*! - * \brief Write a blackoil reservoir state to disk for later inspection with - * visualization tools like ResInsight. This function will not write - * any cell properties (e.g., those requested by RPTRST keyword) - */ - void writeTimeStepWithoutCellProperties( - const SimulatorTimerInterface& timer, - const SimulationDataContainer& reservoirState, - const Opm::WellStateFullyImplicitBlackoil& wellState, - bool substep = false); - - /*! - * \brief Write a blackoil reservoir state to disk for later inspection withS - * visualization tools like ResInsight. This is the function which does - * the actual write to file. - */ - void writeTimeStepSerial(const SimulatorTimerInterface& timer, - const SimulationDataContainer& reservoirState, - const Opm::WellStateFullyImplicitBlackoil& wellState, - const data::Solution& simProps, - bool substep); - - /** \brief return output directory */ - const std::string& outputDirectory() const { return outputDir_; } - - /** \brief return true if output is enabled */ - bool output () const { return output_; } - - /** \brief Whether this process does write to disk */ - bool isIORank () const - { - return parallelOutput_->isIORank(); - } - - void restore(SimulatorTimerInterface& timer, - BlackoilState& state, - WellStateFullyImplicitBlackoilDense& wellState, - const std::string& filename, - const int desiredReportStep); - - - template - void initFromRestartFile(const PhaseUsage& phaseusage, - const Grid& grid, - SimulationDataContainer& simulatorstate, - WellStateFullyImplicitBlackoilDense& wellstate); - - bool isRestart() const; - - bool requireFIPNUM() const; - - protected: - const bool output_; - std::unique_ptr< ParallelDebugOutputInterface > parallelOutput_; - - // Parameters for output. - const std::string outputDir_; - const int output_interval_; - - int lastBackupReportStep_; - - std::ofstream backupfile_; - Opm::PhaseUsage phaseUsage_; - EclipseIO& eclIO_; - const EclipseState& eclipseState_; - - std::unique_ptr< ThreadHandle > asyncOutput_; - }; - - - ////////////////////////////////////////////////////////////// - // - // Implementation - // - ////////////////////////////////////////////////////////////// - template - inline - BlackoilOutputWriterEbos:: - BlackoilOutputWriterEbos(const Grid& grid, - const parameter::ParameterGroup& param, - const Opm::EclipseState& eclipseState, - EclipseIO& eclIO, - const Opm::PhaseUsage &phaseUsage) - : output_( param.getDefault("output", true) ), - parallelOutput_( output_ ? new ParallelDebugOutput< Grid >( grid, eclipseState, phaseUsage.num_phases, phaseUsage ) : 0 ), - outputDir_( output_ ? param.getDefault("output_dir", std::string("output")) : "." ), - output_interval_( output_ ? param.getDefault("output_interval", 1): 0 ), - lastBackupReportStep_( -1 ), - phaseUsage_( phaseUsage ), - eclIO_(eclIO), - eclipseState_(eclipseState), - asyncOutput_() - { - // For output. - if (output_ && parallelOutput_->isIORank() ) { - // Ensure that output dir exists - boost::filesystem::path fpath(outputDir_); - try { - create_directories(fpath); - } - catch (...) { - OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); - } - - // create output thread if enabled and rank is I/O rank - // async output is enabled by default if pthread are enabled -#if HAVE_PTHREAD - const bool asyncOutputDefault = false; -#else - const bool asyncOutputDefault = false; -#endif - if( param.getDefault("async_output", asyncOutputDefault ) ) - { -#if HAVE_PTHREAD - asyncOutput_.reset( new ThreadHandle() ); -#else - OPM_THROW(std::runtime_error,"Pthreads were not found, cannot enable async_output"); -#endif - } - - std::string backupfilename = param.getDefault("backupfile", std::string("") ); - if( ! backupfilename.empty() ) - { - backupfile_.open( backupfilename.c_str() ); - } - } - } - - - template - inline void - BlackoilOutputWriterEbos:: - initFromRestartFile( const PhaseUsage& phaseusage, - const Grid& grid, - SimulationDataContainer& simulatorstate, - WellStateFullyImplicitBlackoilDense& wellstate) - { - std::map solution_keys {{"PRESSURE" , UnitSystem::measure::pressure}, - {"SWAT" , UnitSystem::measure::identity}, - {"SGAS" , UnitSystem::measure::identity}, - {"TEMP" , UnitSystem::measure::temperature}, - {"RS" , UnitSystem::measure::gas_oil_ratio}, - {"RV" , UnitSystem::measure::oil_gas_ratio}}; - - // gives a dummy dynamic_list_econ_limited - DynamicListEconLimited dummy_list_econ_limited; - WellsManager wellsmanager(eclipseState_, - eclipseState_.getInitConfig().getRestartStep(), - Opm::UgGridHelpers::numCells(grid), - Opm::UgGridHelpers::globalCell(grid), - Opm::UgGridHelpers::cartDims(grid), - Opm::UgGridHelpers::dimensions(grid), - Opm::UgGridHelpers::cell2Faces(grid), - Opm::UgGridHelpers::beginFaceCentroids(grid), - dummy_list_econ_limited - // We need to pass the optionaly arguments - // as we get the following error otherwise - // with c++ (Debian 4.9.2-10) 4.9.2 and -std=c++11 - // converting to ‘const std::unordered_set >’ from initializer list would use explicit constructo - , false, - std::vector(), - std::unordered_set()); - - const Wells* wells = wellsmanager.c_wells(); - wellstate.resize(wells, simulatorstate, phaseusage ); //Resize for restart step - auto state = eclIO_.loadRestart(solution_keys); - solutionToSim( state.first, phaseusage, simulatorstate ); - wellsToState( state.second, phaseusage, wellstate ); - } - - - - - - namespace detail { - template - void getOutputDataEbos(data::Solution& output, - const Opm::PhaseUsage& phaseUsage, - const Model& model, - const RestartConfig& restartConfig, - const int reportStepNum, - const bool log) - { - typedef typename Model::FluidSystem FluidSystem; - - //Get the value of each of the keys - std::map outKeywords = restartConfig.getRestartKeywords(reportStepNum); - for (auto& keyValue : outKeywords) { - keyValue.second = restartConfig.getKeyword(keyValue.first, reportStepNum); - } - - const auto& ebosModel = model.ebosSimulator().model(); - - //Get shorthands for water, oil, gas - const int aqua_active = phaseUsage.phase_used[Opm::PhaseUsage::Aqua]; - const int liquid_active = phaseUsage.phase_used[Opm::PhaseUsage::Liquid]; - const int vapour_active = phaseUsage.phase_used[Opm::PhaseUsage::Vapour]; - - // extract everything which can possibly be written to disk - int numCells = ebosModel.numGridDof(); - - std::vector pressureOil(numCells); - std::vector temperature(numCells); - - std::vector satWater(numCells); - std::vector satGas(numCells); - - std::vector bWater(numCells); - std::vector bOil(numCells); - std::vector bGas(numCells); - - std::vector rhoWater(numCells); - std::vector rhoOil(numCells); - std::vector rhoGas(numCells); - - std::vector muWater(numCells); - std::vector muOil(numCells); - std::vector muGas(numCells); - - std::vector krWater(numCells); - std::vector krOil(numCells); - std::vector krGas(numCells); - - std::vector Rs(numCells); - std::vector Rv(numCells); - std::vector RsSat(numCells); - std::vector RvSat(numCells); - - for (int cellIdx = 0; cellIdx < numCells; ++cellIdx) { - const auto& intQuants = *ebosModel.cachedIntensiveQuantities(cellIdx, /*timeIdx=*/0); - const auto& fs = intQuants.fluidState(); - - pressureOil[cellIdx] = fs.pressure(FluidSystem::oilPhaseIdx).value(); - - temperature[cellIdx] = fs.temperature(FluidSystem::oilPhaseIdx).value(); - - if (aqua_active) { - satWater[cellIdx] = fs.saturation(FluidSystem::waterPhaseIdx).value(); - bWater[cellIdx] = fs.invB(FluidSystem::waterPhaseIdx).value(); - rhoWater[cellIdx] = fs.density(FluidSystem::waterPhaseIdx).value(); - muWater[cellIdx] = fs.viscosity(FluidSystem::waterPhaseIdx).value(); - krWater[cellIdx] = intQuants.relativePermeability(FluidSystem::waterPhaseIdx).value(); - } - if (vapour_active) { - satGas[cellIdx] = fs.saturation(FluidSystem::gasPhaseIdx).value(); - bGas[cellIdx] = fs.invB(FluidSystem::gasPhaseIdx).value(); - rhoGas[cellIdx] = fs.density(FluidSystem::gasPhaseIdx).value(); - muGas[cellIdx] = fs.viscosity(FluidSystem::gasPhaseIdx).value(); - krGas[cellIdx] = intQuants.relativePermeability(FluidSystem::gasPhaseIdx).value(); - Rs[cellIdx] = fs.Rs().value(); - Rv[cellIdx] = fs.Rv().value(); - RsSat[cellIdx] = FluidSystem::saturatedDissolutionFactor(fs, - FluidSystem::oilPhaseIdx, - intQuants.pvtRegionIndex(), - /*maxOilSaturation=*/1.0).value(); - RvSat[cellIdx] = FluidSystem::saturatedDissolutionFactor(fs, - FluidSystem::gasPhaseIdx, - intQuants.pvtRegionIndex(), - /*maxOilSaturation=*/1.0).value(); - } - - // oil is always active - bOil[cellIdx] = fs.invB(FluidSystem::oilPhaseIdx).value(); - rhoOil[cellIdx] = fs.density(FluidSystem::oilPhaseIdx).value(); - muOil[cellIdx] = fs.viscosity(FluidSystem::oilPhaseIdx).value(); - krOil[cellIdx] = intQuants.relativePermeability(FluidSystem::oilPhaseIdx).value(); - } - - /** - * Oil Pressures - */ - outKeywords["PRESSURE"] = 0; - output.insert("PRESSURE", - UnitSystem::measure::pressure, - std::move(pressureOil), - data::TargetType::RESTART_SOLUTION); - - /** - * Temperatures - */ - outKeywords["TEMP"] = 0; - output.insert("TEMP", - UnitSystem::measure::temperature, - std::move(temperature), - data::TargetType::RESTART_SOLUTION); - - /** - * Water and gas saturation. - */ - outKeywords["SWAT"] = 0; - outKeywords["SGAS"] = 0; - output.insert("SWAT", - UnitSystem::measure::identity, - std::move(satWater), - data::TargetType::RESTART_SOLUTION); - output.insert("SGAS", - UnitSystem::measure::identity, - std::move(satGas), - data::TargetType::RESTART_SOLUTION); - - /** - * the dissolution factors - */ - outKeywords["RS"] = 0; - outKeywords["RV"] = 0; - output.insert("RS", - UnitSystem::measure::gas_oil_ratio, - std::move(Rs), - data::TargetType::RESTART_SOLUTION); - output.insert("RV", - UnitSystem::measure::oil_gas_ratio, - std::move(Rv), - data::TargetType::RESTART_SOLUTION); - - /** - * Formation volume factors for water, oil, gas - */ - if (aqua_active && outKeywords["BW"] > 0) { - outKeywords["BW"] = 0; - output.insert("BW", - Opm::UnitSystem::measure::water_inverse_formation_volume_factor, - std::move(bWater), - data::TargetType::RESTART_AUXILLARY); - - } - if (liquid_active && outKeywords["BO"] > 0) { - outKeywords["BO"] = 0; - output.insert("BO", - Opm::UnitSystem::measure::oil_inverse_formation_volume_factor, - std::move(bOil), - data::TargetType::RESTART_AUXILLARY); - } - if (vapour_active && outKeywords["BG"] > 0) { - outKeywords["BG"] = 0; - output.insert("BG", - Opm::UnitSystem::measure::gas_inverse_formation_volume_factor, - std::move(bGas), - data::TargetType::RESTART_AUXILLARY); - } - - /** - * Densities for water, oil gas - */ - if (outKeywords["DEN"] > 0) { - outKeywords["DEN"] = 0; - if (aqua_active) { - output.insert("WAT_DEN", - Opm::UnitSystem::measure::density, - std::move(rhoWater), - data::TargetType::RESTART_AUXILLARY); - } - if (liquid_active) { - output.insert("OIL_DEN", - Opm::UnitSystem::measure::density, - std::move(rhoOil), - data::TargetType::RESTART_AUXILLARY); - } - if (vapour_active) { - output.insert("GAS_DEN", - Opm::UnitSystem::measure::density, - std::move(rhoGas), - data::TargetType::RESTART_AUXILLARY); - } - } - - /** - * Viscosities for water, oil gas - */ - if (outKeywords["VISC"] > 0) { - outKeywords["VISC"] = 0; - if (aqua_active) { - output.insert("WAT_VISC", - Opm::UnitSystem::measure::viscosity, - std::move(muWater), - data::TargetType::RESTART_AUXILLARY); - } - if (liquid_active) { - output.insert("OIL_VISC", - Opm::UnitSystem::measure::viscosity, - std::move(muOil), - data::TargetType::RESTART_AUXILLARY); - } - if (vapour_active) { - output.insert("GAS_VISC", - Opm::UnitSystem::measure::viscosity, - std::move(muGas), - data::TargetType::RESTART_AUXILLARY); - } - } - - /** - * Relative permeabilities for water, oil, gas - */ - if (aqua_active && outKeywords["KRW"] > 0) { - outKeywords["KRW"] = 0; - output.insert("WATKR", - Opm::UnitSystem::measure::identity, - std::move(krWater), - data::TargetType::RESTART_AUXILLARY); - } - if (liquid_active && outKeywords["KRO"] > 0) { - outKeywords["KRO"] = 0; - output.insert("OILKR", - Opm::UnitSystem::measure::identity, - std::move(krOil), - data::TargetType::RESTART_AUXILLARY); - } - if (vapour_active && outKeywords["KRG"] > 0) { - outKeywords["KRG"] = 0; - output.insert("GASKR", - Opm::UnitSystem::measure::identity, - std::move(krGas), - data::TargetType::RESTART_AUXILLARY); - } - - /** - * Vaporized and dissolved gas/oil ratio - */ - if (vapour_active && liquid_active && outKeywords["RSSAT"] > 0) { - outKeywords["RSSAT"] = 0; - output.insert("RSSAT", - Opm::UnitSystem::measure::gas_oil_ratio, - std::move(RsSat), - data::TargetType::RESTART_AUXILLARY); - } - if (vapour_active && liquid_active && outKeywords["RVSAT"] > 0) { - outKeywords["RVSAT"] = 0; - output.insert("RVSAT", - Opm::UnitSystem::measure::oil_gas_ratio, - std::move(RvSat), - data::TargetType::RESTART_AUXILLARY); - } - - - /** - * Bubble point and dew point pressures - */ - if (log && vapour_active && - liquid_active && outKeywords["PBPD"] > 0) { - Opm::OpmLog::warning("Bubble/dew point pressure output unsupported", - "Writing bubble points and dew points (PBPD) to file is unsupported, " - "as the simulator does not use these internally."); - } - - //Warn for any unhandled keyword - if (log) { - for (auto& keyValue : outKeywords) { - if (keyValue.second > 0) { - std::string logstring = "Keyword '"; - logstring.append(keyValue.first); - logstring.append("' is unhandled for output to file."); - Opm::OpmLog::warning("Unhandled output keyword", logstring); - } - } - } - - } - - /** - * Checks if the summaryConfig has a keyword with the standardized field, region, or block prefixes. - */ - inline bool hasFRBKeyword(const SummaryConfig& summaryConfig, const std::string keyword) { - std::string field_kw = "F" + keyword; - std::string region_kw = "R" + keyword; - std::string block_kw = "B" + keyword; - return summaryConfig.hasKeyword(field_kw) - || summaryConfig.hasKeyword(region_kw) - || summaryConfig.hasKeyword(block_kw); - } - - - /** - * Returns the data as asked for in the summaryConfig - */ - template - void getSummaryData(data::Solution& output, - const Opm::PhaseUsage& phaseUsage, - const Model& physicalModel, - const SummaryConfig& summaryConfig) { - - const typename Model::FIPData& fip = physicalModel.getFIPData(); - - //Get shorthands for water, oil, gas - const int aqua_active = phaseUsage.phase_used[Opm::PhaseUsage::Aqua]; - const int liquid_active = phaseUsage.phase_used[Opm::PhaseUsage::Liquid]; - const int vapour_active = phaseUsage.phase_used[Opm::PhaseUsage::Vapour]; - - /** - * Now process all of the summary config files - */ - // Water in place - if (aqua_active && hasFRBKeyword(summaryConfig, "WIP")) { - output.insert("WIP", - Opm::UnitSystem::measure::volume, - fip.fip[Model::FIPData::FIP_AQUA], - data::TargetType::SUMMARY ); - } - if (liquid_active) { - const std::vector& oipl = fip.fip[Model::FIPData::FIP_LIQUID]; - const int size = oipl.size(); - - const std::vector& oipg = vapour_active ? fip.fip[Model::FIPData::FIP_VAPORIZED_OIL] : std::vector(size,0.0); - std::vector oip = oipl; - if (vapour_active) { - std::transform(oip.begin(), oip.end(), oipg.begin(), oip.begin(), std::plus()); - } - - //Oil in place (liquid phase only) - if (hasFRBKeyword(summaryConfig, "OIPL")) { - output.insert("OIPL", - Opm::UnitSystem::measure::volume, - oipl, - data::TargetType::SUMMARY ); - } - - //Oil in place (gas phase only) - if (hasFRBKeyword(summaryConfig, "OIPG")) { - output.insert("OIPG", - Opm::UnitSystem::measure::volume, - oipg, - data::TargetType::SUMMARY ); - } - // Oil in place (in liquid and gas phases) - if (hasFRBKeyword(summaryConfig, "OIP")) { - output.insert("OIP", - Opm::UnitSystem::measure::volume, - oip, - data::TargetType::SUMMARY ); - } - - - - } - if (vapour_active) { - const std::vector& gipg = fip.fip[Model::FIPData::FIP_VAPOUR]; - const int size = gipg.size(); - - const std::vector& gipl= liquid_active ? fip.fip[Model::FIPData::FIP_DISSOLVED_GAS] : std::vector(size,0.0); - std::vector gip = gipg; - if (liquid_active) { - std::transform(gip.begin(), gip.end(), gipl.begin(), gip.begin(), std::plus()); - } - - // Gas in place (gas phase only) - if (hasFRBKeyword(summaryConfig, "GIPG")) { - output.insert("GIPG", - Opm::UnitSystem::measure::volume, - gipg, - data::TargetType::SUMMARY ); - } - // Gas in place (liquid phase only) - if (hasFRBKeyword(summaryConfig, "GIPL")) { - output.insert("GIPL", - Opm::UnitSystem::measure::volume, - gipl, - data::TargetType::SUMMARY ); - } - // Gas in place (in both liquid and gas phases) - if (hasFRBKeyword(summaryConfig, "GIP")) { - output.insert("GIP", - Opm::UnitSystem::measure::volume, - gip, - data::TargetType::SUMMARY ); - } - } - // Cell pore volume in reservoir conditions - if (hasFRBKeyword(summaryConfig, "RPV")) { - output.insert("RPV", - Opm::UnitSystem::measure::volume, - fip.fip[Model::FIPData::FIP_PV], - data::TargetType::SUMMARY ); - } - // Pressure averaged value (hydrocarbon pore volume weighted) - if (summaryConfig.hasKeyword("FPRH") || summaryConfig.hasKeyword("RPRH")) { - output.insert("PRH", - Opm::UnitSystem::measure::pressure, - fip.fip[Model::FIPData::FIP_WEIGHTED_PRESSURE], - data::TargetType::SUMMARY ); - } - } - - } - - - - - template - inline void - BlackoilOutputWriterEbos:: - writeTimeStep(const SimulatorTimerInterface& timer, - const SimulationDataContainer& localState, - const WellStateFullyImplicitBlackoil& localWellState, - const Model& physicalModel, - bool substep) - { - data::Solution cellData{}; - const RestartConfig& restartConfig = eclipseState_.getRestartConfig(); - const SummaryConfig& summaryConfig = eclipseState_.getSummaryConfig(); - const int reportStepNum = timer.reportStepNum(); - bool logMessages = output_ && parallelOutput_->isIORank(); - - if( output_ && !parallelOutput_->isParallel() ) - { - - detail::getOutputDataEbos( cellData, phaseUsage_, physicalModel, - restartConfig, reportStepNum, logMessages ); - detail::getSummaryData( cellData, phaseUsage_, physicalModel, summaryConfig ); - } - else - { - if ( logMessages ) - { - std::map rstKeywords = restartConfig.getRestartKeywords(reportStepNum); - std::vector keywords = - { "WIP", "OIPL", "OIPG", "OIP", "GIPG", "GIPL", "GIP", - "RPV", "FRPH", "RPRH"}; - - std::ostringstream str; - str << "Output of restart/summary config not supported in parallel. Requested keywords were "; - std::size_t no_kw = 0; - - auto func = [&] (const char* kw) - { - if ( detail::hasFRBKeyword(summaryConfig, kw) ) - { - str << kw << " "; - ++ no_kw; - } - }; - - std::for_each(keywords.begin(), keywords.end(), func); - - for (auto& keyValue : rstKeywords) - { - str << keyValue.first << " "; - ++ no_kw; - } - - if ( no_kw ) - { - Opm::OpmLog::warning("Unhandled ouput request", str.str()); - } - } - } - - writeTimeStepWithCellProperties(timer, localState, localWellState, cellData, substep); - } -} -#endif diff --git a/opm/autodiff/WellStateFullyImplicitBlackoil.hpp b/opm/autodiff/WellStateFullyImplicitBlackoil.hpp index a22a48750..9daf543ad 100644 --- a/opm/autodiff/WellStateFullyImplicitBlackoil.hpp +++ b/opm/autodiff/WellStateFullyImplicitBlackoil.hpp @@ -20,7 +20,7 @@ #ifndef OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED #define OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED - +#include #include #include #include @@ -176,11 +176,15 @@ namespace Opm } template - void resize(const Wells* wells, const State& state) { + void resize(const Wells* wells, const State& state ) { const WellStateFullyImplicitBlackoil dummy_state{}; // Init with an empty previous state only resizes init(wells, state, dummy_state) ; } + template + void resize(const Wells* wells, const State& state, const PhaseUsage& ) { + resize( wells, state ); + } /// One rate per phase and well connection. std::vector& perfPhaseRates() { return perfphaserates_; } diff --git a/opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.hpp b/opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.hpp index 72b706040..04a99f815 100644 --- a/opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.hpp +++ b/opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.hpp @@ -32,6 +32,7 @@ #include #include #include +#include struct UnstructuredGrid; struct Wells; @@ -76,7 +77,7 @@ namespace Opm { std::vector ads; // Adsorption term. }; - struct SimulatorData { + struct SimulatorData : public Opm::FIPDataEnums { SimulatorData(int num_phases) : rq(num_phases) , rsSat(ADB::null()) @@ -85,21 +86,18 @@ namespace Opm { { } - enum FipId { - FIP_AQUA = Opm::Water, - FIP_LIQUID = Opm::Oil, - FIP_VAPOUR = Opm::Gas, - FIP_DISSOLVED_GAS = 3, - FIP_VAPORIZED_OIL = 4, - FIP_PV = 5, //< Pore volume - FIP_WEIGHTED_PRESSURE = 6 - }; + using Opm::FIPDataEnums::FipId; + using Opm::FIPDataEnums::fipValues; + std::vector rq; ADB rsSat; ADB rvSat; - std::array fip; + std::array fip; }; + typedef Opm::FIPData FIPDataType; + + /// Construct a solver. It will retain references to the /// arguments of this functions, and they are expected to /// remain in scope for the lifetime of the solver. @@ -153,12 +151,17 @@ namespace Opm { return sd_; } + /// Return reservoir simulation data (for output functionality) + FIPDataType getFIPData() const { + return FIPDataType( sd_.fip ); + } + /// Compute fluid in place. /// \param[in] ReservoirState /// \param[in] WellState /// \param[in] FIPNUM for active cells not global cells. - /// \return fluid in place, number of fip regions, each region contains 5 values which are liquid, vapour, water, free gas and dissolved gas. - std::vector > + /// \return fluid in place, number of fip regions, each region contains 5 values which are liquid, vapour, water, free gas and dissolved gas. + std::vector > computeFluidInPlace(const PolymerBlackoilState& x, const std::vector& fipnum); @@ -224,7 +227,7 @@ namespace Opm { void assemble(const double dt, const PolymerBlackoilState& x, - const WellStateFullyImplicitBlackoil& xw, + const WellStateFullyImplicitBlackoil& xw, const std::vector& polymer_inflow); V solveJacobianSystem() const; @@ -298,10 +301,10 @@ namespace Opm { ADB transMult(const ADB& p) const; - + const std::vector phaseCondition() const { return phaseCondition_; } - + void classifyCondition(const PolymerBlackoilState& state); };