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); };