From 9febba59a415340eb8b74a7975060c5b2589c8db Mon Sep 17 00:00:00 2001 From: "Kjell W. Kongsvik" Date: Thu, 17 Mar 2016 13:30:48 +0100 Subject: [PATCH 1/4] Commented out usage of OutputWriter in tutorials --- tutorials/tutorial1.cpp | 9 ++++++--- tutorials/tutorial2.cpp | 18 ++++++++++-------- tutorials/tutorial3.cpp | 14 ++++++++------ tutorials/tutorial4.cpp | 12 +++++++----- 4 files changed, 31 insertions(+), 22 deletions(-) diff --git a/tutorials/tutorial1.cpp b/tutorials/tutorial1.cpp index 428e365c..e4a57243 100644 --- a/tutorials/tutorial1.cpp +++ b/tutorials/tutorial1.cpp @@ -37,7 +37,8 @@ #include #include -#include +// 17.03.2016 Temporarily removed while moving functionality to opm-output +//#include #include #include #include @@ -94,14 +95,16 @@ try /// Opm::writeVtkData() together with the grid /// \snippet tutorial1.cpp data map /// \internal [data map] - Opm::DataMap dm; +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// Opm::DataMap dm; /// \internal [data map] /// \endinternal /// \page tutorial1 /// Call Opm::writeVtkData() to write the output file. /// \snippet tutorial1.cpp write vtk /// \internal [write vtk] - Opm::writeVtkData(*grid.c_grid(), dm, vtkfile); +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// Opm::writeVtkData(*grid.c_grid(), dm, vtkfile); /// \internal [write vtk] /// \endinternal } diff --git a/tutorials/tutorial2.cpp b/tutorials/tutorial2.cpp index d8e1e8c4..445e76a6 100644 --- a/tutorials/tutorial2.cpp +++ b/tutorials/tutorial2.cpp @@ -34,7 +34,8 @@ #include #include -#include +// 17.03.2016 Temporarily removed while moving functionality to opm-output +//#include #include #include #include @@ -186,13 +187,14 @@ try /// (pressure) or a vector per cell (cell_velocity). /// \snippet tutorial2.cpp write output /// \internal [write output] - std::ofstream vtkfile("tutorial2.vtu"); - Opm::DataMap dm; - dm["pressure"] = &state.pressure(); - std::vector cell_velocity; - Opm::estimateCellVelocity(*grid.c_grid(), state.faceflux(), cell_velocity); - dm["velocity"] = &cell_velocity; - Opm::writeVtkData(*grid.c_grid(), dm, vtkfile); +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// std::ofstream vtkfile("tutorial2.vtu"); +// Opm::DataMap dm; +// dm["pressure"] = &state.pressure(); +// std::vector cell_velocity; +// Opm::estimateCellVelocity(*grid.c_grid(), state.faceflux(), cell_velocity); +// dm["velocity"] = &cell_velocity; +// Opm::writeVtkData(*grid.c_grid(), dm, vtkfile); /// \internal [write output] /// \endinternal } diff --git a/tutorials/tutorial3.cpp b/tutorials/tutorial3.cpp index 976839a5..398dbc44 100644 --- a/tutorials/tutorial3.cpp +++ b/tutorials/tutorial3.cpp @@ -29,7 +29,8 @@ #include #include #include -#include +// 17.03.2016 Temporarily removed while moving functionality to opm-output +//#include #include #include #include @@ -313,11 +314,12 @@ try /// \internal [write output] vtkfilename.str(""); vtkfilename << "tutorial3-" << std::setw(3) << std::setfill('0') << i << ".vtu"; - std::ofstream vtkfile(vtkfilename.str().c_str()); - Opm::DataMap dm; - dm["saturation"] = &state.saturation(); - dm["pressure"] = &state.pressure(); - Opm::writeVtkData(grid, dm, vtkfile); +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// std::ofstream vtkfile(vtkfilename.str().c_str()); +// Opm::DataMap dm; +// dm["saturation"] = &state.saturation(); +// dm["pressure"] = &state.pressure(); +// Opm::writeVtkData(grid, dm, vtkfile); } } catch (const std::exception &e) { diff --git a/tutorials/tutorial4.cpp b/tutorials/tutorial4.cpp index f67c55fa..eee0de49 100644 --- a/tutorials/tutorial4.cpp +++ b/tutorials/tutorial4.cpp @@ -29,7 +29,8 @@ #include #include #include -#include +// 17.03.2016 Temporarily removed while moving functionality to opm-output +//#include #include #include #include @@ -434,10 +435,11 @@ try vtkfilename.str(""); vtkfilename << "tutorial4-" << std::setw(3) << std::setfill('0') << i << ".vtu"; std::ofstream vtkfile(vtkfilename.str().c_str()); - Opm::DataMap dm; - dm["saturation"] = &state.saturation(); - dm["pressure"] = &state.pressure(); - Opm::writeVtkData(grid, dm, vtkfile); +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// Opm::DataMap dm; +// dm["saturation"] = &state.saturation(); +// dm["pressure"] = &state.pressure(); +// Opm::writeVtkData(grid, dm, vtkfile); } destroy_wells(wells); From 03d5a5ba25ee9f149e8a63664559b0a1f67a330a Mon Sep 17 00:00:00 2001 From: "Kjell W. Kongsvik" Date: Thu, 17 Mar 2016 14:12:34 +0100 Subject: [PATCH 2/4] Commented out usage of OutputWriter in simulator --- .../SimulatorCompressibleTwophase.cpp | 137 +++++++++--------- .../simulator/SimulatorIncompTwophase.cpp | 136 ++++++++--------- opm/core/simulator/SimulatorOutput.cpp | 12 +- opm/core/simulator/SimulatorOutput.hpp | 6 +- 4 files changed, 152 insertions(+), 139 deletions(-) diff --git a/opm/core/simulator/SimulatorCompressibleTwophase.cpp b/opm/core/simulator/SimulatorCompressibleTwophase.cpp index c0875a2e..c223b665 100644 --- a/opm/core/simulator/SimulatorCompressibleTwophase.cpp +++ b/opm/core/simulator/SimulatorCompressibleTwophase.cpp @@ -34,7 +34,8 @@ #include #include #include -#include +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// #include #include #include @@ -138,70 +139,72 @@ namespace Opm - static void outputStateVtk(const UnstructuredGrid& grid, - const Opm::BlackoilState& state, - const int step, - const std::string& output_dir) - { - // Write data in VTK format. - std::ostringstream vtkfilename; - vtkfilename << output_dir << "/vtk_files"; - boost::filesystem::path fpath(vtkfilename.str()); - try { - create_directories(fpath); - } - catch (...) { - OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); - } - vtkfilename << "/output-" << std::setw(3) << std::setfill('0') << step << ".vtu"; - std::ofstream vtkfile(vtkfilename.str().c_str()); - if (!vtkfile) { - OPM_THROW(std::runtime_error, "Failed to open " << vtkfilename.str()); - } - Opm::DataMap dm; - dm["saturation"] = &state.saturation(); - dm["pressure"] = &state.pressure(); - std::vector cell_velocity; - Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); - dm["velocity"] = &cell_velocity; - Opm::writeVtkData(grid, dm, vtkfile); - } +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// static void outputStateVtk(const UnstructuredGrid& grid, +// const Opm::BlackoilState& state, +// const int step, +// const std::string& output_dir) +// { +// // Write data in VTK format. +// std::ostringstream vtkfilename; +// vtkfilename << output_dir << "/vtk_files"; +// boost::filesystem::path fpath(vtkfilename.str()); +// try { +// create_directories(fpath); +// } +// catch (...) { +// OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); +// } +// vtkfilename << "/output-" << std::setw(3) << std::setfill('0') << step << ".vtu"; +// std::ofstream vtkfile(vtkfilename.str().c_str()); +// if (!vtkfile) { +// OPM_THROW(std::runtime_error, "Failed to open " << vtkfilename.str()); +// } +// Opm::DataMap dm; +// dm["saturation"] = &state.saturation(); +// dm["pressure"] = &state.pressure(); +// std::vector cell_velocity; +// Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); +// dm["velocity"] = &cell_velocity; +// Opm::writeVtkData(grid, dm, vtkfile); +// } - static void outputStateMatlab(const UnstructuredGrid& grid, - const Opm::BlackoilState& state, - const int step, - const std::string& output_dir) - { - Opm::DataMap dm; - dm["saturation"] = &state.saturation(); - dm["pressure"] = &state.pressure(); - dm["surfvolume"] = &state.surfacevol(); - std::vector cell_velocity; - Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); - dm["velocity"] = &cell_velocity; - - // Write data (not grid) in Matlab format - for (Opm::DataMap::const_iterator it = dm.begin(); it != dm.end(); ++it) { - std::ostringstream fname; - fname << output_dir << "/" << it->first; - boost::filesystem::path fpath = fname.str(); - try { - create_directories(fpath); - } - catch (...) { - OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); - } - fname << "/" << std::setw(3) << std::setfill('0') << step << ".txt"; - std::ofstream file(fname.str().c_str()); - if (!file) { - OPM_THROW(std::runtime_error, "Failed to open " << fname.str()); - } - file.precision(15); - const std::vector& d = *(it->second); - std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); - } - } +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// static void outputStateMatlab(const UnstructuredGrid& grid, +// const Opm::BlackoilState& state, +// const int step, +// const std::string& output_dir) +// { +// Opm::DataMap dm; +// dm["saturation"] = &state.saturation(); +// dm["pressure"] = &state.pressure(); +// dm["surfvolume"] = &state.surfacevol(); +// std::vector cell_velocity; +// Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); +// dm["velocity"] = &cell_velocity; +// +// // Write data (not grid) in Matlab format +// for (Opm::DataMap::const_iterator it = dm.begin(); it != dm.end(); ++it) { +// std::ostringstream fname; +// fname << output_dir << "/" << it->first; +// boost::filesystem::path fpath = fname.str(); +// try { +// create_directories(fpath); +// } +// catch (...) { +// OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); +// } +// fname << "/" << std::setw(3) << std::setfill('0') << step << ".txt"; +// std::ofstream file(fname.str().c_str()); +// if (!file) { +// OPM_THROW(std::runtime_error, "Failed to open " << fname.str()); +// } +// file.precision(15); +// const std::vector& d = *(it->second); +// std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); +// } +// } static void outputWaterCut(const Opm::Watercut& watercut, @@ -348,9 +351,9 @@ namespace Opm timer.report(std::cout); if (output_ && (timer.currentStepNum() % output_interval_ == 0)) { if (output_vtk_) { - outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +// outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); } - outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +// outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); } SimulatorReport sreport; @@ -515,9 +518,9 @@ namespace Opm if (output_) { if (output_vtk_) { - outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +// outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); } - outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +// outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); outputWaterCut(watercut, output_dir_); if (wells_) { outputWellReport(wellreport, output_dir_); diff --git a/opm/core/simulator/SimulatorIncompTwophase.cpp b/opm/core/simulator/SimulatorIncompTwophase.cpp index 385e6ca4..a8efd60b 100644 --- a/opm/core/simulator/SimulatorIncompTwophase.cpp +++ b/opm/core/simulator/SimulatorIncompTwophase.cpp @@ -35,7 +35,8 @@ #include #include #include -#include +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// #include #include #include @@ -179,34 +180,36 @@ namespace Opm os.precision(8); } - static void outputStateVtk(const UnstructuredGrid& grid, - const Opm::TwophaseState& state, - const int step, - const std::string& output_dir) - { - // Write data in VTK format. - std::ostringstream vtkfilename; - vtkfilename << output_dir << "/vtk_files"; - boost::filesystem::path fpath(vtkfilename.str()); - try { - create_directories(fpath); - } - catch (...) { - OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); - } - vtkfilename << "/output-" << std::setw(3) << std::setfill('0') << step << ".vtu"; - std::ofstream vtkfile(vtkfilename.str().c_str()); - if (!vtkfile) { - OPM_THROW(std::runtime_error, "Failed to open " << vtkfilename.str()); - } - Opm::DataMap dm; - dm["saturation"] = &state.saturation(); - dm["pressure"] = &state.pressure(); - std::vector cell_velocity; - Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); - dm["velocity"] = &cell_velocity; - Opm::writeVtkData(grid, dm, vtkfile); - } + +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// static void outputStateVtk(const UnstructuredGrid& grid, +// const Opm::TwophaseState& state, +// const int step, +// const std::string& output_dir) +// { +// // Write data in VTK format. +// std::ostringstream vtkfilename; +// vtkfilename << output_dir << "/vtk_files"; +// boost::filesystem::path fpath(vtkfilename.str()); +// try { +// create_directories(fpath); +// } +// catch (...) { +// OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); +// } +// vtkfilename << "/output-" << std::setw(3) << std::setfill('0') << step << ".vtu"; +// std::ofstream vtkfile(vtkfilename.str().c_str()); +// if (!vtkfile) { +// OPM_THROW(std::runtime_error, "Failed to open " << vtkfilename.str()); +// } +// Opm::DataMap dm; +// dm["saturation"] = &state.saturation(); +// dm["pressure"] = &state.pressure(); +// std::vector cell_velocity; +// Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); +// dm["velocity"] = &cell_velocity; +// Opm::writeVtkData(grid, dm, vtkfile); +// } static void outputVectorMatlab(const std::string& name, const std::vector& vec, @@ -230,39 +233,40 @@ namespace Opm std::copy(vec.begin(), vec.end(), std::ostream_iterator(file, "\n")); } - static void outputStateMatlab(const UnstructuredGrid& grid, - const Opm::TwophaseState& state, - const int step, - const std::string& output_dir) - { - Opm::DataMap dm; - dm["saturation"] = &state.saturation(); - dm["pressure"] = &state.pressure(); - std::vector cell_velocity; - Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); - dm["velocity"] = &cell_velocity; - - // Write data (not grid) in Matlab format - for (Opm::DataMap::const_iterator it = dm.begin(); it != dm.end(); ++it) { - std::ostringstream fname; - fname << output_dir << "/" << it->first; - boost::filesystem::path fpath = fname.str(); - try { - create_directories(fpath); - } - catch (...) { - OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); - } - fname << "/" << std::setw(3) << std::setfill('0') << step << ".txt"; - std::ofstream file(fname.str().c_str()); - if (!file) { - OPM_THROW(std::runtime_error, "Failed to open " << fname.str()); - } - file.precision(15); - const std::vector& d = *(it->second); - std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); - } - } +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// static void outputStateMatlab(const UnstructuredGrid& grid, +// const Opm::TwophaseState& state, +// const int step, +// const std::string& output_dir) +// { +// Opm::DataMap dm; +// dm["saturation"] = &state.saturation(); +// dm["pressure"] = &state.pressure(); +// std::vector cell_velocity; +// Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); +// dm["velocity"] = &cell_velocity; +// +// // Write data (not grid) in Matlab format +// for (Opm::DataMap::const_iterator it = dm.begin(); it != dm.end(); ++it) { +// std::ostringstream fname; +// fname << output_dir << "/" << it->first; +// boost::filesystem::path fpath = fname.str(); +// try { +// create_directories(fpath); +// } +// catch (...) { +// OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); +// } +// fname << "/" << std::setw(3) << std::setfill('0') << step << ".txt"; +// std::ofstream file(fname.str().c_str()); +// if (!file) { +// OPM_THROW(std::runtime_error, "Failed to open " << fname.str()); +// } +// file.precision(15); +// const std::vector& d = *(it->second); +// std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); +// } +// } static void outputWaterCut(const Opm::Watercut& watercut, @@ -460,9 +464,9 @@ namespace Opm timer.report(*log_); if (output_ && (timer.currentStepNum() % output_interval_ == 0)) { if (output_vtk_) { - outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +// outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); } - outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +// outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); if (use_reorder_) { // This use of dynamic_cast is not ideal, but should be safe. outputVectorMatlab(std::string("reorder_it"), @@ -620,9 +624,9 @@ namespace Opm if (output_) { if (output_vtk_) { - outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +// outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); } - outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +// outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); if (use_reorder_) { // This use of dynamic_cast is not ideal, but should be safe. outputVectorMatlab(std::string("reorder_it"), diff --git a/opm/core/simulator/SimulatorOutput.cpp b/opm/core/simulator/SimulatorOutput.cpp index 1216a8fd..923c80b3 100644 --- a/opm/core/simulator/SimulatorOutput.cpp +++ b/opm/core/simulator/SimulatorOutput.cpp @@ -21,7 +21,8 @@ // we need complete definitions for these types #include -#include +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// #include #include #include // partial_sum @@ -45,13 +46,15 @@ SimulatorOutputBase::SimulatorOutputBase ( // process parameters into a writer. we don't setup a new chain in // every timestep! - , writer_ (std::move (OutputWriter::create (params, eclipseState, phaseUsage, grid))) +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// , writer_ (std::move (OutputWriter::create (params, eclipseState, phaseUsage, grid))) // always start from the first timestep , next_ (0) { // write the static initialization files, even before simulation starts - writer_->writeInit (*timer); +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// writer_->writeInit (*timer); } // default destructor is OK, just need to be defined @@ -88,7 +91,8 @@ SimulatorOutputBase::writeOutput () { // relay the request to the handlers (setup in the constructor // from parameters) - writer_->writeTimeStep (*timer_, *reservoirState_, *wellState_ , false); +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// writer_->writeTimeStep (*timer_, *reservoirState_, *wellState_ , false); // advance to the next reporting time ++next_; diff --git a/opm/core/simulator/SimulatorOutput.hpp b/opm/core/simulator/SimulatorOutput.hpp index 5d54fb50..c5833aac 100644 --- a/opm/core/simulator/SimulatorOutput.hpp +++ b/opm/core/simulator/SimulatorOutput.hpp @@ -34,7 +34,8 @@ namespace Opm { // forward definitions class Deck; class EclipseState; -class OutputWriter; +// 17.03.2016 Temporarily removed while moving functionality to opm-output +//class OutputWriter; namespace parameter { class ParameterGroup; } class SimulationDataContainer; class SimulatorTimer; @@ -84,7 +85,8 @@ protected: std::shared_ptr wellState_; /// Created locally and destructed together with us - std::unique_ptr writer_; +// 17.03.2016 Temporarily removed while moving functionality to opm-output +// std::unique_ptr writer_; /// Call the writers that were created based on the parameters virtual void writeOutput (); From 8e68a4d8158c2469d97e5cf019d164de4345275b Mon Sep 17 00:00:00 2001 From: "Kjell W. Kongsvik" Date: Thu, 17 Mar 2016 14:15:50 +0100 Subject: [PATCH 3/4] Deleted all files in opm/core/io This removes OutputWriter and eclipse, vtk, vag as this functionality has moved to opm-output. --- CMakeLists_files.cmake | 24 - opm/core/io/OutputWriter.cpp | 103 -- opm/core/io/OutputWriter.hpp | 118 -- opm/core/io/eclipse/CornerpointChopper.hpp | 457 ----- opm/core/io/eclipse/EclipseGridInspector.cpp | 350 ---- opm/core/io/eclipse/EclipseGridInspector.hpp | 103 -- opm/core/io/eclipse/EclipseIOUtil.hpp | 34 - opm/core/io/eclipse/EclipseReader.cpp | 247 --- opm/core/io/eclipse/EclipseReader.hpp | 39 - opm/core/io/eclipse/EclipseUnits.hpp | 63 - .../io/eclipse/EclipseWriteRFTHandler.cpp | 144 -- .../io/eclipse/EclipseWriteRFTHandler.hpp | 77 - opm/core/io/eclipse/EclipseWriter.cpp | 1508 ----------------- opm/core/io/eclipse/EclipseWriter.hpp | 139 -- opm/core/io/eclipse/writeECLData.cpp | 193 --- opm/core/io/eclipse/writeECLData.hpp | 46 - opm/core/io/vag/vag.cpp | 408 ----- opm/core/io/vag/vag.hpp | 165 -- opm/core/io/vtk/writeVtkData.cpp | 319 ---- opm/core/io/vtk/writeVtkData.hpp | 48 - 20 files changed, 4585 deletions(-) delete mode 100644 opm/core/io/OutputWriter.cpp delete mode 100644 opm/core/io/OutputWriter.hpp delete mode 100644 opm/core/io/eclipse/CornerpointChopper.hpp delete mode 100644 opm/core/io/eclipse/EclipseGridInspector.cpp delete mode 100644 opm/core/io/eclipse/EclipseGridInspector.hpp delete mode 100644 opm/core/io/eclipse/EclipseIOUtil.hpp delete mode 100644 opm/core/io/eclipse/EclipseReader.cpp delete mode 100644 opm/core/io/eclipse/EclipseReader.hpp delete mode 100644 opm/core/io/eclipse/EclipseUnits.hpp delete mode 100644 opm/core/io/eclipse/EclipseWriteRFTHandler.cpp delete mode 100644 opm/core/io/eclipse/EclipseWriteRFTHandler.hpp delete mode 100644 opm/core/io/eclipse/EclipseWriter.cpp delete mode 100644 opm/core/io/eclipse/EclipseWriter.hpp delete mode 100644 opm/core/io/eclipse/writeECLData.cpp delete mode 100644 opm/core/io/eclipse/writeECLData.hpp delete mode 100644 opm/core/io/vag/vag.cpp delete mode 100644 opm/core/io/vag/vag.hpp delete mode 100644 opm/core/io/vtk/writeVtkData.cpp delete mode 100644 opm/core/io/vtk/writeVtkData.hpp diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 0a011da8..0f157018 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -45,14 +45,6 @@ list (APPEND MAIN_SOURCE_FILES opm/core/grid/cpgpreprocess/uniquepoints.c opm/core/grid/grid.c opm/core/grid/grid_equal.cpp - opm/core/io/OutputWriter.cpp - opm/core/io/eclipse/EclipseGridInspector.cpp - opm/core/io/eclipse/EclipseReader.cpp - opm/core/io/eclipse/EclipseWriteRFTHandler.cpp - opm/core/io/eclipse/EclipseWriter.cpp - opm/core/io/eclipse/writeECLData.cpp - opm/core/io/vag/vag.cpp - opm/core/io/vtk/writeVtkData.cpp opm/core/linalg/LinearSolverFactory.cpp opm/core/linalg/LinearSolverInterface.cpp opm/core/linalg/LinearSolverIstl.cpp @@ -151,10 +143,6 @@ list (APPEND MAIN_SOURCE_FILES # originally generated with the command: # find tests -name '*.cpp' -a ! -wholename '*/not-unit/*' -printf '\t%p\n' | sort list (APPEND TEST_SOURCE_FILES - tests/test_writenumwells.cpp - tests/test_writeReadRestartFile.cpp - tests/test_EclipseWriter.cpp - tests/test_EclipseWriteRFTHandler.cpp tests/test_compressedpropertyaccess.cpp tests/test_dgbasis.cpp tests/test_cartgrid.cpp @@ -220,7 +208,6 @@ list (APPEND TEST_DATA_FILES tests/satfuncEPS_D.DATA tests/testBlackoilState1.DATA tests/testBlackoilState2.DATA - tests/testBlackoilState3.DATA tests/testPinch1.DATA tests/wells_manager_data.data tests/wells_manager_data_expanded.data @@ -294,17 +281,6 @@ list (APPEND PUBLIC_HEADER_FILES opm/core/grid/cpgpreprocess/geometry.h opm/core/grid/cpgpreprocess/preprocess.h opm/core/grid/cpgpreprocess/uniquepoints.h - opm/core/io/OutputWriter.hpp - opm/core/io/eclipse/CornerpointChopper.hpp - opm/core/io/eclipse/EclipseGridInspector.hpp - opm/core/io/eclipse/EclipseIOUtil.hpp - opm/core/io/eclipse/EclipseReader.hpp - opm/core/io/eclipse/EclipseUnits.hpp - opm/core/io/eclipse/EclipseWriteRFTHandler.hpp - opm/core/io/eclipse/EclipseWriter.hpp - opm/core/io/eclipse/writeECLData.hpp - opm/core/io/vag/vag.hpp - opm/core/io/vtk/writeVtkData.hpp opm/core/linalg/LinearSolverFactory.hpp opm/core/linalg/LinearSolverInterface.hpp opm/core/linalg/LinearSolverIstl.hpp diff --git a/opm/core/io/OutputWriter.cpp b/opm/core/io/OutputWriter.cpp deleted file mode 100644 index a3eb0614..00000000 --- a/opm/core/io/OutputWriter.cpp +++ /dev/null @@ -1,103 +0,0 @@ -#include "OutputWriter.hpp" - -#include -#include -#include -#include - -#include -#include -#include // unique_ptr - -using namespace std; -using namespace Opm; -using namespace Opm::parameter; - -namespace { - -/// Multiplexer over a list of output writers -struct MultiWriter : public OutputWriter { - /// Shorthand for a list of owned output writers - typedef forward_list > writers_t; - typedef writers_t::iterator it_t; - typedef unique_ptr ptr_t; - - /// Adopt a list of writers - MultiWriter (ptr_t writers) : writers_ (std::move (writers)) { } - - /// Forward the call to all writers - virtual void writeInit(const SimulatorTimerInterface &timer) { - for (it_t it = writers_->begin (); it != writers_->end (); ++it) { - (*it)->writeInit (timer); - } - } - - virtual void writeTimeStep(const SimulatorTimerInterface& timer, - const SimulationDataContainer& reservoirState, - const WellState& wellState, - bool isSubstep) { - for (it_t it = writers_->begin (); it != writers_->end(); ++it) { - (*it)->writeTimeStep (timer, reservoirState, wellState, isSubstep); - } - } - -private: - ptr_t writers_; -}; - -/// Psuedo-constructor, can appear in template -template unique_ptr -create (const ParameterGroup& params, - std::shared_ptr eclipseState, - const Opm::PhaseUsage &phaseUsage, - std::shared_ptr grid) { - return unique_ptr (new Format (params, - eclipseState, - phaseUsage, - grid->number_of_cells, - grid->global_cell)); -} - -/// Map between keyword in configuration and the corresponding -/// constructor function (type) that should be called when detected. -/// The writer must have a constructor which takes params and parser. -/// -/// If you want to add more possible writer formats, just add them -/// to the list below! -typedef map (*)( - const ParameterGroup&, - std::shared_ptr eclipseState, - const Opm::PhaseUsage &phaseUsage, - std::shared_ptr )> map_t; -map_t FORMATS = { - { "output_ecl", &create }, -}; - -} // anonymous namespace - -unique_ptr -OutputWriter::create (const ParameterGroup& params, - std::shared_ptr eclipseState, - const Opm::PhaseUsage &phaseUsage, - std::shared_ptr grid) { - // allocate a list which will be filled with writers. this list - // is initially empty (no output). - MultiWriter::ptr_t list (new MultiWriter::writers_t ()); - - // loop through the map and see if we can find the key that is - // specified there - typedef map_t::iterator map_it_t; - for (map_it_t it = FORMATS.begin (); it != FORMATS.end(); ++it) { - // keyword which would indicate that this format should be used - const std::string name (it->first); - - // invoke the constructor for the type if we found the keyword - // and put the pointer to this writer onto the list - if (params.getDefault (name, false)) { - list->push_front (it->second (params, eclipseState, phaseUsage, grid)); - } - } - - // create a multiplexer from the list of formats we found - return unique_ptr (new MultiWriter (std::move (list))); -} diff --git a/opm/core/io/OutputWriter.hpp b/opm/core/io/OutputWriter.hpp deleted file mode 100644 index de9f3460..00000000 --- a/opm/core/io/OutputWriter.hpp +++ /dev/null @@ -1,118 +0,0 @@ -/* - Copyright (c) 2013 Uni Research 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 . -*/ - -#ifndef OPM_OUTPUT_WRITER_HPP -#define OPM_OUTPUT_WRITER_HPP - -#include // unique_ptr, shared_ptr -#include - -struct UnstructuredGrid; - -namespace Opm { - -// forward declaration -class EclipseState; -namespace parameter { class ParameterGroup; } -class SimulationDataContainer; -class WellState; -struct PhaseUsage; - -/*! - * Interface for writing non-compositional (blackoil, two-phase) simulation - * state to files. - * - * Use the create() function to setup a chain of writer based on the - * configuration values, e.g. - * - * \example - * \code{.cpp} - * ParameterGroup params (argc, argv, false); - * auto parser = std::make_shared ( - * params.get ("deck_filename")); - * - * std::unique_ptr writer = - * OutputWriter::create (params, parser); - * - * // before the first timestep - * writer->writeInit (timer); - * - * // after each timestep - * writer->writeTimeStep (timer, state, wellState); - * - * \endcode - */ -class OutputWriter { -public: - /// Allow derived classes to be used in the unique_ptr that is returned - /// from the create() method. (Every class that should be delete'd should - /// have a proper constructor, and if the base class isn't virtual then - /// the compiler won't call the right one when the unique_ptr goes out of - /// scope). - virtual ~OutputWriter () { } - - /** - * Write the static data (grid, PVT curves, etc) to disk. - * - * This routine should be called before the first timestep (i.e. when - * timer.currentStepNum () == 0) - */ - virtual void writeInit(const SimulatorTimerInterface &timer) = 0; - - /*! - * \brief Write a blackoil reservoir state to disk for later inspection with - * visualization tools like ResInsight - * - * \param[in] timer The timer providing time, time step, etc. information - * \param[in] reservoirState The thermodynamic state of the reservoir - * \param[in] wellState The production/injection data for all wells - * - * This routine should be called after the timestep has been advanced, - * i.e. timer.currentStepNum () > 0. - */ - virtual void writeTimeStep(const SimulatorTimerInterface& timer, - const SimulationDataContainer& reservoirState, - const WellState& wellState, - bool isSubstep) = 0; - - /*! - * Create a suitable set of output formats based on configuration. - * - * @param params Configuration properties. This function will setup a - * multiplexer of applicable output formats based on the - * desired configuration values. - * - * @param deck Input deck used to set up the simulation. - * - * @param eclipseState The internalized input deck. - * - * @return Pointer to a multiplexer to all applicable output formats. - * - * @see Opm::share_obj - */ - static std::unique_ptr - create (const parameter::ParameterGroup& params, - std::shared_ptr eclipseState, - const Opm::PhaseUsage &phaseUsage, - std::shared_ptr grid); -}; - -} // namespace Opm - -#endif /* OPM_OUTPUT_WRITER_HPP */ diff --git a/opm/core/io/eclipse/CornerpointChopper.hpp b/opm/core/io/eclipse/CornerpointChopper.hpp deleted file mode 100644 index 57fb09ae..00000000 --- a/opm/core/io/eclipse/CornerpointChopper.hpp +++ /dev/null @@ -1,457 +0,0 @@ -/* - Copyright 2010 SINTEF ICT, Applied Mathematics. - - 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 . -*/ - -#ifndef OPM_CORNERPOINTCHOPPER_HEADER_INCLUDED -#define OPM_CORNERPOINTCHOPPER_HEADER_INCLUDED - -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - -namespace Opm -{ - - class CornerPointChopper - { - public: - CornerPointChopper(const std::string& file) - { - Opm::ParseContext parseContext; - Opm::ParserPtr parser(new Opm::Parser()); - deck_ = parser->parseFile(file , parseContext); - - metricUnits_.reset(Opm::UnitSystem::newMETRIC()); - - const auto& specgridRecord = deck_->getKeyword("SPECGRID").getRecord(0); - dims_[0] = specgridRecord.getItem("NX").get< int >(0); - dims_[1] = specgridRecord.getItem("NY").get< int >(0); - dims_[2] = specgridRecord.getItem("NZ").get< int >(0); - - int layersz = 8*dims_[0]*dims_[1]; - const std::vector& ZCORN = deck_->getKeyword("ZCORN").getRawDoubleData(); - botmax_ = *std::max_element(ZCORN.begin(), ZCORN.begin() + layersz/2); - topmin_ = *std::min_element(ZCORN.begin() + dims_[2]*layersz - layersz/2, - ZCORN.begin() + dims_[2]*layersz); - - abszmax_ = *std::max_element(ZCORN.begin(), ZCORN.end()); - abszmin_ = *std::min_element(ZCORN.begin(), ZCORN.end()); - - std::cout << "Parsed grdecl file with dimensions (" - << dims_[0] << ", " << dims_[1] << ", " << dims_[2] << ")" << std::endl; - } - - - - - const int* dimensions() const - { - return dims_; - } - - - - - const int* newDimensions() const - { - return new_dims_; - } - - - - - const std::pair zLimits() const - { - return std::make_pair(botmax_, topmin_); - } - - const std::pair abszLimits() const - { - return std::make_pair(abszmin_, abszmax_); - } - - - void verifyInscribedShoebox(int imin, int ilen, int imax, - int jmin, int jlen, int jmax, - double zmin, double zlen, double zmax) - { - if (imin < 0) { - std::cerr << "Error! imin < 0 (imin = " << imin << ")\n"; - throw std::runtime_error("Inconsistent user input."); - } - if (ilen > dims_[0]) { - std::cerr << "Error! ilen larger than grid (ilen = " << ilen <<")\n"; - throw std::runtime_error("Inconsistent user input."); - } - if (imax > dims_[0]) { - std::cerr << "Error! imax larger than input grid (imax = " << imax << ")\n"; - throw std::runtime_error("Inconsistent user input."); - } - if (jmin < 0) { - std::cerr << "Error! jmin < 0 (jmin = " << jmin << ")\n"; - throw std::runtime_error("Inconsistent user input."); - } - if (jlen > dims_[1]) { - std::cerr << "Error! jlen larger than grid (jlen = " << jlen <<")\n"; - throw std::runtime_error("Inconsistent user input."); - } - if (jmax > dims_[1]) { - std::cerr << "Error! jmax larger than input grid (jmax = " << jmax << ")\n"; - throw std::runtime_error("Inconsistent user input."); - } - if (zmin < abszmin_) { - std::cerr << "Error! zmin ("<< zmin << ") less than minimum ZCORN value ("<< abszmin_ << ")\n"; - throw std::runtime_error("Inconsistent user input."); - } - if (zmax > abszmax_) { - std::cerr << "Error! zmax ("<< zmax << ") larger than maximal ZCORN value ("<< abszmax_ << ")\n"; - throw std::runtime_error("Inconsistent user input."); - } - if (zlen > (abszmax_ - abszmin_)) { - std::cerr << "Error! zlen ("<< zlen <<") larger than maximal ZCORN (" << abszmax_ << ") minus minimal ZCORN ("<< abszmin_ <<")\n"; - throw std::runtime_error("Inconsistent user input."); - } - } - - void chop(int imin, int imax, int jmin, int jmax, double zmin, double zmax, bool resettoorigin=true) - { - new_dims_[0] = imax - imin; - new_dims_[1] = jmax - jmin; - - // Filter the coord field - const std::vector& COORD = deck_->getKeyword("COORD").getRawDoubleData(); - int num_coord = COORD.size(); - if (num_coord != 6*(dims_[0] + 1)*(dims_[1] + 1)) { - std::cerr << "Error! COORD size (" << COORD.size() << ") not consistent with SPECGRID\n"; - throw std::runtime_error("Inconsistent COORD and SPECGRID."); - } - int num_new_coord = 6*(new_dims_[0] + 1)*(new_dims_[1] + 1); - double x_correction = COORD[6*((dims_[0] + 1)*jmin + imin)]; - double y_correction = COORD[6*((dims_[0] + 1)*jmin + imin) + 1]; - new_COORD_.resize(num_new_coord, 1e100); - for (int j = jmin; j < jmax + 1; ++j) { - for (int i = imin; i < imax + 1; ++i) { - int pos = (dims_[0] + 1)*j + i; - int new_pos = (new_dims_[0] + 1)*(j-jmin) + (i-imin); - // Copy all 6 coordinates for a pillar. - std::copy(COORD.begin() + 6*pos, COORD.begin() + 6*(pos + 1), new_COORD_.begin() + 6*new_pos); - if (resettoorigin) { - // Substract lowest x value from all X-coords, similarly for y, and truncate in z-direction - new_COORD_[6*new_pos] -= x_correction; - new_COORD_[6*new_pos + 1] -= y_correction; - new_COORD_[6*new_pos + 2] = 0; - new_COORD_[6*new_pos + 3] -= x_correction; - new_COORD_[6*new_pos + 4] -= y_correction; - new_COORD_[6*new_pos + 5] = zmax-zmin; - } - } - } - - // Get the z limits, check if they must be changed to make a shoe-box. - // This means that zmin must be greater than or equal to the highest - // coordinate of the bottom surface, while zmax must be less than or - // equal to the lowest coordinate of the top surface. - int layersz = 8*dims_[0]*dims_[1]; - const std::vector& ZCORN = deck_->getKeyword("ZCORN").getRawDoubleData(); - int num_zcorn = ZCORN.size(); - if (num_zcorn != layersz*dims_[2]) { - std::cerr << "Error! ZCORN size (" << ZCORN.size() << ") not consistent with SPECGRID\n"; - throw std::runtime_error("Inconsistent ZCORN and SPECGRID."); - } - - zmin = std::max(zmin, botmax_); - zmax = std::min(zmax, topmin_); - if (zmin >= zmax) { - std::cerr << "Error: zmin >= zmax (zmin = " << zmin << ", zmax = " << zmax << ")\n"; - throw std::runtime_error("zmin >= zmax"); - } - std::cout << "Chopping subsample, i: (" << imin << "--" << imax << ") j: (" << jmin << "--" << jmax << ") z: (" << zmin << "--" << zmax << ")" << std::endl; - - // We must find the maximum and minimum k value for the given z limits. - // First, find the first layer with a z-coordinate strictly above zmin. - int kmin = -1; - for (int k = 0; k < dims_[2]; ++k) { - double layer_max = *std::max_element(ZCORN.begin() + k*layersz, ZCORN.begin() + (k + 1)*layersz); - if (layer_max > zmin) { - kmin = k; - break; - } - } - // Then, find the last layer with a z-coordinate strictly below zmax. - int kmax = -1; - for (int k = dims_[2]; k > 0; --k) { - double layer_min = *std::min_element(ZCORN.begin() + (k - 1)*layersz, ZCORN.begin() + k*layersz); - if (layer_min < zmax) { - kmax = k; - break; - } - } - new_dims_[2] = kmax - kmin; - - // Filter the ZCORN field, build mapping from new to old cells. - double z_origin_correction = 0.0; - if (resettoorigin) { - z_origin_correction = zmin; - } - new_ZCORN_.resize(8*new_dims_[0]*new_dims_[1]*new_dims_[2], 1e100); - new_to_old_cell_.resize(new_dims_[0]*new_dims_[1]*new_dims_[2], -1); - int cellcount = 0; - int delta[3] = { 1, 2*dims_[0], 4*dims_[0]*dims_[1] }; - int new_delta[3] = { 1, 2*new_dims_[0], 4*new_dims_[0]*new_dims_[1] }; - for (int k = kmin; k < kmax; ++k) { - for (int j = jmin; j < jmax; ++j) { - for (int i = imin; i < imax; ++i) { - new_to_old_cell_[cellcount++] = dims_[0]*dims_[1]*k + dims_[0]*j + i; - int old_ix = 2*(i*delta[0] + j*delta[1] + k*delta[2]); - int new_ix = 2*((i-imin)*new_delta[0] + (j-jmin)*new_delta[1] + (k-kmin)*new_delta[2]); - int old_indices[8] = { old_ix, old_ix + delta[0], - old_ix + delta[1], old_ix + delta[1] + delta[0], - old_ix + delta[2], old_ix + delta[2] + delta[0], - old_ix + delta[2] + delta[1], old_ix + delta[2] + delta[1] + delta[0] }; - int new_indices[8] = { new_ix, new_ix + new_delta[0], - new_ix + new_delta[1], new_ix + new_delta[1] + new_delta[0], - new_ix + new_delta[2], new_ix + new_delta[2] + new_delta[0], - new_ix + new_delta[2] + new_delta[1], new_ix + new_delta[2] + new_delta[1] + new_delta[0] }; - for (int cc = 0; cc < 8; ++cc) { - new_ZCORN_[new_indices[cc]] = std::min(zmax, std::max(zmin, ZCORN[old_indices[cc]])) - z_origin_correction; - } - } - } - } - - filterIntegerField("ACTNUM", new_ACTNUM_); - filterDoubleField("PORO", new_PORO_); - filterDoubleField("NTG", new_NTG_); - filterDoubleField("SWCR", new_SWCR_); - filterDoubleField("SOWCR", new_SOWCR_); - filterDoubleField("PERMX", new_PERMX_); - filterDoubleField("PERMY", new_PERMY_); - filterDoubleField("PERMZ", new_PERMZ_); - filterIntegerField("SATNUM", new_SATNUM_); - } - - /// Return a sub-deck with fields corresponding to the selected subset. - Opm::DeckConstPtr subDeck() - { - Opm::DeckPtr subDeck(new Opm::Deck); - - Opm::DeckKeyword specGridKw("SPECGRID"); - Opm::DeckRecord specGridRecord; - - auto nxItem = Opm::DeckItem::make< int >("NX"); - auto nyItem = Opm::DeckItem::make< int >("NY"); - auto nzItem = Opm::DeckItem::make< int >("NZ"); - auto numresItem = Opm::DeckItem::make< int >("NUMRES"); - auto coordTypeItem = Opm::DeckItem::make< std::string >("COORD_TYPE"); - - nxItem.push_back(new_dims_[0]); - nyItem.push_back(new_dims_[1]); - nzItem.push_back(new_dims_[2]); - numresItem.push_back(1); - coordTypeItem.push_back("F"); - - specGridRecord.addItem(std::move(nxItem)); - specGridRecord.addItem(std::move(nyItem)); - specGridRecord.addItem(std::move(nzItem)); - specGridRecord.addItem(std::move(numresItem)); - specGridRecord.addItem(std::move(coordTypeItem)); - - specGridKw.addRecord(std::move(specGridRecord)); - - subDeck->addKeyword(std::move(specGridKw)); - addDoubleKeyword_(subDeck, "COORD", /*dimension=*/"Length", new_COORD_); - addDoubleKeyword_(subDeck, "ZCORN", /*dimension=*/"Length", new_ZCORN_); - addIntKeyword_(subDeck, "ACTNUM", new_ACTNUM_); - addDoubleKeyword_(subDeck, "PORO", /*dimension=*/"1", new_PORO_); - addDoubleKeyword_(subDeck, "NTG", /*dimension=*/"1", new_NTG_); - addDoubleKeyword_(subDeck, "SWCR", /*dimension=*/"1", new_SWCR_); - addDoubleKeyword_(subDeck, "SOWCR", /*dimension=*/"1", new_SOWCR_); - addDoubleKeyword_(subDeck, "PERMX", /*dimension=*/"Permeability", new_PERMX_); - addDoubleKeyword_(subDeck, "PERMY", /*dimension=*/"Permeability", new_PERMY_); - addDoubleKeyword_(subDeck, "PERMZ", /*dimension=*/"Permeability", new_PERMZ_); - addIntKeyword_(subDeck, "SATNUM", new_SATNUM_); - return subDeck; - } - void writeGrdecl(const std::string& filename) - { - // Output new versions of SPECGRID, COORD, ZCORN, ACTNUM, PERMX, PORO, SATNUM. - std::ofstream out(filename.c_str()); - if (!out) { - std::cerr << "Could not open file " << filename << "\n"; - throw std::runtime_error("Could not open output file."); - } - out << "SPECGRID\n" << new_dims_[0] << ' ' << new_dims_[1] << ' ' << new_dims_[2] - << " 1 F\n/\n\n"; - - out.precision(15); - out.setf(std::ios::scientific); - - outputField(out, new_COORD_, "COORD", /* nl = */ 3); - outputField(out, new_ZCORN_, "ZCORN", /* nl = */ 4); - outputField(out, new_ACTNUM_, "ACTNUM"); - outputField(out, new_PORO_, "PORO", 4); - if (hasNTG()) {outputField(out, new_NTG_, "NTG", 4);} - if (hasSWCR()) {outputField(out, new_SWCR_, "SWCR", 4);} - if (hasSOWCR()) {outputField(out, new_SOWCR_, "SOWCR", 4);} - outputField(out, new_PERMX_, "PERMX", 4); - outputField(out, new_PERMY_, "PERMY", 4); - outputField(out, new_PERMZ_, "PERMZ", 4); - outputField(out, new_SATNUM_, "SATNUM"); - } - bool hasNTG() const {return !new_NTG_.empty(); } - bool hasSWCR() const {return !new_SWCR_.empty(); } - bool hasSOWCR() const {return !new_SOWCR_.empty(); } - - private: - Opm::DeckConstPtr deck_; - std::shared_ptr metricUnits_; - - double botmax_; - double topmin_; - double abszmin_; - double abszmax_; - std::vector new_COORD_; - std::vector new_ZCORN_; - std::vector new_ACTNUM_; - std::vector new_PORO_; - std::vector new_NTG_; - std::vector new_SWCR_; - std::vector new_SOWCR_; - std::vector new_PERMX_; - std::vector new_PERMY_; - std::vector new_PERMZ_; - std::vector new_SATNUM_; - int dims_[3]; - int new_dims_[3]; - std::vector new_to_old_cell_; - - void addDoubleKeyword_(Opm::DeckPtr subDeck, - const std::string& keywordName, - const std::string& dimensionString, - const std::vector& data) - { - if (data.empty()) - return; - - Opm::DeckKeyword dataKw(keywordName); - Opm::DeckRecord dataRecord; - auto dataItem = Opm::DeckItem::make< double >("DATA"); - - for (size_t i = 0; i < data.size(); ++i) { - dataItem.push_back(data[i]); - } - - std::shared_ptr dimension = metricUnits_->parse(dimensionString); - dataItem.push_backDimension(/*active=*/dimension, /*default=*/dimension); - - dataRecord.addItem(std::move(dataItem)); - dataKw.addRecord(std::move(dataRecord)); - subDeck->addKeyword(std::move(dataKw)); - } - - void addIntKeyword_(Opm::DeckPtr subDeck, - const std::string& keywordName, - const std::vector& data) - { - if (data.empty()) - return; - - Opm::DeckKeyword dataKw(keywordName); - Opm::DeckRecord dataRecord; - auto dataItem = Opm::DeckItem::make< int >("DATA"); - - for (size_t i = 0; i < data.size(); ++i) { - dataItem.push_back(data[i]); - } - - dataRecord.addItem(std::move(dataItem)); - dataKw.addRecord(std::move(dataRecord)); - subDeck->addKeyword(std::move(dataKw)); - } - - template - void outputField(std::ostream& os, - const std::vector& field, - const std::string& keyword, - const typename std::vector::size_type nl = 20) - { - if (field.empty()) return; - - os << keyword << '\n'; - - typedef typename std::vector::size_type sz_t; - - const sz_t n = field.size(); - for (sz_t i = 0; i < n; ++i) { - os << field[i] - << (((i + 1) % nl == 0) ? '\n' : ' '); - } - if (n % nl != 0) { - os << '\n'; - } - os << "/\n\n"; - } - - - - template - void filterField(const std::vector& field, - std::vector& output_field) - { - int sz = new_to_old_cell_.size(); - output_field.resize(sz); - for (int i = 0; i < sz; ++i) { - output_field[i] = field[new_to_old_cell_[i]]; - } - } - - void filterDoubleField(const std::string& keyword, std::vector& output_field) - { - if (deck_->hasKeyword(keyword)) { - const std::vector& field = deck_->getKeyword(keyword).getRawDoubleData(); - filterField(field, output_field); - } - } - - void filterIntegerField(const std::string& keyword, std::vector& output_field) - { - if (deck_->hasKeyword(keyword)) { - const std::vector& field = deck_->getKeyword(keyword).getIntData(); - filterField(field, output_field); - } - } - - }; - -} - - - - -#endif // OPM_CORNERPOINTCHOPPER_HEADER_INCLUDED diff --git a/opm/core/io/eclipse/EclipseGridInspector.cpp b/opm/core/io/eclipse/EclipseGridInspector.cpp deleted file mode 100644 index ad7666f5..00000000 --- a/opm/core/io/eclipse/EclipseGridInspector.cpp +++ /dev/null @@ -1,350 +0,0 @@ -//=========================================================================== -// -// File: EclipseGridInspector.C -// -// Created: Mon Jun 2 12:17:51 2008 -// -// Author: Atgeirr F Rasmussen -// -// $Date$ -// -// $Revision$ -// -// Revision: $Id: EclipseGridInspector.C,v 1.2 2008/08/18 14:16:13 atgeirr Exp $ -// -//=========================================================================== - -/* - Copyright 2009, 2010 SINTEF ICT, Applied Mathematics. - Copyright 2009, 2010 Statoil ASA. - - 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 . -*/ - -#if HAVE_CONFIG_H -#include "config.h" -#endif -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - -namespace Opm -{ -EclipseGridInspector::EclipseGridInspector(Opm::DeckConstPtr deck) - : deck_(deck) -{ - init_(); -} - -void EclipseGridInspector::init_() -{ - if (!deck_->hasKeyword("COORD")) { - OPM_THROW(std::runtime_error, "Needed field \"COORD\" is missing in file"); - } - if (!deck_->hasKeyword("ZCORN")) { - OPM_THROW(std::runtime_error, "Needed field \"ZCORN\" is missing in file"); - } - - if (deck_->hasKeyword("SPECGRID")) { - const auto& specgridRecord = - deck_->getKeyword("SPECGRID").getRecord(0); - logical_gridsize_[0] = specgridRecord.getItem("NX").get< int >(0); - logical_gridsize_[1] = specgridRecord.getItem("NY").get< int >(0); - logical_gridsize_[2] = specgridRecord.getItem("NZ").get< int >(0); - } else if (deck_->hasKeyword("DIMENS")) { - const auto& dimensRecord = - deck_->getKeyword("DIMENS").getRecord(0); - logical_gridsize_[0] = dimensRecord.getItem("NX").get< int >(0); - logical_gridsize_[1] = dimensRecord.getItem("NY").get< int >(0); - logical_gridsize_[2] = dimensRecord.getItem("NZ").get< int >(0); - } else { - OPM_THROW(std::runtime_error, "Found neither SPECGRID nor DIMENS in file. At least one is needed."); - } - -} - -/** - Return the dip slopes for the cell relative to xy-plane in x- and y- direction. - Dip slope is average rise in positive x-direction over cell length in x-direction. - Similarly for y. - - Current implementation is for vertical pillars, but is not difficult to fix. - - @returns a std::pair with x-dip in first component and y-dip in second. -*/ -std::pair EclipseGridInspector::cellDips(int i, int j, int k) const -{ - checkLogicalCoords(i, j, k); - const std::vector& pillc = - deck_->getKeyword("COORD").getSIDoubleData(); - int num_pillars = (logical_gridsize_[0] + 1)*(logical_gridsize_[1] + 1); - if (6*num_pillars != int(pillc.size())) { - throw std::runtime_error("Wrong size of COORD field."); - } - const std::vector& z = - deck_->getKeyword("ZCORN").getSIDoubleData(); - int num_cells = logical_gridsize_[0]*logical_gridsize_[1]*logical_gridsize_[2]; - if (8*num_cells != int(z.size())) { - throw std::runtime_error("Wrong size of ZCORN field"); - } - - // Pick ZCORN-value for all 8 corners of the given cell - std::array cellz = cellZvals(i, j, k); - - // Compute rise in positive x-direction for all four edges (and then find mean) - // Current implementation is for regularly placed and vertical pillars! - int numxpill = logical_gridsize_[0] + 1; - int pix = i + j*numxpill; - double cell_xlength = pillc[6*(pix + 1)] - pillc[6*pix]; - flush(std::cout); - double xrise[4] = { (cellz[1] - cellz[0])/cell_xlength, // LLL -> HLL - (cellz[3] - cellz[2])/cell_xlength, // LHL -> HHL - (cellz[5] - cellz[4])/cell_xlength, // LLH -> HLH - (cellz[7] - cellz[6])/cell_xlength}; // LHH -> HHH - - double cell_ylength = pillc[6*(pix + numxpill) + 1] - pillc[6*pix + 1]; - double yrise[4] = { (cellz[2] - cellz[0])/cell_ylength, // LLL -> LHL - (cellz[3] - cellz[1])/cell_ylength, // HLL -> HHL - (cellz[6] - cellz[4])/cell_ylength, // LLH -> LHH - (cellz[7] - cellz[5])/cell_ylength}; // HLH -> HHH - - - // Now ignore those edges that touch the global top or bottom surface - // of the entire grdecl model. This is to avoid bias, as these edges probably - // don't follow an overall dip for the model if it exists. - int x_edges = 4; - int y_edges = 4; - std::array gridlimits = getGridLimits(); - double zmin = gridlimits[4]; - double zmax = gridlimits[5]; - // LLL -> HLL - if ((cellz[1] == zmin) || (cellz[0] == zmin)) { - xrise[0] = 0; x_edges--; - } - // LHL -> HHL - if ((cellz[2] == zmin) || (cellz[3] == zmin)) { - xrise[1] = 0; x_edges--; - } - // LLH -> HLH - if ((cellz[4] == zmax) || (cellz[5] == zmax)) { - xrise[2] = 0; x_edges--; - } - // LHH -> HHH - if ((cellz[6] == zmax) || (cellz[7] == zmax)) { - xrise[3] = 0; x_edges--; - } - // LLL -> LHL - if ((cellz[0] == zmin) || (cellz[2] == zmin)) { - yrise[0] = 0; y_edges--; - } - // HLL -> HHL - if ((cellz[1] == zmin) || (cellz[3] == zmin)) { - yrise[1] = 0; y_edges--; - } - // LLH -> LHH - if ((cellz[6] == zmax) || (cellz[4] == zmax)) { - yrise[2] = 0; y_edges--; - } - // HLH -> HHH - if ((cellz[7] == zmax) || (cellz[5] == zmax)) { - yrise[3] = 0; y_edges--; - } - - return std::make_pair( (xrise[0] + xrise[1] + xrise[2] + xrise[3])/x_edges, - (yrise[0] + yrise[1] + yrise[2] + yrise[3])/y_edges); -} -/** - Wrapper for cellDips(i, j, k). -*/ -std::pair EclipseGridInspector::cellDips(int cell_idx) const -{ - std::array idxs = cellIdxToLogicalCoords(cell_idx); - return cellDips(idxs[0], idxs[1], idxs[2]); -} - -std::array EclipseGridInspector::cellIdxToLogicalCoords(int cell_idx) const -{ - - int i,j,k; // Position of cell in cell hierarchy - int horIdx = (cell_idx+1) - int(std::floor(((double)(cell_idx+1))/((double)(logical_gridsize_[0]*logical_gridsize_[1]))))*logical_gridsize_[0]*logical_gridsize_[1]; // index in the corresponding horizon - if (horIdx == 0) { - horIdx = logical_gridsize_[0]*logical_gridsize_[1]; - } - i = horIdx - int(std::floor(((double)horIdx)/((double)logical_gridsize_[0])))*logical_gridsize_[0]; - if (i == 0) { - i = logical_gridsize_[0]; - } - j = (horIdx-i)/logical_gridsize_[0]+1; - k = ((cell_idx+1)-logical_gridsize_[0]*(j-1)-1)/(logical_gridsize_[0]*logical_gridsize_[1])+1; - - std::array a = {{i-1, j-1, k-1}}; - return a; //std::array {{i-1, j-1, k-1}}; -} - -double EclipseGridInspector::cellVolumeVerticalPillars(int i, int j, int k) const -{ - // Checking parameters and obtaining values from parser. - checkLogicalCoords(i, j, k); - const std::vector& pillc = - deck_->getKeyword("COORD").getSIDoubleData(); - int num_pillars = (logical_gridsize_[0] + 1)*(logical_gridsize_[1] + 1); - if (6*num_pillars != int(pillc.size())) { - throw std::runtime_error("Wrong size of COORD field."); - } - const std::vector& z = - deck_->getKeyword("ZCORN").getSIDoubleData(); - int num_cells = logical_gridsize_[0]*logical_gridsize_[1]*logical_gridsize_[2]; - if (8*num_cells != int(z.size())) { - throw std::runtime_error("Wrong size of ZCORN field"); - } - - // Computing the base area as half the 2d cross product of the diagonals. - int numxpill = logical_gridsize_[0] + 1; - int pix = i + j*numxpill; - double px[4] = { pillc[6*pix], - pillc[6*(pix + 1)], - pillc[6*(pix + numxpill)], - pillc[6*(pix + numxpill + 1)] }; - double py[4] = { pillc[6*pix + 1], - pillc[6*(pix + 1) + 1], - pillc[6*(pix + numxpill) + 1], - pillc[6*(pix + numxpill + 1) + 1] }; - double diag1[2] = { px[3] - px[0], py[3] - py[0] }; - double diag2[2] = { px[2] - px[1], py[2] - py[1] }; - double area = 0.5*(diag1[0]*diag2[1] - diag1[1]*diag2[0]); - - // Computing the average of the z-differences along each pillar. - int delta[3] = { 1, - 2*logical_gridsize_[0], - 4*logical_gridsize_[0]*logical_gridsize_[1] }; - int ix = 2*(i*delta[0] + j*delta[1] + k*delta[2]); - double cellz[8] = { z[ix], z[ix + delta[0]], - z[ix + delta[1]], z[ix + delta[1] + delta[0]], - z[ix + delta[2]], z[ix + delta[2] + delta[0]], - z[ix + delta[2] + delta[1]], z[ix + delta[2] + delta[1] + delta[0]] }; - double diffz[4] = { cellz[4] - cellz[0], - cellz[5] - cellz[1], - cellz[6] - cellz[2], - cellz[7] - cellz[3] }; - double averzdiff = 0.25*std::accumulate(diffz, diffz + 4, 0.0); - return averzdiff*area; -} - - -double EclipseGridInspector::cellVolumeVerticalPillars(int cell_idx) const -{ - std::array idxs = cellIdxToLogicalCoords(cell_idx); - return cellVolumeVerticalPillars(idxs[0], idxs[1], idxs[2]); -} - -void EclipseGridInspector::checkLogicalCoords(int i, int j, int k) const -{ - if (i < 0 || i >= logical_gridsize_[0]) - throw std::runtime_error("First coordinate out of bounds"); - if (j < 0 || j >= logical_gridsize_[1]) - throw std::runtime_error("Second coordinate out of bounds"); - if (k < 0 || k >= logical_gridsize_[2]) - throw std::runtime_error("Third coordinate out of bounds"); -} - - -std::array EclipseGridInspector::getGridLimits() const -{ - if (! (deck_->hasKeyword("COORD") && deck_->hasKeyword("ZCORN") && deck_->hasKeyword("SPECGRID")) ) { - throw std::runtime_error("EclipseGridInspector: Grid does not have SPECGRID, COORD, and ZCORN, can't find dimensions."); - } - - std::vector coord = deck_->getKeyword("COORD").getSIDoubleData(); - std::vector zcorn = deck_->getKeyword("ZCORN").getSIDoubleData(); - - double xmin = +DBL_MAX; - double xmax = -DBL_MAX; - double ymin = +DBL_MAX; - double ymax = -DBL_MAX; - - - int pillars = (logical_gridsize_[0]+1) * (logical_gridsize_[1]+1); - - for (int pillarindex = 0; pillarindex < pillars; ++pillarindex) { - if (coord[pillarindex * 6 + 0] > xmax) - xmax = coord[pillarindex * 6 + 0]; - if (coord[pillarindex * 6 + 0] < xmin) - xmin = coord[pillarindex * 6 + 0]; - if (coord[pillarindex * 6 + 1] > ymax) - ymax = coord[pillarindex * 6 + 1]; - if (coord[pillarindex * 6 + 1] < ymin) - ymin = coord[pillarindex * 6 + 1]; - if (coord[pillarindex * 6 + 3] > xmax) - xmax = coord[pillarindex * 6 + 3]; - if (coord[pillarindex * 6 + 3] < xmin) - xmin = coord[pillarindex * 6 + 3]; - if (coord[pillarindex * 6 + 4] > ymax) - ymax = coord[pillarindex * 6 + 4]; - if (coord[pillarindex * 6 + 4] < ymin) - ymin = coord[pillarindex * 6 + 4]; - } - - std::array gridlimits = {{ xmin, xmax, ymin, ymax, - *min_element(zcorn.begin(), zcorn.end()), - *max_element(zcorn.begin(), zcorn.end()) }}; - return gridlimits; -} - - - -std::array EclipseGridInspector::gridSize() const -{ - std::array retval = {{ logical_gridsize_[0], - logical_gridsize_[1], - logical_gridsize_[2] }}; - return retval; -} - - -std::array EclipseGridInspector::cellZvals(int i, int j, int k) const -{ - // Get the zcorn field. - const std::vector& z = deck_->getKeyword("ZCORN").getSIDoubleData(); - int num_cells = logical_gridsize_[0]*logical_gridsize_[1]*logical_gridsize_[2]; - if (8*num_cells != int(z.size())) { - throw std::runtime_error("Wrong size of ZCORN field"); - } - - // Make the coordinate array. - int delta[3] = { 1, - 2*logical_gridsize_[0], - 4*logical_gridsize_[0]*logical_gridsize_[1] }; - int ix = 2*(i*delta[0] + j*delta[1] + k*delta[2]); - std::array cellz = {{ z[ix], z[ix + delta[0]], - z[ix + delta[1]], z[ix + delta[1] + delta[0]], - z[ix + delta[2]], z[ix + delta[2] + delta[0]], - z[ix + delta[2] + delta[1]], z[ix + delta[2] + delta[1] + delta[0]] }}; - return cellz; -} - - -} // namespace Opm diff --git a/opm/core/io/eclipse/EclipseGridInspector.hpp b/opm/core/io/eclipse/EclipseGridInspector.hpp deleted file mode 100644 index df0dba0c..00000000 --- a/opm/core/io/eclipse/EclipseGridInspector.hpp +++ /dev/null @@ -1,103 +0,0 @@ -//=========================================================================== -// -// File: EclipseGridInspector.h -// -// Created: Mon Jun 2 09:46:08 2008 -// -// Author: Atgeirr F Rasmussen -// -// $Date$ -// -// Revision: $Id: EclipseGridInspector.h,v 1.2 2008/08/18 14:16:12 atgeirr Exp $ -// -//=========================================================================== - -/* - Copyright 2009, 2010 SINTEF ICT, Applied Mathematics. - Copyright 2009, 2010 Statoil ASA. - - 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 . -*/ - -#ifndef OPM_ECLIPSEGRIDINSPECTOR_HEADER -#define OPM_ECLIPSEGRIDINSPECTOR_HEADER - -#include -#include - -#include - -namespace Opm -{ - -/** - @brief A class for inspecting the contents of an eclipse file. - - Given an Eclipse deck, this class may be used to answer certain - queries about its contents. - - @author Atgeirr F. Rasmussen - @date 2008/06/02 09:46:08 -*/ -class EclipseGridInspector -{ -public: - /// Constructor taking a parser as argument. - /// The parser must already have read an Eclipse file. - EclipseGridInspector(Opm::DeckConstPtr deck); - - /// Assuming that the pillars are vertical, compute the - /// volume of the cell given by logical coordinates (i, j, k). - double cellVolumeVerticalPillars(int i, int j, int k) const; - - /// Assuming that the pillars are vertical, compute the - /// volume of the cell given by the cell index - double cellVolumeVerticalPillars(int cell_idx) const; - - /// Compute the average dip in x- and y-direction of the - /// cell tops and bottoms relative to the xy-plane - std::pair cellDips(int i, int j, int k) const; - std::pair cellDips(int cell_idx) const; - - // Convert global cell index to logical ijk-coordinates - std::array cellIdxToLogicalCoords(int cell_idx) const; - - /// Returns a vector with the outer limits of grid (in the grid's unit). - /// The vector contains [xmin, xmax, ymin, ymax, zmin, zmax], as - /// read from COORDS and ZCORN - std::array getGridLimits() const; - - /// Returns the extent of the logical cartesian grid - /// as number of cells in the (i, j, k) directions. - std::array gridSize() const; - - /// Returns the eight z-values associated with a given cell. - /// The ordering is such that i runs fastest. That is, with - /// L = low and H = high: - /// {LLL, HLL, LHL, HHL, LLH, HLH, LHH, HHH }. - std::array cellZvals(int i, int j, int k) const; - -private: - Opm::DeckConstPtr deck_; - int logical_gridsize_[3]; - void init_(); - void checkLogicalCoords(int i, int j, int k) const; -}; - -} // namespace Opm - -#endif // OPM_ECLIPSEGRIDINSPECTOR_HEADER - diff --git a/opm/core/io/eclipse/EclipseIOUtil.hpp b/opm/core/io/eclipse/EclipseIOUtil.hpp deleted file mode 100644 index 1b8d3c20..00000000 --- a/opm/core/io/eclipse/EclipseIOUtil.hpp +++ /dev/null @@ -1,34 +0,0 @@ -#ifndef ECLIPSE_IO_UTIL_HPP -#define ECLIPSE_IO_UTIL_HPP - -#include -#include -#include - -namespace Opm -{ -namespace EclipseIOUtil -{ - - template - void addToStripedData(const std::vector& data, std::vector& result, size_t offset, size_t stride) { - int dataindex = 0; - for (size_t index = offset; index < result.size(); index += stride) { - result[index] = data[dataindex]; - ++dataindex; - } - } - - - template - void extractFromStripedData(const std::vector& data, std::vector& result, size_t offset, size_t stride) { - for (size_t index = offset; index < data.size(); index += stride) { - result.push_back(data[index]); - } - } - - -} //namespace EclipseIOUtil -} //namespace Opm - -#endif //ECLIPSE_IO_UTIL_HPP diff --git a/opm/core/io/eclipse/EclipseReader.cpp b/opm/core/io/eclipse/EclipseReader.cpp deleted file mode 100644 index ab7cf90b..00000000 --- a/opm/core/io/eclipse/EclipseReader.cpp +++ /dev/null @@ -1,247 +0,0 @@ -/* - Copyright 2015 Statoil ASA. - - 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 "EclipseReader.hpp" -#include -#include -#include -#include -#include - -#include -#include -#include -#include -#include - -#include - -#include - -namespace Opm -{ - - static void restoreTemperatureData(const ecl_file_type* file, - EclipseStateConstPtr eclipse_state, - int numcells, - SimulationDataContainer& simulator_state) { - const char* temperature = "TEMP"; - - if (ecl_file_has_kw(file , temperature)) { - ecl_kw_type* temperature_kw = ecl_file_iget_named_kw(file, temperature, 0); - - if (ecl_kw_get_size(temperature_kw) != numcells) { - throw std::runtime_error("Read of restart file: Could not restore temperature data, length of data from file not equal number of cells"); - } - - float* temperature_data = ecl_kw_get_float_ptr(temperature_kw); - - // factor and offset from the temperature values given in the deck to Kelvin - double scaling = eclipse_state->getDeckUnitSystem().parse("Temperature")->getSIScaling(); - double offset = eclipse_state->getDeckUnitSystem().parse("Temperature")->getSIOffset(); - - for (size_t index = 0; index < simulator_state.temperature().size(); ++index) { - simulator_state.temperature()[index] = unit::convert::from((double)temperature_data[index] - offset, scaling); - } - } else { - throw std::runtime_error("Read of restart file: File does not contain TEMP data\n"); - } - } - - - void restorePressureData(const ecl_file_type* file, - EclipseStateConstPtr eclipse_state, - int numcells, - SimulationDataContainer& simulator_state) { - const char* pressure = "PRESSURE"; - - if (ecl_file_has_kw(file , pressure)) { - - ecl_kw_type* pressure_kw = ecl_file_iget_named_kw(file, pressure, 0); - - if (ecl_kw_get_size(pressure_kw) != numcells) { - throw std::runtime_error("Read of restart file: Could not restore pressure data, length of data from file not equal number of cells"); - } - - float* pressure_data = ecl_kw_get_float_ptr(pressure_kw); - const double deck_pressure_unit = (eclipse_state->getDeckUnitSystem().getType() == UnitSystem::UNIT_TYPE_METRIC) ? Opm::unit::barsa : Opm::unit::psia; - for (size_t index = 0; index < simulator_state.pressure().size(); ++index) { - simulator_state.pressure()[index] = unit::convert::from((double)pressure_data[index], deck_pressure_unit); - } - } else { - throw std::runtime_error("Read of restart file: File does not contain PRESSURE data\n"); - } - } - - - static void restoreSaturation(const ecl_file_type* file_type, - const PhaseUsage& phaseUsage, - int numcells, - SimulationDataContainer& simulator_state) { - - float* sgas_data = NULL; - float* swat_data = NULL; - - if (phaseUsage.phase_used[BlackoilPhases::Aqua]) { - const char* swat = "SWAT"; - if (ecl_file_has_kw(file_type, swat)) { - ecl_kw_type* swat_kw = ecl_file_iget_named_kw(file_type , swat, 0); - swat_data = ecl_kw_get_float_ptr(swat_kw); - std::vector swat_datavec(&swat_data[0], &swat_data[numcells]); - EclipseIOUtil::addToStripedData(swat_datavec, simulator_state.saturation(), phaseUsage.phase_pos[BlackoilPhases::Aqua], phaseUsage.num_phases); - } else { - std::string error_str = "Restart file is missing SWAT data!\n"; - throw std::runtime_error(error_str); - } - } - - if (phaseUsage.phase_used[BlackoilPhases::Vapour]) { - const char* sgas = "SGAS"; - if (ecl_file_has_kw(file_type, sgas)) { - ecl_kw_type* sgas_kw = ecl_file_iget_named_kw(file_type , sgas, 0); - sgas_data = ecl_kw_get_float_ptr(sgas_kw); - std::vector sgas_datavec(&sgas_data[0], &sgas_data[numcells]); - EclipseIOUtil::addToStripedData(sgas_datavec, simulator_state.saturation(), phaseUsage.phase_pos[BlackoilPhases::Vapour], phaseUsage.num_phases); - } else { - std::string error_str = "Restart file is missing SGAS data!\n"; - throw std::runtime_error(error_str); - } - } - } - - - static void restoreRSandRV(const ecl_file_type* file_type, - SimulationConfigConstPtr sim_config, - int numcells, - SimulationDataContainer& simulator_state) { - - if (sim_config->hasDISGAS()) { - const char* RS = "RS"; - if (ecl_file_has_kw(file_type, RS)) { - ecl_kw_type* rs_kw = ecl_file_iget_named_kw(file_type, RS, 0); - float* rs_data = ecl_kw_get_float_ptr(rs_kw); - auto& rs = simulator_state.getCellData( BlackoilState::GASOILRATIO ); - for (int i = 0; i < ecl_kw_get_size( rs_kw ); i++) { - rs[i] = rs_data[i]; - } - } else { - throw std::runtime_error("Restart file is missing RS data!\n"); - } - } - - if (sim_config->hasVAPOIL()) { - const char* RV = "RV"; - if (ecl_file_has_kw(file_type, RV)) { - ecl_kw_type* rv_kw = ecl_file_iget_named_kw(file_type, RV, 0); - float* rv_data = ecl_kw_get_float_ptr(rv_kw); - auto& rv = simulator_state.getCellData( BlackoilState::RV ); - for (int i = 0; i < ecl_kw_get_size( rv_kw ); i++) { - rv[i] = rv_data[i]; - } - } else { - throw std::runtime_error("Restart file is missing RV data!\n"); - } - } - } - - - static void restoreSOLUTION(const std::string& restart_filename, - int reportstep, - bool unified, - EclipseStateConstPtr eclipseState, - int numcells, - const PhaseUsage& phaseUsage, - SimulationDataContainer& simulator_state) - { - const char* filename = restart_filename.c_str(); - ecl_file_type* file_type = ecl_file_open(filename, 0); - - if (file_type) { - bool block_selected = unified ? ecl_file_select_rstblock_report_step(file_type , reportstep) : true; - - if (block_selected) { - restorePressureData(file_type, eclipseState, numcells, simulator_state); - restoreTemperatureData(file_type, eclipseState, numcells, simulator_state); - restoreSaturation(file_type, phaseUsage, numcells, simulator_state); - if (simulator_state.hasCellData( BlackoilState::RV )) { - SimulationConfigConstPtr sim_config = eclipseState->getSimulationConfig(); - restoreRSandRV(file_type, sim_config, numcells, simulator_state ); - } - } else { - std::string error_str = "Restart file " + restart_filename + " does not contain data for report step " + std::to_string(reportstep) + "!\n"; - throw std::runtime_error(error_str); - } - ecl_file_close(file_type); - } else { - std::string error_str = "Restart file " + restart_filename + " not found!\n"; - throw std::runtime_error(error_str); - } - } - - - static void restoreOPM_XWELKeyword(const std::string& restart_filename, int reportstep, bool unified, WellState& wellstate) - { - const char * keyword = "OPM_XWEL"; - const char* filename = restart_filename.c_str(); - ecl_file_type* file_type = ecl_file_open(filename, 0); - - if (file_type != NULL) { - - bool block_selected = unified ? ecl_file_select_rstblock_report_step(file_type , reportstep) : true; - - if (block_selected) { - ecl_kw_type* xwel = ecl_file_iget_named_kw(file_type , keyword, 0); - const double* xwel_data = ecl_kw_get_double_ptr(xwel); - std::copy_n(xwel_data + wellstate.getRestartTemperatureOffset(), wellstate.temperature().size(), wellstate.temperature().begin()); - std::copy_n(xwel_data + wellstate.getRestartBhpOffset(), wellstate.bhp().size(), wellstate.bhp().begin()); - std::copy_n(xwel_data + wellstate.getRestartPerfPressOffset(), wellstate.perfPress().size(), wellstate.perfPress().begin()); - std::copy_n(xwel_data + wellstate.getRestartPerfRatesOffset(), wellstate.perfRates().size(), wellstate.perfRates().begin()); - std::copy_n(xwel_data + wellstate.getRestartWellRatesOffset(), wellstate.wellRates().size(), wellstate.wellRates().begin()); - } else { - std::string error_str = "Restart file " + restart_filename + " does not contain data for report step " + std::to_string(reportstep) + "!\n"; - throw std::runtime_error(error_str); - } - ecl_file_close(file_type); - } else { - std::string error_str = "Restart file " + restart_filename + " not found!\n"; - throw std::runtime_error(error_str); - } - } - - - - void init_from_restart_file(EclipseStateConstPtr eclipse_state, - int numcells, - const PhaseUsage& phase_usage, - SimulationDataContainer& simulator_state, - WellState& wellstate) { - - InitConfigConstPtr initConfig = eclipse_state->getInitConfig(); - IOConfigConstPtr ioConfig = eclipse_state->getIOConfig(); - int restart_step = initConfig->getRestartStep(); - const std::string& restart_file_root = initConfig->getRestartRootName(); - bool output = false; - const std::string& restart_file_name = ioConfig->getRestartFileName(restart_file_root, restart_step, output); - - Opm::restoreSOLUTION(restart_file_name, restart_step, ioConfig->getUNIFIN(), eclipse_state, numcells, phase_usage, simulator_state); - Opm::restoreOPM_XWELKeyword(restart_file_name, restart_step, ioConfig->getUNIFIN(), wellstate); - } - - -} // namespace Opm diff --git a/opm/core/io/eclipse/EclipseReader.hpp b/opm/core/io/eclipse/EclipseReader.hpp deleted file mode 100644 index 8db43572..00000000 --- a/opm/core/io/eclipse/EclipseReader.hpp +++ /dev/null @@ -1,39 +0,0 @@ -#ifndef ECLIPSEREADER_HPP -#define ECLIPSEREADER_HPP - -#include - -#include -#include -#include - - - -namespace Opm -{ -/// -/// \brief init_from_restart_file -/// Reading from the restart file, information stored under the OPM_XWEL keyword and SOLUTION data is in this method filled into -/// an instance of a wellstate object and a SimulatorState object. -/// \param grid -/// UnstructuredGrid reference -/// \param pu -/// PhaseUsage reference -/// \param simulator_state -/// An instance of a SimulatorState object -/// \param wellstate -/// An instance of a WellState object, with correct size for each of the 5 contained std::vector objects -/// - - class SimulationDataContainer; - - void init_from_restart_file(EclipseStateConstPtr eclipse_state, - int numcells, - const PhaseUsage& pu, - SimulationDataContainer& simulator_state, - WellState& wellstate); - - -} - -#endif // ECLIPSEREADER_HPP diff --git a/opm/core/io/eclipse/EclipseUnits.hpp b/opm/core/io/eclipse/EclipseUnits.hpp deleted file mode 100644 index 866ba0cd..00000000 --- a/opm/core/io/eclipse/EclipseUnits.hpp +++ /dev/null @@ -1,63 +0,0 @@ -/* - Copyright 2010 SINTEF ICT, Applied Mathematics. - - 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 . -*/ - -#ifndef OPM_ECLIPSEUNITS_HEADER_INCLUDED -#define OPM_ECLIPSEUNITS_HEADER_INCLUDED - -namespace Opm -{ - struct EclipseUnits - { - double length; - double time; - double density; - double polymer_density; - double pressure; - double compressibility; - double viscosity; - double permeability; - double liqvol_s; - double liqvol_r; - double gasvol_s; - double gasvol_r; - double transmissibility; - - void setToOne() - { - length = 1.0; - time = 1.0; - density = 1.0; - polymer_density = 1.0; - pressure = 1.0; - compressibility = 1.0; - viscosity = 1.0; - permeability = 1.0; - liqvol_s = 1.0; - liqvol_r = 1.0; - gasvol_s = 1.0; - gasvol_r = 1.0; - transmissibility = 1.0; - } - }; - - -} // namespace Opm - - -#endif // OPM_ECLIPSEUNITS_HEADER_INCLUDED diff --git a/opm/core/io/eclipse/EclipseWriteRFTHandler.cpp b/opm/core/io/eclipse/EclipseWriteRFTHandler.cpp deleted file mode 100644 index 093efe4c..00000000 --- a/opm/core/io/eclipse/EclipseWriteRFTHandler.cpp +++ /dev/null @@ -1,144 +0,0 @@ -/* - Copyright 2015 Statoil ASA. - - 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 - -#include -#include -#include -#include -#include - - -#include -#include -#include - -#include - - -#include -#include - - - -namespace Opm { -namespace EclipseWriterDetails { - - EclipseWriteRFTHandler::EclipseWriteRFTHandler(const int * compressedToCartesianCellIdx, size_t numCells, size_t cartesianSize) { - initGlobalToActiveIndex(compressedToCartesianCellIdx, numCells, cartesianSize); - } - - void EclipseWriteRFTHandler::writeTimeStep(const std::string& filename, - const ert_ecl_unit_enum ecl_unit, - const SimulatorTimerInterface& simulatorTimer, - std::vector& wells, - EclipseGridConstPtr eclipseGrid, - std::vector& pressure, - std::vector& swat, - std::vector& sgas) { - - - - std::vector rft_nodes; - for (std::vector::const_iterator ci = wells.begin(); ci != wells.end(); ++ci) { - WellConstPtr well = *ci; - if ((well->getRFTActive(simulatorTimer.reportStepNum())) || (well->getPLTActive(simulatorTimer.reportStepNum()))) { - ecl_rft_node_type * ecl_node = createEclRFTNode(well, - simulatorTimer, - eclipseGrid, - pressure, - swat, - sgas); - - // TODO: replace this silenced warning with an appropriate - // use of the OpmLog facilities. - // if (well->getPLTActive(simulatorTimer.reportStepNum())) { - // std::cerr << "PLT not supported, writing RFT data" << std::endl; - // } - - rft_nodes.push_back(ecl_node); - } - } - - - if (rft_nodes.size() > 0) { - ecl_rft_file_update(filename.c_str(), rft_nodes.data(), rft_nodes.size(), ecl_unit); - } - - //Cleanup: The ecl_rft_file_update method takes care of freeing the ecl_rft_nodes that it receives. - // Each ecl_rft_node is again responsible for freeing it's cells. - } - - - - - ecl_rft_node_type * EclipseWriteRFTHandler::createEclRFTNode(WellConstPtr well, - const SimulatorTimerInterface& simulatorTimer, - EclipseGridConstPtr eclipseGrid, - const std::vector& pressure, - const std::vector& swat, - const std::vector& sgas) { - - - const std::string& well_name = well->name(); - size_t timestep = (size_t)simulatorTimer.reportStepNum(); - time_t recording_date = simulatorTimer.currentPosixTime(); - double days = Opm::unit::convert::to(simulatorTimer.simulationTimeElapsed(), Opm::unit::day); - - std::string type = "RFT"; - ecl_rft_node_type * ecl_rft_node = ecl_rft_node_alloc_new(well_name.c_str(), type.c_str(), recording_date, days); - - CompletionSetConstPtr completionsSet = well->getCompletions(timestep); - for (size_t index = 0; index < completionsSet->size(); ++index) { - CompletionConstPtr completion = completionsSet->get(index); - size_t i = (size_t)completion->getI(); - size_t j = (size_t)completion->getJ(); - size_t k = (size_t)completion->getK(); - - size_t global_index = eclipseGrid->getGlobalIndex(i,j,k); - int active_index = globalToActiveIndex_[global_index]; - - if (active_index > -1) { - double depth = eclipseGrid->getCellDepth(i,j,k); - double completion_pressure = pressure.size() > 0 ? pressure[active_index] : 0.0; - double saturation_water = swat.size() > 0 ? swat[active_index] : 0.0; - double saturation_gas = sgas.size() > 0 ? sgas[active_index] : 0.0; - - ecl_rft_cell_type * ecl_rft_cell = ecl_rft_cell_alloc_RFT( i ,j, k , depth, completion_pressure, saturation_water, saturation_gas); - ecl_rft_node_append_cell( ecl_rft_node , ecl_rft_cell); - } - } - - return ecl_rft_node; - } - - - void EclipseWriteRFTHandler::initGlobalToActiveIndex(const int * compressedToCartesianCellIdx, size_t numCells, size_t cartesianSize) { - globalToActiveIndex_.resize(cartesianSize, -1); - for (size_t active_index = 0; active_index < numCells; ++active_index) { - //If compressedToCartesianCellIdx is NULL, assume no compressed to cartesian mapping, set global equal to active index - int global_index = (NULL != compressedToCartesianCellIdx) ? compressedToCartesianCellIdx[active_index] : active_index; - globalToActiveIndex_[global_index] = active_index; - } - } - - -}//namespace EclipseWriterDetails -}//namespace Opm diff --git a/opm/core/io/eclipse/EclipseWriteRFTHandler.hpp b/opm/core/io/eclipse/EclipseWriteRFTHandler.hpp deleted file mode 100644 index ea178746..00000000 --- a/opm/core/io/eclipse/EclipseWriteRFTHandler.hpp +++ /dev/null @@ -1,77 +0,0 @@ -/* - Copyright 2015 Statoil ASA. - - 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 . -*/ - -#ifndef OPM_ECLIPSE_WRITE_RFT_HANDLER_HPP -#define OPM_ECLIPSE_WRITE_RFT_HANDLER_HPP - -#include - -#include - -#include -#include - - -namespace Opm { - class EclipseGrid; - class Well; - -namespace EclipseWriterDetails { - - - class EclipseWriteRFTHandler { - - public: - EclipseWriteRFTHandler(const int * compressedToCartesianCellIdx, size_t numCells, size_t cartesianSize); - - - void writeTimeStep(const std::string& filename, - const ert_ecl_unit_enum ecl_unit, - const SimulatorTimerInterface& simulatorTimer, - std::vector>& wells, - std::shared_ptr< const EclipseGrid > eclipseGrid, - std::vector& pressure, - std::vector& swat, - std::vector& sgas); - - - - private: - - ecl_rft_node_type * createEclRFTNode(std::shared_ptr< const Well > well, - const SimulatorTimerInterface& simulatorTimer, - std::shared_ptr< const EclipseGrid > eclipseGrid, - const std::vector& pressure, - const std::vector& swat, - const std::vector& sgas); - - void initGlobalToActiveIndex(const int * compressedToCartesianCellIdx, size_t numCells, size_t cartesianSize); - - std::vector globalToActiveIndex_; - - }; - - - - -}//namespace EclipseWriterDetails -}//namespace Opm - - -#endif //OPM_ECLIPSE_WRITE_RFT_HANDLER_HPP diff --git a/opm/core/io/eclipse/EclipseWriter.cpp b/opm/core/io/eclipse/EclipseWriter.cpp deleted file mode 100644 index 75909acb..00000000 --- a/opm/core/io/eclipse/EclipseWriter.cpp +++ /dev/null @@ -1,1508 +0,0 @@ -/* - Copyright (c) 2013-2015 Andreas Lauser - Copyright (c) 2013 SINTEF ICT, Applied Mathematics. - Copyright (c) 2013 Uni Research AS - Copyright (c) 2015 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 "EclipseWriter.hpp" - -#include - - -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include // WellType - -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#include // to_upper_copy -#include -#include // path - -#include // mktime -#include -#include // unique_ptr -#include // move - -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#define OPM_XWEL "OPM_XWEL" - -// namespace start here since we don't want the ERT headers in it -namespace Opm { -namespace EclipseWriterDetails { -/// Names of the saturation property for each phase. The order of these -/// names are critical; they must be the same as the BlackoilPhases enum -static const char* saturationKeywordNames[] = { "SWAT", "SOIL", "SGAS" }; - -// throw away the data for all non-active cells and reorder to the Cartesian logic of -// eclipse -void restrictAndReorderToActiveCells(std::vector &data, - int numCells, - const int* compressedToCartesianCellIdx) -{ - if (!compressedToCartesianCellIdx) - // if there is no active -> global mapping, all cells - // are considered active - return; - - std::vector eclData; - eclData.reserve( numCells ); - - // activate those cells that are actually there - for (int i = 0; i < numCells; ++i) { - eclData.push_back( data[ compressedToCartesianCellIdx[i] ] ); - } - data.swap( eclData ); -} - -// convert the units of an array -void convertFromSiTo(std::vector &siValues, double toSiConversionFactor, double toSiOffset = 0) -{ - for (size_t curIdx = 0; curIdx < siValues.size(); ++curIdx) { - siValues[curIdx] = unit::convert::to(siValues[curIdx] - toSiOffset, toSiConversionFactor); - } -} - -// extract a sub-array of a larger one which represents multiple -// striped ones -void extractFromStripedData(std::vector &data, - int offset, - int stride) -{ - size_t tmpIdx = 0; - for (size_t curIdx = offset; curIdx < data.size(); curIdx += stride) { - assert(tmpIdx <= curIdx); - data[tmpIdx] = data[curIdx]; - ++tmpIdx; - } - // shirk the result - data.resize(tmpIdx); -} - -/// Convert OPM phase usage to ERT bitmask -int ertPhaseMask(const PhaseUsage uses) -{ - return (uses.phase_used[BlackoilPhases::Liquid] ? ECL_OIL_PHASE : 0) - | (uses.phase_used[BlackoilPhases::Aqua] ? ECL_WATER_PHASE : 0) - | (uses.phase_used[BlackoilPhases::Vapour] ? ECL_GAS_PHASE : 0); -} - - -/** - * Eclipse "keyword" (i.e. named data) for a vector. - */ -template -class Keyword : private boost::noncopyable -{ -public: - // Default constructor - Keyword() - : ertHandle_(0) - {} - - /// Initialization from double-precision array. - Keyword(const std::string& name, - const std::vector& data) - : ertHandle_(0) - { set(name, data); } - - /// Initialization from double-precision array. - Keyword(const std::string& name, - const std::vector& data) - : ertHandle_(0) - { set(name, data); } - - Keyword(const std::string& name, - const std::vector& data) - : ertHandle_(0) - {set(name, data); } - - ~Keyword() - { - if (ertHandle_) - ecl_kw_free(ertHandle_); - } - - template - void set(const std::string name, const std::vector& data) - { - - if(ertHandle_) { - ecl_kw_free(ertHandle_); - } - - - ertHandle_ = ecl_kw_alloc(name.c_str(), - data.size(), - ertType_()); - - // number of elements to take - const int numEntries = data.size(); - - // fill it with values - - T* target = static_cast(ecl_kw_get_ptr(ertHandle())); - for (int i = 0; i < numEntries; ++i) { - target[i] = static_cast(data[i]); - } - } - - void set(const std::string name, const std::vector& data) - { - if(ertHandle_) { - ecl_kw_free(ertHandle_); - } - - - ertHandle_ = ecl_kw_alloc(name.c_str(), - data.size(), - ertType_()); - - // number of elements to take - const int numEntries = data.size(); - for (int i = 0; i < numEntries; ++i) { - ecl_kw_iset_char_ptr( ertHandle_, i, data[i]); - } - - } - - ecl_kw_type *ertHandle() const - { return ertHandle_; } - -private: - static ecl_type_enum ertType_() - { - if (std::is_same::value) - { return ECL_FLOAT_TYPE; } - if (std::is_same::value) - { return ECL_DOUBLE_TYPE; } - if (std::is_same::value) - { return ECL_INT_TYPE; } - if (std::is_same::value) - { return ECL_CHAR_TYPE; } - - - OPM_THROW(std::logic_error, - "Unhandled type for data elements in EclipseWriterDetails::Keyword"); - } - - ecl_kw_type *ertHandle_; -}; - -/** - * Pointer to memory that holds the name to an Eclipse output file. - */ -class FileName : private boost::noncopyable -{ -public: - FileName(const std::string& outputDir, - const std::string& baseName, - ecl_file_enum type, - int writeStepIdx, - bool formatted) - { - ertHandle_ = ecl_util_alloc_filename(outputDir.c_str(), - baseName.c_str(), - type, - formatted, - writeStepIdx); - } - - ~FileName() - { std::free(ertHandle_); } - - const char *ertHandle() const - { return ertHandle_; } - -private: - char *ertHandle_; -}; - -class Restart : private boost::noncopyable -{ -public: - static const int NIWELZ = 11; //Number of data elements per well in IWEL array in restart file - static const int NZWELZ = 3; //Number of 8-character words per well in ZWEL array restart file - static const int NICONZ = 15; //Number of data elements per completion in ICON array restart file - - /** - * The constants NIWELZ and NZWELZ referes to the number of elements per - * well that we write to the IWEL and ZWEL eclipse restart file data - * arrays. The constant NICONZ refers to the number of elements per - * completion in the eclipse restart file ICON data array.These numbers are - * written to the INTEHEAD header. - - * The elements are added in the method addRestartFileIwelData(...) and and - * addRestartFileIconData(...), respectively. We write as many elements - * that we need to be able to view the restart file in Resinsight. The - * restart file will not be possible to restart from with Eclipse, we write - * to little information to be able to do this. - * - * Observe that all of these values are our "current-best-guess" for how - * many numbers are needed; there might very well be third party - * applications out there which have a hard expectation for these values. - */ - - - Restart(const std::string& outputDir, - const std::string& baseName, - int writeStepIdx, - IOConfigConstPtr ioConfig) - { - - ecl_file_enum type_of_restart_file = ioConfig->getUNIFOUT() ? ECL_UNIFIED_RESTART_FILE : ECL_RESTART_FILE; - - restartFileName_ = ecl_util_alloc_filename(outputDir.c_str(), - baseName.c_str(), - type_of_restart_file, - ioConfig->getFMTOUT(), // use formatted instead of binary output? - writeStepIdx); - - if ((writeStepIdx > 0) && (ECL_UNIFIED_RESTART_FILE == type_of_restart_file)) { - restartFileHandle_ = ecl_rst_file_open_append(restartFileName_); - } - else { - restartFileHandle_ = ecl_rst_file_open_write(restartFileName_); - } - } - - template - void add_kw(const Keyword& kw) - { ecl_rst_file_add_kw(restartFileHandle_, kw.ertHandle()); } - - - void addRestartFileIwelData(std::vector& iwel_data, size_t currentStep , WellConstPtr well, size_t offset) const { - CompletionSetConstPtr completions = well->getCompletions( currentStep ); - - iwel_data[ offset + IWEL_HEADI_ITEM ] = well->getHeadI() + 1; - iwel_data[ offset + IWEL_HEADJ_ITEM ] = well->getHeadJ() + 1; - iwel_data[ offset + IWEL_CONNECTIONS_ITEM ] = completions->size(); - iwel_data[ offset + IWEL_GROUP_ITEM ] = 1; - - { - WellType welltype = well->isProducer(currentStep) ? PRODUCER : INJECTOR; - int ert_welltype = EclipseWriter::eclipseWellTypeMask(welltype, well->getInjectionProperties(currentStep).injectorType); - iwel_data[ offset + IWEL_TYPE_ITEM ] = ert_welltype; - } - - iwel_data[ offset + IWEL_STATUS_ITEM ] = EclipseWriter::eclipseWellStatusMask(well->getStatus(currentStep)); - } - - void addRestartFileXwelData(const WellState& wellState, std::vector& xwel_data) const { - std::copy_n(wellState.bhp().begin(), wellState.bhp().size(), xwel_data.begin() + wellState.getRestartBhpOffset()); - std::copy_n(wellState.perfPress().begin(), wellState.perfPress().size(), xwel_data.begin() + wellState.getRestartPerfPressOffset()); - std::copy_n(wellState.perfRates().begin(), wellState.perfRates().size(), xwel_data.begin() + wellState.getRestartPerfRatesOffset()); - std::copy_n(wellState.temperature().begin(), wellState.temperature().size(), xwel_data.begin() + wellState.getRestartTemperatureOffset()); - std::copy_n(wellState.wellRates().begin(), wellState.wellRates().size(), xwel_data.begin() + wellState.getRestartWellRatesOffset()); - } - - void addRestartFileIconData(std::vector& icon_data, CompletionSetConstPtr completions , size_t wellICONOffset) const { - for (size_t i = 0; i < completions->size(); ++i) { - CompletionConstPtr completion = completions->get(i); - size_t iconOffset = wellICONOffset + i * Opm::EclipseWriterDetails::Restart::NICONZ; - icon_data[ iconOffset + ICON_IC_ITEM] = 1; - - icon_data[ iconOffset + ICON_I_ITEM] = completion->getI() + 1; - icon_data[ iconOffset + ICON_J_ITEM] = completion->getJ() + 1; - icon_data[ iconOffset + ICON_K_ITEM] = completion->getK() + 1; - - { - WellCompletion::StateEnum completion_state = completion->getState(); - if (completion_state == WellCompletion::StateEnum::OPEN) { - icon_data[ iconOffset + ICON_STATUS_ITEM ] = 1; - } else { - icon_data[ iconOffset + ICON_STATUS_ITEM ] = 0; - } - } - icon_data[ iconOffset + ICON_DIRECTION_ITEM] = (int)completion->getDirection(); - } - } - - - ~Restart() - { - free(restartFileName_); - ecl_rst_file_close(restartFileHandle_); - } - - void writeHeader(const SimulatorTimerInterface& /*timer*/, - int writeStepIdx, - ecl_rsthead_type * rsthead_data) - { - - ecl_util_set_date_values( rsthead_data->sim_time , &rsthead_data->day , &rsthead_data->month , &rsthead_data->year ); - - ecl_rst_file_fwrite_header(restartFileHandle_, - writeStepIdx, - rsthead_data); - - } - - ecl_rst_file_type *ertHandle() const - { return restartFileHandle_; } - -private: - char *restartFileName_; - ecl_rst_file_type *restartFileHandle_; -}; - -/** - * The Solution class wraps the actions that must be done to the restart file while - * writing solution variables; it is not a handle on its own. - */ -class Solution : private boost::noncopyable -{ -public: - Solution(Restart& restartHandle) - : restartHandle_(&restartHandle) - { ecl_rst_file_start_solution(restartHandle_->ertHandle()); } - - ~Solution() - { ecl_rst_file_end_solution(restartHandle_->ertHandle()); } - - template - void add(const Keyword& kw) - { ecl_rst_file_add_kw(restartHandle_->ertHandle(), kw.ertHandle()); } - - ecl_rst_file_type *ertHandle() const - { return restartHandle_->ertHandle(); } - -private: - Restart* restartHandle_; -}; - -/// Supported well types. Enumeration doesn't let us get all the members, -/// so we must have an explicit array. -static WellType WELL_TYPES[] = { INJECTOR, PRODUCER }; - -class WellReport; - -class Summary : private boost::noncopyable -{ -public: - Summary(const std::string& outputDir, - const std::string& baseName, - const SimulatorTimerInterface& timer, - bool time_in_days , - int nx, - int ny, - int nz) - { - boost::filesystem::path casePath(outputDir); - casePath /= boost::to_upper_copy(baseName); - - // convert the start time to seconds since 1970-1-1@00:00:00 - boost::posix_time::ptime startTime - = timer.startDateTime(); - tm t = boost::posix_time::to_tm(startTime); - double secondsSinceEpochStart = std::mktime(&t); - - ertHandle_ = ecl_sum_alloc_writer(casePath.string().c_str(), - false, /* formatted */ - true, /* unified */ - ":", /* join string */ - secondsSinceEpochStart, - time_in_days, - nx, - ny, - nz); - } - - ~Summary() - { ecl_sum_free(ertHandle_); } - - typedef std::unique_ptr SummaryReportVar; - typedef std::vector SummaryReportVarCollection; - - Summary& addWell(SummaryReportVar var) - { - summaryReportVars_.push_back(std::move(var)); - return *this; - } - - // no inline implementation of these two methods since they depend - // on the classes defined in the following. - - // add rate variables for each of the well in the input file - void addAllWells(Opm::EclipseStateConstPtr eclipseState, - const PhaseUsage& uses); - void writeTimeStep(int writeStepIdx, - const SimulatorTimerInterface& timer, - const WellState& wellState); - - ecl_sum_type *ertHandle() const - { return ertHandle_; } - -private: - ecl_sum_type *ertHandle_; - - Opm::EclipseStateConstPtr eclipseState_; - SummaryReportVarCollection summaryReportVars_; -}; - -class SummaryTimeStep : private boost::noncopyable -{ -public: - SummaryTimeStep(Summary& summaryHandle, - int writeStepIdx, - const SimulatorTimerInterface &timer) - { - ertHandle_ = ecl_sum_add_tstep(summaryHandle.ertHandle(), - writeStepIdx, - timer.simulationTimeElapsed()); - } - - // no destructor in this class as ERT takes care of freeing the - // handle as part of freeing the solution handle! - - ecl_sum_tstep_type *ertHandle() const - { return ertHandle_; }; - -private: - ecl_sum_tstep_type *ertHandle_; -}; - - -/** - * Initialization file which contains static properties (such as - * porosity and permeability) for the simulation field. - */ -class Init : private boost::noncopyable -{ -public: - Init(const std::string& outputDir, - const std::string& baseName, - int writeStepIdx, - IOConfigConstPtr ioConfig) : egridFileName_(outputDir, - baseName, - ECL_EGRID_FILE, - writeStepIdx, - ioConfig->getFMTOUT()) - { - bool formatted = ioConfig->getFMTOUT(); - - FileName initFileName(outputDir, - baseName, - ECL_INIT_FILE, - writeStepIdx, - formatted); - - ertHandle_ = fortio_open_writer(initFileName.ertHandle(), - formatted, - ECL_ENDIAN_FLIP); - } - - ~Init() - { fortio_fclose(ertHandle_); } - - void writeHeader(int numCells, - const int* compressedToCartesianCellIdx, - const SimulatorTimerInterface& timer, - Opm::EclipseStateConstPtr eclipseState, - const PhaseUsage uses) - { - auto dataField = eclipseState->getDoubleGridProperty("PORO")->getData(); - restrictAndReorderToActiveCells(dataField, numCells, compressedToCartesianCellIdx); - - auto eclGrid = eclipseState->getEclipseGridCopy(); - - // update the ACTNUM array using the processed cornerpoint grid - std::vector actnumData(eclGrid->getCartesianSize(), 1); - if (compressedToCartesianCellIdx) { - std::fill(actnumData.begin(), actnumData.end(), 0); - for (int cellIdx = 0; cellIdx < numCells; ++cellIdx) { - int cartesianCellIdx = compressedToCartesianCellIdx[cellIdx]; - actnumData[cartesianCellIdx] = 1; - } - } - - eclGrid->resetACTNUM(&actnumData[0]); - - - // finally, write the grid to disk - IOConfigConstPtr ioConfig = eclipseState->getIOConfigConst(); - if (ioConfig->getWriteEGRIDFile()) { - if (eclipseState->getDeckUnitSystem().getType() == UnitSystem::UNIT_TYPE_METRIC){ - eclGrid->fwriteEGRID(egridFileName_.ertHandle(), true); - }else{ - eclGrid->fwriteEGRID(egridFileName_.ertHandle(), false); - } - } - - - if (ioConfig->getWriteINITFile()) { - Keyword poro_kw("PORO", dataField); - ecl_init_file_fwrite_header(ertHandle(), - eclGrid->c_ptr(), - poro_kw.ertHandle(), - ertPhaseMask(uses), - timer.currentPosixTime()); - } - } - - void writeKeyword(const std::string& keywordName, const std::vector &data) - { - Keyword kw(keywordName, data); - ecl_kw_fwrite(kw.ertHandle(), ertHandle()); - } - - fortio_type *ertHandle() const - { return ertHandle_; } - -private: - fortio_type *ertHandle_; - FileName egridFileName_; -}; - - - - -/** - * Summary variable that reports a characteristics of a well. - */ -class WellReport : private boost::noncopyable -{ -protected: - WellReport(const Summary& summary, /* section to add to */ - Opm::EclipseStateConstPtr eclipseState, - Opm::WellConstPtr& well, - PhaseUsage uses, /* phases present */ - BlackoilPhases::PhaseIndex phase, /* oil, water or gas */ - WellType type, /* prod. or inj. */ - char aggregation, /* rate or total */ - std::string unit) - // save these for when we update the value in a timestep - : eclipseState_(eclipseState) - , well_(well) - , phaseUses_(uses) - , phaseIdx_(phase) - { - // producers can be seen as negative injectors - if (type == INJECTOR) - sign_ = +1.0; - else - sign_ = -1.0; - ertHandle_ = ecl_sum_add_var(summary.ertHandle(), - varName_(phase, - type, - aggregation).c_str(), - well_->name().c_str(), - /*num=*/ 0, - unit.c_str(), - /*defaultValue=*/ 0.); - } - -public: - /// Retrieve the value which the monitor is supposed to write to the summary file - /// according to the state of the well. - virtual double retrieveValue(const int writeStepIdx, - const SimulatorTimerInterface& timer, - const WellState& wellState, - const std::map& nameToIdxMap) = 0; - - smspec_node_type *ertHandle() const - { return ertHandle_; } - -protected: - - - void updateTimeStepWellIndex_(const std::map& nameToIdxMap) - { - const std::string& wellName = well_->name(); - - const auto wellIdxIt = nameToIdxMap.find(wellName); - if (wellIdxIt == nameToIdxMap.end()) { - timeStepWellIdx_ = -1; - flatIdx_ = -1; - return; - } - - timeStepWellIdx_ = wellIdxIt->second; - flatIdx_ = timeStepWellIdx_*phaseUses_.num_phases + phaseUses_.phase_pos[phaseIdx_]; - } - - // return m^3/s of injected or produced fluid - double rate(const WellState& wellState) - { - double value = 0; - if (wellState.wellRates().size() > 0) { - assert(int(wellState.wellRates().size()) > flatIdx_); - value = sign_ * wellState.wellRates()[flatIdx_]; - } - return value; - } - - double bhp(const WellState& wellState) - { - if (wellState.bhp().size() > 0) { - // Note that 'flatIdx_' is used here even though it is meant - // to give a (well,phase) pair. - const int numPhases = wellState.wellRates().size() / wellState.bhp().size(); - - return wellState.bhp()[flatIdx_/numPhases]; - } - return 0.0; - } - - /// Get the index associated a well name - int wellIndex_(Opm::EclipseStateConstPtr eclipseState) - { - const Opm::ScheduleConstPtr schedule = eclipseState->getSchedule(); - - const std::string& wellName = well_->name(); - const auto& wells = schedule->getWells(); - for (size_t wellIdx = 0; wellIdx < wells.size(); ++wellIdx) { - if (wells[wellIdx]->name() == wellName) { - return wellIdx; - } - } - - OPM_THROW(std::runtime_error, - "Well '" << wellName << "' is not present in deck"); - } - - /// Compose the name of the summary variable, e.g. "WOPR" for - /// well oil production rate. - std::string varName_(BlackoilPhases::PhaseIndex phase, - WellType type, - char aggregation) - { - std::string name; - name += 'W'; // well - if (aggregation == 'B') { - name += "BHP"; - } else { - switch (phase) { - case BlackoilPhases::Aqua: name += 'W'; break; /* water */ - case BlackoilPhases::Vapour: name += 'G'; break; /* gas */ - case BlackoilPhases::Liquid: name += 'O'; break; /* oil */ - default: - OPM_THROW(std::runtime_error, - "Unknown phase used in blackoil reporting"); - } - switch (type) { - case WellType::INJECTOR: name += 'I'; break; - case WellType::PRODUCER: name += 'P'; break; - default: - OPM_THROW(std::runtime_error, - "Unknown well type used in blackoil reporting"); - } - name += aggregation; /* rate ('R') or total ('T') */ - } - return name; - } - - smspec_node_type *ertHandle_; - - Opm::EclipseStateConstPtr eclipseState_; - Opm::WellConstPtr well_; - - PhaseUsage phaseUses_; - BlackoilPhases::PhaseIndex phaseIdx_; - - int timeStepWellIdx_; - - /// index into a (flattened) wellsOfTimeStep*phases matrix - int flatIdx_; - - /// natural sign of the rate - double sign_; -}; - -/// Monitors the rate given by a well. -class WellRate : public WellReport -{ -public: - WellRate(const Summary& summary, - Opm::EclipseStateConstPtr eclipseState, - Opm::WellConstPtr well, - PhaseUsage uses, - BlackoilPhases::PhaseIndex phase, - WellType type, - UnitSystem::UnitType unitType) - : WellReport(summary, - eclipseState, - well, - uses, - phase, - type, - 'R', - handleUnit_(phase, unitType)) - { - } - - virtual double retrieveValue(const int /* writeStepIdx */, - const SimulatorTimerInterface& timer, - const WellState& wellState, - const std::map& wellNameToIdxMap) - { - // find the index for the quantity in the wellState - this->updateTimeStepWellIndex_(wellNameToIdxMap); - if (this->flatIdx_ < 0) { - // well not active in current time step - return 0.0; - } - - if (well_->getStatus(timer.reportStepNum()) == WellCommon::SHUT) { - // well is shut in the current time step - return 0.0; - } - - // TODO: Why only positive rates? - using namespace Opm::unit; - return convert::to(std::max(0., rate(wellState)), - targetRateToSiConversionFactor_); - } - -private: - const std::string handleUnit_(BlackoilPhases::PhaseIndex phase, UnitSystem::UnitType unitType) { - using namespace Opm::unit; - if (phase == BlackoilPhases::Liquid || phase == BlackoilPhases::Aqua) { - if (unitType == UnitSystem::UNIT_TYPE_FIELD) { - unitName_ = "STB/DAY"; - targetRateToSiConversionFactor_ = stb/day; // m^3/s -> STB/day - } - else if (unitType == UnitSystem::UNIT_TYPE_METRIC) { - unitName_ = "SM3/DAY"; - targetRateToSiConversionFactor_ = cubic(meter)/day; // m^3/s -> m^3/day - } - else - OPM_THROW(std::logic_error, "Deck uses unexpected unit system"); - } - else if (phase == BlackoilPhases::Vapour) { - if (unitType == UnitSystem::UNIT_TYPE_FIELD) { - unitName_ = "MSCF/DAY"; - targetRateToSiConversionFactor_ = 1000*cubic(feet)/day; // m^3/s -> MSCF^3/day - } - else if (unitType == UnitSystem::UNIT_TYPE_METRIC) { - unitName_ = "SM3/DAY"; - targetRateToSiConversionFactor_ = cubic(meter)/day; // m^3/s -> m^3/day - } - else - OPM_THROW(std::logic_error, "Deck uses unexpected unit system"); - } - else - OPM_THROW(std::logic_error, - "Unexpected phase " << phase); - return unitName_; - } - - const char* unitName_; - double targetRateToSiConversionFactor_; -}; - -/// Monitors the total production in a well. -class WellTotal : public WellReport -{ -public: - WellTotal(const Summary& summary, - Opm::EclipseStateConstPtr eclipseState, - Opm::WellConstPtr well, - PhaseUsage uses, - BlackoilPhases::PhaseIndex phase, - WellType type, - UnitSystem::UnitType unitType) - : WellReport(summary, - eclipseState, - well, - uses, - phase, - type, - 'T', - handleUnit_(phase, unitType)) - // nothing produced when the reporting starts - , total_(0.) - { } - - virtual double retrieveValue(const int writeStepIdx, - const SimulatorTimerInterface& timer, - const WellState& wellState, - const std::map& wellNameToIdxMap) - { - if (writeStepIdx == 0) { - // We are at the initial state. - // No step has been taken yet. - return 0.0; - } - - if (well_->getStatus(timer.reportStepNum()) == WellCommon::SHUT) { - // well is shut in the current time step - return 0.0; - } - - // find the index for the quantity in the wellState - this->updateTimeStepWellIndex_(wellNameToIdxMap); - if (this->flatIdx_ < 0) { - // well not active in current time step - return 0.0; - } - - // due to using an Euler method as time integration scheme, the well rate is the - // average for the time step. For more complicated time stepping schemes, the - // integral of the rate is not simply multiplying two numbers... - const double intg = timer.stepLengthTaken() * rate(wellState); - - // add this timesteps production to the total - total_ += intg; - // report the new production total - return unit::convert::to(total_, targetRateToSiConversionFactor_); - } - -private: - const std::string handleUnit_(BlackoilPhases::PhaseIndex phase, UnitSystem::UnitType unitType) { - using namespace Opm::unit; - if (phase == BlackoilPhases::Liquid || phase == BlackoilPhases::Aqua) { - if (unitType == UnitSystem::UNIT_TYPE_FIELD) { - unitName_ = "STB"; - targetRateToSiConversionFactor_ = stb; // m^3 -> STB - } - else if (unitType == UnitSystem::UNIT_TYPE_METRIC) { - unitName_ = "SM3"; - targetRateToSiConversionFactor_ = cubic(meter); // m^3 -> m^3 - } - else - OPM_THROW(std::logic_error, "Deck uses unexpected unit system"); - } - else if (phase == BlackoilPhases::Vapour) { - if (unitType == UnitSystem::UNIT_TYPE_FIELD) { - unitName_ = "MSCF"; - targetRateToSiConversionFactor_ = 1000*cubic(feet); // m^3 -> MSCF^3 - } - else if (unitType == UnitSystem::UNIT_TYPE_METRIC) { - unitName_ = "SM3"; - targetRateToSiConversionFactor_ = cubic(meter); // m^3 -> m^3 - } - else - OPM_THROW(std::logic_error, "Deck uses unexpected unit system"); - } - else - OPM_THROW(std::logic_error, - "Unexpected phase " << phase); - return unitName_; - } - - const char* unitName_; - double targetRateToSiConversionFactor_; - - /// Aggregated value of the course of the simulation - double total_; -}; - -/// Monitors the bottom hole pressure in a well. -class WellBhp : public WellReport -{ -public: - WellBhp(const Summary& summary, - Opm::EclipseStateConstPtr eclipseState, - Opm::WellConstPtr well, - PhaseUsage uses, - BlackoilPhases::PhaseIndex phase, - WellType type, - UnitSystem::UnitType unitType) - : WellReport(summary, - eclipseState, - well, - uses, - phase, - type, - 'B', - handleUnit_(unitType)) - { } - - virtual double retrieveValue(const int /* writeStepIdx */, - const SimulatorTimerInterface& timer, - const WellState& wellState, - const std::map& wellNameToIdxMap) - { - // find the index for the quantity in the wellState - this->updateTimeStepWellIndex_(wellNameToIdxMap); - if (this->flatIdx_ < 0) { - // well not active in current time step - return 0.0; - } - if (well_->getStatus(timer.reportStepNum()) == WellCommon::SHUT) { - // well is shut in the current time step - return 0.0; - } - - return unit::convert::to(bhp(wellState), targetRateToSiConversionFactor_); - } - -private: - const std::string handleUnit_(UnitSystem::UnitType unitType) { - using namespace Opm::unit; - - if (unitType == UnitSystem::UNIT_TYPE_FIELD) { - unitName_ = "PSIA"; - targetRateToSiConversionFactor_ = psia; // Pa -> PSI - } - else if (unitType == UnitSystem::UNIT_TYPE_METRIC) { - unitName_ = "BARSA"; - targetRateToSiConversionFactor_ = barsa; // Pa -> bar - } - else - OPM_THROW(std::logic_error, - "Unexpected unit type " << unitType); - - return unitName_; - } - - const char* unitName_; - double targetRateToSiConversionFactor_; -}; - -// no inline implementation of this since it depends on the -// WellReport type being completed first -void Summary::writeTimeStep(int writeStepIdx, - const SimulatorTimerInterface& timer, - const WellState& wellState) -{ - // create a name -> well index map - const Opm::ScheduleConstPtr schedule = eclipseState_->getSchedule(); - const auto& timeStepWells = schedule->getWells(timer.reportStepNum()); - std::map wellNameToIdxMap; - int openWellIdx = 0; - for (size_t tsWellIdx = 0; tsWellIdx < timeStepWells.size(); ++tsWellIdx) { - if (timeStepWells[tsWellIdx]->getStatus(timer.reportStepNum()) != WellCommon::SHUT ) { - wellNameToIdxMap[timeStepWells[tsWellIdx]->name()] = openWellIdx; - openWellIdx++; - } - } - - // internal view; do not move this code out of Summary! - SummaryTimeStep tstep(*this, writeStepIdx, timer); - // write all the variables - for (auto varIt = summaryReportVars_.begin(); varIt != summaryReportVars_.end(); ++varIt) { - ecl_sum_tstep_iset(tstep.ertHandle(), - smspec_node_get_params_index((*varIt)->ertHandle()), - (*varIt)->retrieveValue(writeStepIdx, timer, wellState, wellNameToIdxMap)); - } - - // write the summary file to disk - ecl_sum_fwrite(ertHandle()); -} - -void Summary::addAllWells(Opm::EclipseStateConstPtr eclipseState, - const PhaseUsage& uses) -{ - eclipseState_ = eclipseState; - auto deckUnitType = eclipseState_->getDeckUnitSystem().getType(); - - // TODO: Only create report variables that are requested with keywords - // (e.g. "WOPR") in the input files, and only for those wells that are - // mentioned in those keywords - Opm::ScheduleConstPtr schedule = eclipseState->getSchedule(); - const auto& wells = schedule->getWells(); - const int numWells = schedule->numWells(); - for (int phaseIdx = 0; phaseIdx != BlackoilPhases::MaxNumPhases; ++phaseIdx) { - const BlackoilPhases::PhaseIndex ertPhaseIdx = - static_cast (phaseIdx); - // don't bother with reporting for phases that aren't there - if (!uses.phase_used[phaseIdx]) { - continue; - } - size_t numWellTypes = sizeof(WELL_TYPES) / sizeof(WELL_TYPES[0]); - for (size_t wellTypeIdx = 0; wellTypeIdx < numWellTypes; ++wellTypeIdx) { - const WellType wellType = WELL_TYPES[wellTypeIdx]; - for (int wellIdx = 0; wellIdx != numWells; ++wellIdx) { - // W{O,G,W}{I,P}R - addWell(std::unique_ptr ( - new WellRate(*this, - eclipseState, - wells[wellIdx], - uses, - ertPhaseIdx, - wellType, - deckUnitType))); - // W{O,G,W}{I,P}T - addWell(std::unique_ptr ( - new WellTotal(*this, - eclipseState, - wells[wellIdx], - uses, - ertPhaseIdx, - wellType, - deckUnitType))); - } - } - } - - // Add BHP monitors - for (int wellIdx = 0; wellIdx != numWells; ++wellIdx) { - // In the call below: uses, phase and the well type arguments - // are not used, except to set up an index that stores the - // well indirectly. For details see the implementation of the - // WellReport constructor, and the method - // WellReport::bhp(). - BlackoilPhases::PhaseIndex ertPhaseIdx = BlackoilPhases::Liquid; - if (!uses.phase_used[BlackoilPhases::Liquid]) { - ertPhaseIdx = BlackoilPhases::Vapour; - } - addWell(std::unique_ptr ( - new WellBhp(*this, - eclipseState, - wells[wellIdx], - uses, - ertPhaseIdx, - WELL_TYPES[0], - deckUnitType))); - } -} -} // end namespace EclipseWriterDetails - - - -/** - * Convert opm-core WellType and InjectorType to eclipse welltype - */ -int EclipseWriter::eclipseWellTypeMask(WellType wellType, WellInjector::TypeEnum injectorType) -{ - int ert_well_type = IWEL_UNDOCUMENTED_ZERO; - - if (PRODUCER == wellType) { - ert_well_type = IWEL_PRODUCER; - } else if (INJECTOR == wellType) { - switch (injectorType) { - case WellInjector::WATER: - ert_well_type = IWEL_WATER_INJECTOR; - break; - case WellInjector::GAS: - ert_well_type = IWEL_GAS_INJECTOR; - break; - case WellInjector::OIL : - ert_well_type = IWEL_OIL_INJECTOR; - break; - default: - ert_well_type = IWEL_UNDOCUMENTED_ZERO; - } - } - - return ert_well_type; -} - - -/** - * Convert opm-core WellStatus to eclipse format: > 0 open, <= 0 shut - */ -int EclipseWriter::eclipseWellStatusMask(WellCommon::StatusEnum wellStatus) -{ - int well_status = 0; - - if (wellStatus == WellCommon::OPEN) { - well_status = 1; - } - return well_status; -} - - - -/** - * Convert opm-core UnitType to eclipse format: ert_ecl_unit_enum - */ -ert_ecl_unit_enum -EclipseWriter::convertUnitTypeErtEclUnitEnum(UnitSystem::UnitType unit) -{ - switch (unit) { - case UnitSystem::UNIT_TYPE_METRIC: - return ERT_ECL_METRIC_UNITS; - - case UnitSystem::UNIT_TYPE_FIELD: - return ERT_ECL_FIELD_UNITS; - - case UnitSystem::UNIT_TYPE_LAB: - return ERT_ECL_LAB_UNITS; - } - - throw std::invalid_argument("unhandled enum value"); -} - - -void EclipseWriter::writeInit(const SimulatorTimerInterface &timer) -{ - // if we don't want to write anything, this method becomes a - // no-op... - if (!enableOutput_) { - return; - } - - writeStepIdx_ = 0; - reportStepIdx_ = -1; - - EclipseWriterDetails::Init fortio(outputDir_, baseName_, /*stepIdx=*/0, eclipseState_->getIOConfigConst()); - fortio.writeHeader(numCells_, - compressedToCartesianCellIdx_, - timer, - eclipseState_, - phaseUsage_); - - IOConfigConstPtr ioConfig = eclipseState_->getIOConfigConst(); - - if (ioConfig->getWriteINITFile()) { - if (eclipseState_->hasDeckDoubleGridProperty("PERMX")) { - auto data = eclipseState_->getDoubleGridProperty("PERMX")->getData(); - EclipseWriterDetails::convertFromSiTo(data, Opm::prefix::milli * Opm::unit::darcy); - EclipseWriterDetails::restrictAndReorderToActiveCells(data, gridToEclipseIdx_.size(), gridToEclipseIdx_.data()); - fortio.writeKeyword("PERMX", data); - } - if (eclipseState_->hasDeckDoubleGridProperty("PERMY")) { - auto data = eclipseState_->getDoubleGridProperty("PERMY")->getData(); - EclipseWriterDetails::convertFromSiTo(data, Opm::prefix::milli * Opm::unit::darcy); - EclipseWriterDetails::restrictAndReorderToActiveCells(data, gridToEclipseIdx_.size(), gridToEclipseIdx_.data()); - fortio.writeKeyword("PERMY", data); - } - if (eclipseState_->hasDeckDoubleGridProperty("PERMZ")) { - auto data = eclipseState_->getDoubleGridProperty("PERMZ")->getData(); - EclipseWriterDetails::convertFromSiTo(data, Opm::prefix::milli * Opm::unit::darcy); - EclipseWriterDetails::restrictAndReorderToActiveCells(data, gridToEclipseIdx_.size(), gridToEclipseIdx_.data()); - fortio.writeKeyword("PERMZ", data); - } - } - - /* Create summary object (could not do it at construction time, - since it requires knowledge of the start time). */ - { - auto eclGrid = eclipseState_->getEclipseGrid(); - auto deckUnitType = eclipseState_->getDeckUnitSystem().getType(); - bool time_in_days = true; - - if (deckUnitType == UnitSystem::UNIT_TYPE_LAB) - time_in_days = false; - - summary_.reset(new EclipseWriterDetails::Summary(outputDir_, - baseName_, - timer, - time_in_days, - eclGrid->getNX(), - eclGrid->getNY(), - eclGrid->getNZ())); - summary_->addAllWells(eclipseState_, phaseUsage_); - } -} - -// implementation of the writeTimeStep method -void EclipseWriter::writeTimeStep(const SimulatorTimerInterface& timer, - const SimulationDataContainer& reservoirState, - const WellState& wellState, - bool isSubstep) -{ - - // if we don't want to write anything, this method becomes a - // no-op... - if (!enableOutput_) { - return; - } - - - std::vector pressure = reservoirState.pressure(); - EclipseWriterDetails::convertFromSiTo(pressure, deckToSiPressure_); - EclipseWriterDetails::restrictAndReorderToActiveCells(pressure, gridToEclipseIdx_.size(), gridToEclipseIdx_.data()); - - std::vector saturation_water; - std::vector saturation_gas; - - if (phaseUsage_.phase_used[BlackoilPhases::Aqua]) { - saturation_water = reservoirState.saturation(); - EclipseWriterDetails::extractFromStripedData(saturation_water, - /*offset=*/phaseUsage_.phase_pos[BlackoilPhases::Aqua], - /*stride=*/phaseUsage_.num_phases); - EclipseWriterDetails::restrictAndReorderToActiveCells(saturation_water, gridToEclipseIdx_.size(), gridToEclipseIdx_.data()); - } - - - if (phaseUsage_.phase_used[BlackoilPhases::Vapour]) { - saturation_gas = reservoirState.saturation(); - EclipseWriterDetails::extractFromStripedData(saturation_gas, - /*offset=*/phaseUsage_.phase_pos[BlackoilPhases::Vapour], - /*stride=*/phaseUsage_.num_phases); - EclipseWriterDetails::restrictAndReorderToActiveCells(saturation_gas, gridToEclipseIdx_.size(), gridToEclipseIdx_.data()); - } - - - - IOConfigConstPtr ioConfig = eclipseState_->getIOConfigConst(); - - - // Write restart file - if(!isSubstep && ioConfig->getWriteRestartFile(timer.reportStepNum())) - { - const size_t ncwmax = eclipseState_->getSchedule()->getMaxNumCompletionsForWells(timer.reportStepNum()); - const size_t numWells = eclipseState_->getSchedule()->numWells(timer.reportStepNum()); - std::vector wells_ptr = eclipseState_->getSchedule()->getWells(timer.reportStepNum()); - - std::vector zwell_data( numWells * Opm::EclipseWriterDetails::Restart::NZWELZ , ""); - std::vector iwell_data( numWells * Opm::EclipseWriterDetails::Restart::NIWELZ , 0 ); - std::vector icon_data( numWells * ncwmax * Opm::EclipseWriterDetails::Restart::NICONZ , 0 ); - - EclipseWriterDetails::Restart restartHandle(outputDir_, baseName_, timer.reportStepNum(), ioConfig); - - const size_t sz = wellState.bhp().size() + wellState.perfPress().size() + wellState.perfRates().size() + wellState.temperature().size() + wellState.wellRates().size(); - std::vector xwell_data( sz , 0 ); - - restartHandle.addRestartFileXwelData(wellState, xwell_data); - - for (size_t iwell = 0; iwell < wells_ptr.size(); ++iwell) { - WellConstPtr well = wells_ptr[iwell]; - { - size_t wellIwelOffset = Opm::EclipseWriterDetails::Restart::NIWELZ * iwell; - restartHandle.addRestartFileIwelData(iwell_data, timer.reportStepNum(), well , wellIwelOffset); - } - { - size_t wellIconOffset = ncwmax * Opm::EclipseWriterDetails::Restart::NICONZ * iwell; - restartHandle.addRestartFileIconData(icon_data, well->getCompletions( timer.reportStepNum() ), wellIconOffset); - } - zwell_data[ iwell * Opm::EclipseWriterDetails::Restart::NZWELZ ] = well->name().c_str(); - } - - - { - ecl_rsthead_type rsthead_data = {}; - rsthead_data.sim_time = timer.currentPosixTime(); - rsthead_data.nactive = numCells_; - rsthead_data.nx = cartesianSize_[0]; - rsthead_data.ny = cartesianSize_[1]; - rsthead_data.nz = cartesianSize_[2]; - rsthead_data.nwells = numWells; - rsthead_data.niwelz = EclipseWriterDetails::Restart::NIWELZ; - rsthead_data.nzwelz = EclipseWriterDetails::Restart::NZWELZ; - rsthead_data.niconz = EclipseWriterDetails::Restart::NICONZ; - rsthead_data.ncwmax = ncwmax; - rsthead_data.phase_sum = Opm::EclipseWriterDetails::ertPhaseMask(phaseUsage_); - rsthead_data.sim_days = Opm::unit::convert::to(timer.simulationTimeElapsed(), Opm::unit::day); //data for doubhead - - restartHandle.writeHeader(timer, - timer.reportStepNum(), - &rsthead_data); - } - - - restartHandle.add_kw(EclipseWriterDetails::Keyword(IWEL_KW, iwell_data)); - restartHandle.add_kw(EclipseWriterDetails::Keyword(ZWEL_KW, zwell_data)); - restartHandle.add_kw(EclipseWriterDetails::Keyword(OPM_XWEL, xwell_data)); - restartHandle.add_kw(EclipseWriterDetails::Keyword(ICON_KW, icon_data)); - - - EclipseWriterDetails::Solution sol(restartHandle); - sol.add(EclipseWriterDetails::Keyword("PRESSURE", pressure)); - - - // write the cell temperature - std::vector temperature = reservoirState.temperature(); - EclipseWriterDetails::convertFromSiTo(temperature, deckToSiTemperatureFactor_, deckToSiTemperatureOffset_); - EclipseWriterDetails::restrictAndReorderToActiveCells(temperature, gridToEclipseIdx_.size(), gridToEclipseIdx_.data()); - sol.add(EclipseWriterDetails::Keyword("TEMP", temperature)); - - - if (phaseUsage_.phase_used[BlackoilPhases::Aqua]) { - sol.add(EclipseWriterDetails::Keyword(EclipseWriterDetails::saturationKeywordNames[BlackoilPhases::PhaseIndex::Aqua], saturation_water)); - } - - - if (phaseUsage_.phase_used[BlackoilPhases::Vapour]) { - sol.add(EclipseWriterDetails::Keyword(EclipseWriterDetails::saturationKeywordNames[BlackoilPhases::PhaseIndex::Vapour], saturation_gas)); - } - - - // Write RS - Dissolved GOR - if (reservoirState.hasCellData( BlackoilState::GASOILRATIO )) { - const std::vector& rs = reservoirState.getCellData( BlackoilState::GASOILRATIO ); - sol.add(EclipseWriterDetails::Keyword("RS", rs)); - } - - // Write RV - Volatilized oil/gas ratio - if (reservoirState.hasCellData( BlackoilState::RV )) { - const std::vector& rv = reservoirState.getCellData( BlackoilState::RV ); - sol.add(EclipseWriterDetails::Keyword("RV", rv)); - } - } - - - //Write RFT data for current timestep to RFT file - std::shared_ptr eclipseWriteRFTHandler = std::make_shared( - compressedToCartesianCellIdx_, - numCells_, - eclipseState_->getEclipseGrid()->getCartesianSize()); - - // Write RFT file. - { - char * rft_filename = ecl_util_alloc_filename(outputDir_.c_str(), - baseName_.c_str(), - ECL_RFT_FILE, - ioConfig->getFMTOUT(), - 0); - auto unit_type = eclipseState_->getDeckUnitSystem().getType(); - ert_ecl_unit_enum ecl_unit = convertUnitTypeErtEclUnitEnum(unit_type); - std::vector wells = eclipseState_->getSchedule()->getWells(timer.reportStepNum()); - eclipseWriteRFTHandler->writeTimeStep(rft_filename, - ecl_unit, - timer, - wells, - eclipseState_->getEclipseGrid(), - pressure, - saturation_water, - saturation_gas); - free( rft_filename ); - } - - /* Summary variables (well reporting) */ - // TODO: instead of writing the header (smspec) every time, it should - // only be written when there is a change in the well configuration - // (first timestep, in practice), and reused later. but how to do this - // without keeping the complete summary in memory (which will then - // accumulate all the timesteps)? - // - // Note: The answer to the question above is still not settled, but now we do keep - // the complete summary in memory, as a member variable in the EclipseWriter class, - // instead of creating a temporary EclipseWriterDetails::Summary in this function - // every time it is called. This has been changed so that the final summary file - // will contain data from the whole simulation, instead of just the last step. - summary_->writeTimeStep(writeStepIdx_, timer, wellState); - - ++writeStepIdx_; - // store current report index - reportStepIdx_ = timer.reportStepNum(); -} - - -EclipseWriter::EclipseWriter(const parameter::ParameterGroup& params, - Opm::EclipseStateConstPtr eclipseState, - const Opm::PhaseUsage &phaseUsage, - int numCells, - const int* compressedToCartesianCellIdx) - : eclipseState_(eclipseState) - , numCells_(numCells) - , compressedToCartesianCellIdx_(compressedToCartesianCellIdx) - , gridToEclipseIdx_(numCells, int(-1) ) - , phaseUsage_(phaseUsage) -{ - const auto eclGrid = eclipseState->getEclipseGrid(); - cartesianSize_[0] = eclGrid->getNX(); - cartesianSize_[1] = eclGrid->getNY(); - cartesianSize_[2] = eclGrid->getNZ(); - - if( compressedToCartesianCellIdx ) { - // if compressedToCartesianCellIdx available then - // compute mapping to eclipse order - std::map< int , int > indexMap; - for (int cellIdx = 0; cellIdx < numCells; ++cellIdx) { - int cartesianCellIdx = compressedToCartesianCellIdx[cellIdx]; - indexMap[ cartesianCellIdx ] = cellIdx; - } - - int idx = 0; - for( auto it = indexMap.begin(), end = indexMap.end(); it != end; ++it ) { - gridToEclipseIdx_[ idx++ ] = (*it).second; - } - } - else { - // if not compressedToCartesianCellIdx was given use identity - for (int cellIdx = 0; cellIdx < numCells; ++cellIdx) { - gridToEclipseIdx_[ cellIdx ] = cellIdx; - } - } - - - // factor from the pressure values given in the deck to Pascals - deckToSiPressure_ = - eclipseState->getDeckUnitSystem().parse("Pressure")->getSIScaling(); - - // factor and offset from the temperature values given in the deck to Kelvin - deckToSiTemperatureFactor_ = - eclipseState->getDeckUnitSystem().parse("Temperature")->getSIScaling(); - deckToSiTemperatureOffset_ = - eclipseState->getDeckUnitSystem().parse("Temperature")->getSIOffset(); - - init(params); -} - -void EclipseWriter::init(const parameter::ParameterGroup& params) -{ - // get the base name from the name of the deck - using boost::filesystem::path; - path deckPath(params.get ("deck_filename")); - if (boost::to_upper_copy(path(deckPath.extension()).string()) == ".DATA") { - baseName_ = path(deckPath.stem()).string(); - } - else { - baseName_ = path(deckPath.filename()).string(); - } - - // make uppercase of everything (or otherwise we'll get uppercase - // of some of the files (.SMSPEC, .UNSMRY) and not others - baseName_ = boost::to_upper_copy(baseName_); - - // retrieve the value of the "output" parameter - enableOutput_ = params.getDefault("output", /*defaultValue=*/true); - - // store in current directory if not explicitly set - outputDir_ = params.getDefault("output_dir", "."); - - // set the index of the first time step written to 0... - writeStepIdx_ = 0; - reportStepIdx_ = -1; - - if (enableOutput_) { - // make sure that the output directory exists, if not try to create it - if (!boost::filesystem::exists(outputDir_)) { - std::cout << "Trying to create directory \"" << outputDir_ << "\" for the simulation output\n"; - boost::filesystem::create_directories(outputDir_); - } - - if (!boost::filesystem::is_directory(outputDir_)) { - OPM_THROW(std::runtime_error, - "The path specified as output directory '" << outputDir_ - << "' is not a directory"); - } - } -} - -// default destructor is OK, just need to be defined -EclipseWriter::~EclipseWriter() -{ } - -} // namespace Opm diff --git a/opm/core/io/eclipse/EclipseWriter.hpp b/opm/core/io/eclipse/EclipseWriter.hpp deleted file mode 100644 index 8841e5d6..00000000 --- a/opm/core/io/eclipse/EclipseWriter.hpp +++ /dev/null @@ -1,139 +0,0 @@ -/* - Copyright (c) 2013 Andreas Lauser - Copyright (c) 2013 Uni Research AS - Copyright (c) 2014 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 . -*/ - -#ifndef OPM_ECLIPSE_WRITER_HPP -#define OPM_ECLIPSE_WRITER_HPP - -#include -#include -#include // WellType -#include - -#include -#include -#include - -#include - -#include -#include -#include -#include - -struct UnstructuredGrid; - -namespace Opm { - -// forward declarations -namespace EclipseWriterDetails { -class Summary; -} - -class SimulationDataContainer; -class WellState; - -namespace parameter { class ParameterGroup; } - -/*! - * \brief A class to write the reservoir state and the well state of a - * blackoil simulation to disk using the Eclipse binary format. - * - * This class only writes files if the 'write_output' parameter is set - * to 1. It needs the ERT libraries to write to disk, so if the - * 'write_output' parameter is set but ERT is not available, all - * methods throw a std::runtime_error. - */ -class EclipseWriter : public OutputWriter -{ -public: - /*! - * \brief Sets the common attributes required to write eclipse - * binary files using ERT. - */ - EclipseWriter(const parameter::ParameterGroup& params, - Opm::EclipseStateConstPtr eclipseState, - const Opm::PhaseUsage &phaseUsage, - int numCells, - const int* compressedToCartesianCellIdx); - - /** - * We need a destructor in the compilation unit to avoid the - * EclipseSummary being a complete type here. - */ - virtual ~EclipseWriter (); - - /** - * Write the static eclipse data (grid, PVT curves, etc) to disk. - */ - virtual void writeInit(const SimulatorTimerInterface &timer); - - /*! - * \brief Write a reservoir state and summary information to disk. - * - * - * The reservoir state can be inspected with visualization tools like - * ResInsight. - * - * The summary information can then be visualized using tools from - * ERT or ECLIPSE. Note that calling this method is only - * meaningful after the first time step has been completed. - * - * \param[in] timer The timer providing time step and time information - * \param[in] reservoirState The thermodynamic state of the reservoir - * \param[in] wellState The production/injection data for all wells - */ - virtual void writeTimeStep(const SimulatorTimerInterface& timer, - const SimulationDataContainer& reservoirState, - const WellState& wellState, - bool isSubstep); - - - static int eclipseWellTypeMask(WellType wellType, WellInjector::TypeEnum injectorType); - static int eclipseWellStatusMask(WellCommon::StatusEnum wellStatus); - static ert_ecl_unit_enum convertUnitTypeErtEclUnitEnum(UnitSystem::UnitType unit); - -private: - Opm::EclipseStateConstPtr eclipseState_; - int numCells_; - std::array cartesianSize_; - const int* compressedToCartesianCellIdx_; - std::vector< int > gridToEclipseIdx_; - double deckToSiPressure_; - double deckToSiTemperatureFactor_; - double deckToSiTemperatureOffset_; - bool enableOutput_; - int writeStepIdx_; - int reportStepIdx_; - std::string outputDir_; - std::string baseName_; - PhaseUsage phaseUsage_; // active phases in the input deck - std::shared_ptr summary_; - - void init(const parameter::ParameterGroup& params); -}; - -typedef std::shared_ptr EclipseWriterPtr; -typedef std::shared_ptr EclipseWriterConstPtr; - -} // namespace Opm - - -#endif // OPM_ECLIPSE_WRITER_HPP diff --git a/opm/core/io/eclipse/writeECLData.cpp b/opm/core/io/eclipse/writeECLData.cpp deleted file mode 100644 index 495090f1..00000000 --- a/opm/core/io/eclipse/writeECLData.cpp +++ /dev/null @@ -1,193 +0,0 @@ -/* - Copyright 2012 SINTEF ICT, Applied Mathematics. - Copyright 2012 Statoil ASA. - - 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 . -*/ - -#if HAVE_CONFIG_H -#include "config.h" -#endif - -#include -#include -#include -#include - -#include - -#ifdef HAVE_ERT // This one goes almost to the bottom of the file - -#include -#include -#include - - -namespace Opm -{ - - static ecl_kw_type * ecl_kw_wrapper( const UnstructuredGrid& grid, - const std::string& kw_name , - const std::vector * data , - int offset , - int stride ) { - - if (stride <= 0) - OPM_THROW(std::runtime_error, "Vector strides must be positive. Got stride = " << stride); - if ((stride * std::vector::size_type(grid.number_of_cells)) != data->size()) - OPM_THROW(std::runtime_error, "Internal mismatch grid.number_of_cells: " << grid.number_of_cells << " data size: " << data->size() / stride); - { - ecl_kw_type * ecl_kw = ecl_kw_alloc( kw_name.c_str() , grid.number_of_cells , ECL_FLOAT_TYPE ); - for (int i=0; i < grid.number_of_cells; i++) - ecl_kw_iset_float( ecl_kw , i , (*data)[i*stride + offset]); - return ecl_kw; - } - } - - - /* - This function will write the data solution data in the DataMap - @data as an ECLIPSE restart file, in addition to the solution - fields the ECLIPSE restart file will have a minimum (hopefully - sufficient) amount of header information. - - The ECLIPSE restart files come in two varietes; unified restart - files which have all the report steps lumped together in one large - chunk and non-unified restart files which are one file for each - report step. In addition the files can be either formatted - (i.e. ASCII) or unformatted (i.e. binary). - - The writeECLData() function has two hardcoded settings: - 'file_type' and 'fmt_file' which regulate which type of files the - should be created. The extension of the files follow a convention: - - Unified, formatted : .FUNRST - Unified, unformatted : .UNRST - Multiple, formatted : .Fnnnn - Multiple, unformatted : .Xnnnn - - For the multiple files the 'nnnn' part is the report number, - formatted with '%04d' format specifier. The writeECLData() - function will use the ecl_util_alloc_filename() function to create - an ECLIPSE filename according to this conventions. - */ - - void writeECLData(const UnstructuredGrid& grid, - const DataMap& data, - const int current_step, - const double current_time, - const boost::posix_time::ptime& current_date_time, - const std::string& output_dir, - const std::string& base_name) { - - ecl_file_enum file_type = ECL_UNIFIED_RESTART_FILE; // Alternatively ECL_RESTART_FILE for multiple restart files. - bool fmt_file = false; - - char * filename = ecl_util_alloc_filename(output_dir.c_str() , base_name.c_str() , file_type , fmt_file , current_step ); - int phases = ECL_OIL_PHASE + ECL_WATER_PHASE; - time_t date = 0; - int nx = grid.cartdims[0]; - int ny = grid.cartdims[1]; - int nz = grid.cartdims[2]; - int nactive = grid.number_of_cells; - ecl_rst_file_type * rst_file; - - { - using namespace boost::posix_time; - ptime t0( boost::gregorian::date(1970 , 1 ,1) ); - time_duration::sec_type seconds = (current_date_time - t0).total_seconds(); - - date = time_t( seconds ); - } - - if (current_step > 0 && file_type == ECL_UNIFIED_RESTART_FILE) - rst_file = ecl_rst_file_open_append( filename ); - else - rst_file = ecl_rst_file_open_write( filename ); - - { - ecl_rsthead_type rsthead_data = {}; - - const int num_wells = 0; - const int niwelz = 0; - const int nzwelz = 0; - const int niconz = 0; - const int ncwmax = 0; - - rsthead_data.nx = nx; - rsthead_data.ny = ny; - rsthead_data.nz = nz; - rsthead_data.nwells = num_wells; - rsthead_data.niwelz = niwelz; - rsthead_data.nzwelz = nzwelz; - rsthead_data.niconz = niconz; - rsthead_data.ncwmax = ncwmax; - rsthead_data.nactive = nactive; - rsthead_data.phase_sum = phases; - rsthead_data.sim_time = date; - - rsthead_data.sim_days = Opm::unit::convert::to(current_time, Opm::unit::day); //Data for doubhead - - ecl_rst_file_fwrite_header( rst_file , current_step , &rsthead_data); - } - - ecl_rst_file_start_solution( rst_file ); - - { - DataMap::const_iterator i = data.find("pressure"); - if (i != data.end()) { - ecl_kw_type * pressure_kw = ecl_kw_wrapper( grid , "PRESSURE" , i->second , 0 , 1); - ecl_rst_file_add_kw( rst_file , pressure_kw ); - ecl_kw_free( pressure_kw ); - } - } - - { - DataMap::const_iterator i = data.find("saturation"); - if (i != data.end()) { - if (int(i->second->size()) != 2 * grid.number_of_cells) { - OPM_THROW(std::runtime_error, "writeECLData() requires saturation field to have two phases."); - } - ecl_kw_type * swat_kw = ecl_kw_wrapper( grid , "SWAT" , i->second , 0 , 2); - ecl_rst_file_add_kw( rst_file , swat_kw ); - ecl_kw_free( swat_kw ); - } - } - - ecl_rst_file_end_solution( rst_file ); - ecl_rst_file_close( rst_file ); - free(filename); - } -} - -#else // that is, we have not defined HAVE_ERT - -namespace Opm -{ - - void writeECLData(const UnstructuredGrid&, - const DataMap&, - const int, - const double, - const boost::posix_time::ptime&, - const std::string&, - const std::string&) - { - OPM_THROW(std::runtime_error, "Cannot call writeECLData() without ERT library support. Reconfigure opm-core with ERT support and recompile."); - } -} - -#endif diff --git a/opm/core/io/eclipse/writeECLData.hpp b/opm/core/io/eclipse/writeECLData.hpp deleted file mode 100644 index e962573a..00000000 --- a/opm/core/io/eclipse/writeECLData.hpp +++ /dev/null @@ -1,46 +0,0 @@ -/* - Copyright 2012 SINTEF ICT, Applied Mathematics. - Copyright 2012 Statoil ASA. - - 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 . -*/ - -#ifndef OPM_WRITEECLDATA_HEADER_INCLUDED -#define OPM_WRITEECLDATA_HEADER_INCLUDED - - -#include -#include - -#include - -struct UnstructuredGrid; - -namespace Opm -{ - - // ECLIPSE output for general grids. - void writeECLData(const UnstructuredGrid& grid, - const DataMap& data, - const int current_step, - const double current_time, - const boost::posix_time::ptime& current_date_time, - const std::string& output_dir, - const std::string& base_name); - -} - -#endif diff --git a/opm/core/io/vag/vag.cpp b/opm/core/io/vag/vag.cpp deleted file mode 100644 index fa3f4b48..00000000 --- a/opm/core/io/vag/vag.cpp +++ /dev/null @@ -1,408 +0,0 @@ -/*=========================================================================== -// -// File: vag.cpp -// -// Created: 2012-06-08 15:45:53+0200 -// -// Authors: Knut-Andreas Lie -// Halvor M. Nilsen -// Atgeirr F. Rasmussen -// Xavier Raynaud -// Bård Skaflestad -// -//==========================================================================*/ - - -/* - Copyright 2012 SINTEF ICT, Applied Mathematics. - Copyright 2012 Statoil ASA. - - 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 -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -namespace Opm -{ - void readPosStruct(std::istream& is,int n,PosStruct& pos_struct){ - using namespace std; - //PosStruct pos_struct; - pos_struct.pos.resize(n+1); - pos_struct.pos[0]=0; - for(int i=0;i< n;++i){ - int number; - is >> number ; - //cout <> value; - // cout << value << " "; - pos_struct.value.push_back(value); - } - //cout << endl; - } - if(int(pos_struct.value.size()) != pos_struct.pos[n]){ - cerr << "Failed to read pos structure" << endl; - cerr << "pos_struct.value.size()" << pos_struct.value.size() << endl; - cerr << "pos_struct.pos[n+1]" << pos_struct.pos[n] << endl; - } - } - void writePosStruct(std::ostream& os,PosStruct& pos_struct){ - using namespace std; - //PosStruct pos_struct; - if(pos_struct.pos.size()==0){ - return; - } - int n=pos_struct.pos.size()-1; - pos_struct.pos.resize(n+1); - pos_struct.pos[0]=0; - for(int i=0;i< n;++i){ - int number=pos_struct.pos[i+1]-pos_struct.pos[i]; - os << number << " "; - for(int j=0;j< number;++j){ - os << pos_struct.value[pos_struct.pos[i]+j] << " "; - } - os << endl; - } - } - void readVagGrid(std::istream& is,Opm::VAG& vag_grid){ - using namespace std; - using namespace Opm; - while (!is.eof()) { - string keyword; - is >> keyword; - //cout << keyword<< endl; - if(keyword == "Number"){ - string stmp; - is >> stmp; - if(stmp == "of"){ - string entity; - is >> entity; - getline(is,stmp); - int number; - is >> number; - if(entity=="vertices"){ - vag_grid.number_of_vertices=number; - }else if((entity=="volumes") || (entity=="control")){ - vag_grid.number_of_volumes=number; - }else if(entity=="faces"){ - vag_grid.number_of_faces=number; - }else if(entity=="edges"){ - vag_grid.number_of_edges=number; - } - cout << "Found Number of: " << entity <<" " << number << endl; - } else { - cerr << "Wrong format: Not of after Number" << endl; - return; - } - }else{ - // read geometry defined by vertices - if(keyword=="Vertices"){ - int number; - is >> number; - vag_grid.vertices.resize(3*number);// assume 3d data - readVector(is,vag_grid.vertices); - } - // here starts the reding of all pos structures - else if(keyword=="Volumes->Faces" || keyword=="Volumes->faces"){ - //vag_grid.volumes_to_faces= - int number; - is >> number; - readPosStruct(is,number,vag_grid.volumes_to_faces); - cout << "Volumes->Faces: Number of " << number << endl; - }else if(keyword=="Faces->edges" || keyword=="Faces->Edges" || keyword=="Faces->Edgess"){ - int number; - is >> number; - //vag_grid.volumes_to_faces= - readPosStruct(is,number,vag_grid.faces_to_edges); - cout << "Faces->edges: Number of " << number << endl; - }else if(keyword=="Faces->Vertices" || keyword=="Faces->vertices"){ - int number; - is >> number; - //vag_grid.volumes_to_faces= - readPosStruct(is,number,vag_grid.faces_to_vertices); - cout << "Faces->Vertices: Number of " << number << endl; - }else if(keyword=="Volumes->Vertices" || keyword=="Volumes->Verticess"){ - int number; - is >> number; - //vag_grid.volumes_to_faces= - readPosStruct(is,number,vag_grid.volumes_to_vertices); - cout << "Volumes->Vertices: Number of " << number << endl; - } - - // read simple mappings - else if(keyword=="Edge" || keyword=="Edges"){ - int number; - is >> number; - vag_grid.edges.resize(2*number); - readVector(is,vag_grid.edges); - cout << "Edges: Number of " << number << endl; - }else if(keyword=="Faces->Volumes" || keyword=="Faces->Control"){ - int number; - if(keyword=="Faces->Control"){ - string vol; - is >> vol; - } - is >> number; - vag_grid.faces_to_volumes.resize(2*number); - readVector(is,vag_grid.faces_to_volumes); - cout << "Faces->Volumes: Number of " << number << endl; - } - // read material - else if(keyword=="Material"){ - string snum; - is >> snum; - int number; - is >> number; - cout << "Material number " << number << endl; - // we read all the rest into doubles - while(!is.eof()){ - double value; - is >> value; - //cout << value << endl; - vag_grid.material.push_back(value); - } - }else{ - //cout << "keyword; - } - //cout << "Found" << keyword << "Number of " << number << endl; - } - } - } - void vagToUnstructuredGrid(Opm::VAG& vag_grid,UnstructuredGrid& grid){ - using namespace std; - using namespace Opm; - cout << "Converting grid" << endl; - cout << "Warning:: orignial grid may not be edge confomal" << endl; - cout << " inverse mappings from edges will be wrong" << endl; - grid.dimensions=3; - grid.number_of_cells=vag_grid.number_of_volumes; - grid.number_of_faces=vag_grid.number_of_faces; - grid.number_of_nodes=vag_grid.number_of_vertices; - - // fill face_nodes - for(int i=0;i< int(vag_grid.faces_to_vertices.pos.size());++i){ - grid.face_nodepos[i] = vag_grid.faces_to_vertices.pos[i]; - } - for(int i=0;i< int(vag_grid.faces_to_vertices.value.size());++i){ - grid.face_nodes[i] = vag_grid.faces_to_vertices.value[i]-1; - } - // fill cell_face - for(int i=0;i< int(vag_grid.volumes_to_faces.pos.size());++i){ - grid.cell_facepos[i] = vag_grid.volumes_to_faces.pos[i]; - } - for(int i=0;i< int(vag_grid.volumes_to_faces.value.size());++i){ - grid.cell_faces[i] = vag_grid.volumes_to_faces.value[i]-1; - } - // fill face_cells - for(int i=0;i< int(vag_grid.faces_to_volumes.size());++i){ - grid.face_cells[i] = vag_grid.faces_to_volumes[i]-1; - } - - // fill node_cordinates. This is the only geometry given in the vag - for(int i=0;i< int(vag_grid.vertices.size());++i){ - grid.node_coordinates[i] = vag_grid.vertices[i]; - } - // informations in edges, faces_to_eges, faces_to_vertices, volume_to_vertices and materials - // is not used - cout << "Computing geometry" << endl; - compute_geometry(&grid); - - } - - void unstructuredGridToVag(UnstructuredGrid& grid,Opm::VAG& vag_grid){ - using namespace std; - using namespace Opm; - cout << "Converting grid" << endl; - // grid.dimensions=3; - vag_grid.number_of_volumes=grid.number_of_cells; - vag_grid.number_of_faces=grid.number_of_faces; - vag_grid.number_of_vertices=grid.number_of_nodes; - - // resizing vectors - vag_grid.vertices.resize(grid.number_of_nodes*3); - vag_grid.faces_to_vertices.pos.resize(grid.number_of_faces+1); - vag_grid.faces_to_vertices.value.resize(grid.face_nodepos[grid.number_of_faces]); - vag_grid.faces_to_volumes.resize(2*grid.number_of_faces); - vag_grid.volumes_to_faces.pos.resize(grid.number_of_cells+1); - vag_grid.volumes_to_faces.value.resize(grid.cell_facepos[grid.number_of_cells]);//not known - - - - - // fill face_nodes - for(int i=0;i< int(vag_grid.faces_to_vertices.pos.size());++i){ - vag_grid.faces_to_vertices.pos[i] = grid.face_nodepos[i]; - } - - for(int i=0;i< int(vag_grid.faces_to_vertices.value.size());++i){ - vag_grid.faces_to_vertices.value[i] = grid.face_nodes[i] +1; - } - - // fill cell_face - for(int i=0;i< int(vag_grid.volumes_to_faces.pos.size());++i){ - vag_grid.volumes_to_faces.pos[i] = grid.cell_facepos[i]; - } - for(int i=0;i< int(vag_grid.volumes_to_faces.value.size());++i){ - vag_grid.volumes_to_faces.value[i] = grid.cell_faces[i] +1; - } - // fill face_cells - for(int i=0;i< int(vag_grid.faces_to_volumes.size());++i){ - vag_grid.faces_to_volumes[i] = grid.face_cells[i] +1; - } - - // fill node_cordinates. This is the only geometry given in the vag - for(int i=0;i< int(vag_grid.vertices.size());++i){ - vag_grid.vertices[i] = grid.node_coordinates[i]; - } - - - // The missing field need to be constructed - // gennerate volume to vertice mapping - std::vector< std::set > volumes_to_vertices(grid.number_of_cells); - for(int i=0;i < grid.number_of_cells; ++i){ - int nlf=grid.cell_facepos[i+1]-grid.cell_facepos[i]; - std::set nodes; - for(int j=0; j < nlf; ++j){ - int face = grid.cell_faces[grid.cell_facepos[i]+j]; - int nlv = grid.face_nodepos[face+1]-grid.face_nodepos[face]; - for(int k=0; k< nlv; ++k){ - int node = grid.face_nodes[grid.face_nodepos[face]+k]+1; - nodes.insert(node); - } - } - volumes_to_vertices[i]=nodes; - } - - // fill volume to vertice map - vag_grid.volumes_to_vertices.pos.resize(grid.number_of_cells+1); - vag_grid.volumes_to_vertices.value.resize(0); - vag_grid.volumes_to_vertices.pos[0]=0; - for(int i=0;i < grid.number_of_cells;++i){ - int nv=volumes_to_vertices[i].size(); - vag_grid.volumes_to_vertices.pos[i+1]=vag_grid.volumes_to_vertices.pos[i]+nv; - std::set::iterator it; - for(it=volumes_to_vertices[i].begin();it!=volumes_to_vertices[i].end();++it){ - vag_grid.volumes_to_vertices.value.push_back(*it); - } - } - - std::set< std::set > edges; - std::vector< std::vector< std::set > > faces_spares; - int nfe=0; - faces_spares.resize(grid.number_of_faces); - for(int i=0;i < grid.number_of_faces;++i){ - int ne=grid.face_nodepos[i+1]-grid.face_nodepos[i]; - nfe=nfe+ne; - - for(int j=0; j < ne-1;++j){ - int node1=grid.face_nodes[grid.face_nodepos[i]+j]+1; - int node2=grid.face_nodes[grid.face_nodepos[i]+j+1]+1; - std::set spair; - spair.insert(node1); - spair.insert(node2); - edges.insert(spair); - faces_spares[i].push_back(spair); - } - // add end segment - { - std::set spair; - int node1=grid.face_nodes[grid.face_nodepos[i]+ne-1]+1; - int node2=grid.face_nodes[grid.face_nodepos[i]]+1; - spair.insert(node1); - spair.insert(node2); - edges.insert(spair); - faces_spares[i].push_back(spair); - } - } - - // make edge numbering and fill edges - std::map, int> edge_map; - std::set< std::set >::iterator it; - vag_grid.edges.resize(0); - int k=0; - for(it=edges.begin(); it!=edges.end();++it){ - edge_map.insert(std::pair< std::set , int >(*it,k)); - k=k+1; - std::set::iterator sit; - for(sit=(*it).begin();sit!=(*it).end();++sit){ - vag_grid.edges.push_back(*sit); - } - } - // fill face_to_egdes - vag_grid.number_of_edges=edges.size(); - vag_grid.faces_to_edges.pos.resize(vag_grid.number_of_faces+1); - for(int i=0;i < grid.number_of_faces;++i){ - int ne=grid.face_nodepos[i+1]-grid.face_nodepos[i]; - vag_grid.faces_to_edges.pos[i+1]=vag_grid.faces_to_edges.pos[i]+ne; - for(int j=0;jfaces " << vag_grid.volumes_to_faces.pos.size()-1 << endl; - writePosStruct(os, vag_grid.volumes_to_faces); - os << "Volumes->Vertices " << vag_grid.volumes_to_vertices.pos.size()-1 << endl; - writePosStruct(os, vag_grid.volumes_to_vertices); - os << "Faces->edges " << vag_grid.faces_to_edges.pos.size()-1 << endl; - writePosStruct(os, vag_grid.faces_to_edges); - os << "Faces->vertices " << vag_grid.faces_to_vertices.pos.size()-1 << endl; - writePosStruct(os, vag_grid.faces_to_vertices); - os << "Faces->Control volumes " << floor(vag_grid.faces_to_volumes.size()/2) << endl; - writeVector(os,vag_grid.faces_to_volumes,2); - os << "Edges " << floor(vag_grid.edges.size()/2) << endl; - writeVector(os,vag_grid.edges,2); - /* - assert(vag_grid.material.size()%vag_grid.number_of_volumes==0); - int lines= floor(vag_grid.material.size()/vag_grid.number_of_volumes); - os << "Material number " << 1 << endl; - writeVector(os,vag_grid.material,lines); - */ - - } - - -} diff --git a/opm/core/io/vag/vag.hpp b/opm/core/io/vag/vag.hpp deleted file mode 100644 index ab84af77..00000000 --- a/opm/core/io/vag/vag.hpp +++ /dev/null @@ -1,165 +0,0 @@ -/*=========================================================================== -// -// File: vag.hpp -// -// Created: 2012-06-08 15:46:23+0200 -// -// Authors: Knut-Andreas Lie -// Halvor M. Nilsen -// Atgeirr F. Rasmussen -// Xavier Raynaud -// Bård Skaflestad -// -//==========================================================================*/ - - -/* - Copyright 2012 SINTEF ICT, Applied Mathematics. - Copyright 2012 Statoil ASA. - - 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 . -*/ - -#ifndef OPM_VAG_HPP_HEADER -#define OPM_VAG_HPP_HEADER - - - -#include -#include -#include -#include -#include - -namespace Opm -{ - /** - Struct to hold mapping from the natural numbers less than pos.size()-1 to - a set of integers. value(pos(i):pos(i+1)-1) hold the integers corresponding to i. - pos(end)-1==value.size(); - */ - struct PosStruct{ - std::vector pos; - std::vector value; - }; - /** - Structure to represent the unstructured vag grid format. The format is only for - 3D grids. - */ - struct VAG{ - int number_of_vertices; - int number_of_volumes; - int number_of_faces; - int number_of_edges; - /** Vertices. The coordinates of vertice i is [vetices[3*i:3*i+2]*/ - std::vector vertices; - /** Mapping from volumes to faces */ - PosStruct volumes_to_faces; - /** Mapping from volumes to vertices */ - PosStruct volumes_to_vertices; - /** Mapping from faces to edges */ - PosStruct faces_to_edges; - /** Mapping from faces to vertices */ - PosStruct faces_to_vertices; - /** The edge i is given by the nodes edges[2*i:2*i+1] */ - std::vector edges; - /** The two neigbours of the face i is faces_to_volumes[2*i:2*i+1] */ - std::vector faces_to_volumes; - /** A vector containing information of each volume. The size is n*number_of_volumes. - For each i this is the information: - material[n*i] is the volume number and should be transformed to integer - material[n*i+1] is a tag and should be transformed to integer - material[n*i+2:n*(i+1)-1] represent propertices. - */ - std::vector material; - }; - /** - Function the vag grid format and make a vag_grid struct. This structure - is intended to be converted to a grid. - \param[in] is is is stream of the file. - \param[out] vag_grid is a reference to a vag_grid struct. - */ - void readVagGrid(std::istream& is,Opm::VAG& vag_grid); - /** - Function to write vag format. - \param[out] is is is stream of the file. - \param[in] vag_grid is a reference to a vag_grid struct. - - */ - void writeVagFormat(std::ostream& os,Opm::VAG& vag_grid); - /** - Function to read a vector of some type from a stream. - \param[in] os is is stream of the file. - \param[out] vag_grid is a resized and filled vector containing the quantiy read. - */ - template - void readVector(std::istream& is,std::vector& vec){ - using namespace std; - for(int i=0;i< int(vec.size());++i){ - is >> vec[i]; - } - } - /** - Function to write a vector of some type from a stream. - \param[in] os is is stream of the file. - \param[out] vag_grid is a resized and filled vector containing the quantiy read. - \param[in] n number of doubles on each line. - */ - template - void writeVector(std::ostream& os,std::vector& vec,int n){ - typedef typename std::vector::size_type sz_t; - - const sz_t nn = n; - - for (sz_t i = 0; i < vec.size(); ++i) { - os << vec[i] << (((i % nn) == 0) ? '\n' : ' '); - } - - if ((vec.size() % nn) != 0) { - os << '\n'; - } - } - - /** - Read pos struct type mapping from a stream - \param[in] is is stream - \param[in] n number of lines to read - \param[out] pos_struct reference to PosStruct - */ - void readPosStruct(std::istream& is,int n,PosStruct& pos_struct); - /** - Read pos struct type mapping from a stream - \param[in] os is stream to write to - \param[in] pos_struct to write - */ - void writePosStruct(std::ostream& os,PosStruct& pos_struct); - - /** - Fill a UnstructuredGrid from a vag_grid. - \param[out] vag_grid s is a valid vag_grid struct. - \param[in] grid is a grid with have allocated correct size to each pointer. - */ - void vagToUnstructuredGrid(Opm::VAG& vag_grid,UnstructuredGrid& grid); - - /** - Fill a vag_grid from UnstructuredGrid - \param[out] vag_grid s is a valid vag_grid struct. - \param[in] grid is a grid with have allocated correct size to each pointer. - */ - void unstructuredGridToVag(UnstructuredGrid& grid, Opm::VAG& vag_grid); -} -#endif /* OPM_VAG_HPP_HEADER */ - diff --git a/opm/core/io/vtk/writeVtkData.cpp b/opm/core/io/vtk/writeVtkData.cpp deleted file mode 100644 index 2e5a8fa4..00000000 --- a/opm/core/io/vtk/writeVtkData.cpp +++ /dev/null @@ -1,319 +0,0 @@ -/* - Copyright 2012 SINTEF ICT, Applied Mathematics. - - 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 -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - - - -namespace Opm -{ - - void writeVtkData(const std::array& dims, - const std::array& cell_size, - const DataMap& data, - std::ostream& os) - { - // Dimension is hardcoded in the prototype and the next two lines, - // but the rest is flexible (allows dimension == 2 or 3). - int dimension = 3; - int num_cells = dims[0]*dims[1]*dims[2]; - - assert(dimension == 2 || dimension == 3); - assert(num_cells == dims[0]*dims[1]* (dimension == 2 ? 1 : dims[2])); - - os << "# vtk DataFile Version 2.0\n"; - os << "Structured Grid\n \n"; - os << "ASCII \n"; - os << "DATASET STRUCTURED_POINTS\n"; - - os << "DIMENSIONS " - << dims[0] + 1 << " " - << dims[1] + 1 << " "; - if (dimension == 3) { - os << dims[2] + 1; - } else { - os << 1; - } - os << "\n"; - - os << "ORIGIN " << 0.0 << " " << 0.0 << " " << 0.0 << "\n"; - - os << "SPACING " << cell_size[0] << " " << cell_size[1]; - if (dimension == 3) { - os << " " << cell_size[2]; - } else { - os << " " << 0.0; - } - os << "\n"; - - os << "\nCELL_DATA " << num_cells << '\n'; - for (DataMap::const_iterator dit = data.begin(); dit != data.end(); ++dit) { - std::string name = dit->first; - os << "SCALARS " << name << " float" << '\n'; - os << "LOOKUP_TABLE " << name << "_table " << '\n'; - const std::vector& field = *(dit->second); - // We always print only the first data item for every - // cell, using 'stride'. - // This is a hack to get water saturation nicely. - // \TODO: Extend to properly printing vector data. - const int stride = field.size()/num_cells; - const int num_per_line = 5; - for (int c = 0; c < num_cells; ++c) { - os << field[stride*c] << ' '; - if (c % num_per_line == num_per_line - 1 - || c == num_cells - 1) { - os << '\n'; - } - } - } - } - - typedef std::map PMap; - - - struct Tag - { - Tag(const std::string& tag, const PMap& props, std::ostream& os) - : name_(tag), os_(os) - { - indent(os); - os << "<" << tag; - for (PMap::const_iterator it = props.begin(); it != props.end(); ++it) { - os << " " << it->first << "=\"" << it->second << "\""; - } - os << ">\n"; - ++indent_; - } - Tag(const std::string& tag, std::ostream& os) - : name_(tag), os_(os) - { - indent(os); - os << "<" << tag << ">\n"; - ++indent_; - } - ~Tag() - { - --indent_; - indent(os_); - os_ << "\n"; - } - static void indent(std::ostream& os) - { - for (int i = 0; i < indent_; ++i) { - os << " "; - } - } - private: - static int indent_; - std::string name_; - std::ostream& os_; - }; - - int Tag::indent_ = 0; - - - void writeVtkData(const UnstructuredGrid& grid, - const DataMap& data, - std::ostream& os) - { - if (grid.dimensions != 3) { - OPM_THROW(std::runtime_error, "Vtk output for 3d grids only"); - } - os.precision(12); - os << "\n"; - PMap pm; - pm["type"] = "UnstructuredGrid"; - Tag vtkfiletag("VTKFile", pm, os); - Tag ugtag("UnstructuredGrid", os); - int num_pts = grid.number_of_nodes; - int num_cells = grid.number_of_cells; - pm.clear(); - pm["NumberOfPoints"] = boost::lexical_cast(num_pts); - pm["NumberOfCells"] = boost::lexical_cast(num_cells); - Tag piecetag("Piece", pm, os); - { - Tag pointstag("Points", os); - pm.clear(); - pm["type"] = "Float64"; - pm["Name"] = "Coordinates"; - pm["NumberOfComponents"] = "3"; - pm["format"] = "ascii"; - Tag datag("DataArray", pm, os); - for (int i = 0; i < num_pts; ++i) { - Tag::indent(os); - os << grid.node_coordinates[3*i + 0] << ' ' - << grid.node_coordinates[3*i + 1] << ' ' - << grid.node_coordinates[3*i + 2] << '\n'; - } - } - { - Tag cellstag("Cells", os); - pm.clear(); - pm["type"] = "Int32"; - pm["NumberOfComponents"] = "1"; - pm["format"] = "ascii"; - std::vector cell_numpts; - cell_numpts.reserve(num_cells); - { - pm["Name"] = "connectivity"; - Tag t("DataArray", pm, os); - int hf = 0; - for (int c = 0; c < num_cells; ++c) { - std::set cell_pts; - for (; hf < grid.cell_facepos[c+1]; ++hf) { - int f = grid.cell_faces[hf]; - const int* fnbeg = grid.face_nodes + grid.face_nodepos[f]; - const int* fnend = grid.face_nodes + grid.face_nodepos[f+1]; - cell_pts.insert(fnbeg, fnend); - } - cell_numpts.push_back(cell_pts.size()); - Tag::indent(os); - std::copy(cell_pts.begin(), cell_pts.end(), - std::ostream_iterator(os, " ")); - os << '\n'; - } - } - { - pm["Name"] = "offsets"; - Tag t("DataArray", pm, os); - int offset = 0; - const int num_per_line = 10; - for (int c = 0; c < num_cells; ++c) { - if (c % num_per_line == 0) { - Tag::indent(os); - } - offset += cell_numpts[c]; - os << offset << ' '; - if (c % num_per_line == num_per_line - 1 - || c == num_cells - 1) { - os << '\n'; - } - } - } - std::vector cell_foffsets; - cell_foffsets.reserve(num_cells); - { - pm["Name"] = "faces"; - Tag t("DataArray", pm, os); - const int* fp = grid.cell_facepos; - int offset = 0; - for (int c = 0; c < num_cells; ++c) { - Tag::indent(os); - os << fp[c+1] - fp[c] << '\n'; - ++offset; - for (int hf = fp[c]; hf < fp[c+1]; ++hf) { - int f = grid.cell_faces[hf]; - const int* np = grid.face_nodepos; - int f_num_pts = np[f+1] - np[f]; - Tag::indent(os); - os << f_num_pts << ' '; - ++offset; - std::copy(grid.face_nodes + np[f], - grid.face_nodes + np[f+1], - std::ostream_iterator(os, " ")); - os << '\n'; - offset += f_num_pts; - } - cell_foffsets.push_back(offset); - } - } - { - pm["Name"] = "faceoffsets"; - Tag t("DataArray", pm, os); - const int num_per_line = 10; - for (int c = 0; c < num_cells; ++c) { - if (c % num_per_line == 0) { - Tag::indent(os); - } - os << cell_foffsets[c] << ' '; - if (c % num_per_line == num_per_line - 1 - || c == num_cells - 1) { - os << '\n'; - } - } - } - { - pm["type"] = "UInt8"; - pm["Name"] = "types"; - Tag t("DataArray", pm, os); - const int num_per_line = 10; - for (int c = 0; c < num_cells; ++c) { - if (c % num_per_line == 0) { - Tag::indent(os); - } - os << "42 "; - if (c % num_per_line == num_per_line - 1 - || c == num_cells - 1) { - os << '\n'; - } - } - } - } - { - pm.clear(); - if (data.find("saturation") != data.end()) { - pm["Scalars"] = "saturation"; - } else if (data.find("pressure") != data.end()) { - pm["Scalars"] = "pressure"; - } - Tag celldatatag("CellData", pm, os); - pm.clear(); - pm["NumberOfComponents"] = "1"; - pm["format"] = "ascii"; - pm["type"] = "Float64"; - for (DataMap::const_iterator dit = data.begin(); dit != data.end(); ++dit) { - pm["Name"] = dit->first; - const std::vector& field = *(dit->second); - const int num_comps = field.size()/grid.number_of_cells; - pm["NumberOfComponents"] = boost::lexical_cast(num_comps); - Tag ptag("DataArray", pm, os); - const int num_per_line = num_comps == 1 ? 5 : num_comps; - for (int item = 0; item < num_cells*num_comps; ++item) { - if (item % num_per_line == 0) { - Tag::indent(os); - } - double value = field[item]; - if (std::fabs(value) < std::numeric_limits::min()) { - // Avoiding denormal numbers to work around - // bug in Paraview. - value = 0.0; - } - os << value << ' '; - if (item % num_per_line == num_per_line - 1 - || item == num_cells - 1) { - os << '\n'; - } - } - } - } - } - -} // namespace Opm - diff --git a/opm/core/io/vtk/writeVtkData.hpp b/opm/core/io/vtk/writeVtkData.hpp deleted file mode 100644 index 226750dd..00000000 --- a/opm/core/io/vtk/writeVtkData.hpp +++ /dev/null @@ -1,48 +0,0 @@ -/* - Copyright 2012 SINTEF ICT, Applied Mathematics. - - 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 . -*/ - -#ifndef OPM_WRITEVTKDATA_HEADER_INCLUDED -#define OPM_WRITEVTKDATA_HEADER_INCLUDED - - -#include -#include -#include -#include -#include -#include - -struct UnstructuredGrid; - -namespace Opm -{ - - /// Vtk output for cartesian grids. - void writeVtkData(const std::array& dims, - const std::array& cell_size, - const DataMap& data, - std::ostream& os); - - /// Vtk output for general grids. - void writeVtkData(const UnstructuredGrid& grid, - const DataMap& data, - std::ostream& os); -} // namespace Opm - -#endif // OPM_WRITEVTKDATA_HEADER_INCLUDED From 77a6c4229a48e32adbcb25ae3d29dcc9697d5dd3 Mon Sep 17 00:00:00 2001 From: "Kjell W. Kongsvik" Date: Fri, 18 Mar 2016 13:17:39 +0100 Subject: [PATCH 4/4] Replace comment with "ifdef DISABLE_OUTPUT" --- .../SimulatorCompressibleTwophase.cpp | 148 ++++++++++-------- .../simulator/SimulatorIncompTwophase.cpp | 147 +++++++++-------- opm/core/simulator/SimulatorOutput.cpp | 17 +- opm/core/simulator/SimulatorOutput.hpp | 9 +- tutorials/tutorial1.cpp | 12 +- tutorials/tutorial2.cpp | 20 ++- tutorials/tutorial3.cpp | 16 +- tutorials/tutorial4.cpp | 14 +- 8 files changed, 219 insertions(+), 164 deletions(-) diff --git a/opm/core/simulator/SimulatorCompressibleTwophase.cpp b/opm/core/simulator/SimulatorCompressibleTwophase.cpp index c223b665..5acad2c2 100644 --- a/opm/core/simulator/SimulatorCompressibleTwophase.cpp +++ b/opm/core/simulator/SimulatorCompressibleTwophase.cpp @@ -35,7 +35,9 @@ #include #include // 17.03.2016 Temporarily removed while moving functionality to opm-output -// #include +#ifdef DISABLE_OUTPUT +#include +#endif #include #include @@ -140,71 +142,75 @@ namespace Opm // 17.03.2016 Temporarily removed while moving functionality to opm-output -// static void outputStateVtk(const UnstructuredGrid& grid, -// const Opm::BlackoilState& state, -// const int step, -// const std::string& output_dir) -// { -// // Write data in VTK format. -// std::ostringstream vtkfilename; -// vtkfilename << output_dir << "/vtk_files"; -// boost::filesystem::path fpath(vtkfilename.str()); -// try { -// create_directories(fpath); -// } -// catch (...) { -// OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); -// } -// vtkfilename << "/output-" << std::setw(3) << std::setfill('0') << step << ".vtu"; -// std::ofstream vtkfile(vtkfilename.str().c_str()); -// if (!vtkfile) { -// OPM_THROW(std::runtime_error, "Failed to open " << vtkfilename.str()); -// } -// Opm::DataMap dm; -// dm["saturation"] = &state.saturation(); -// dm["pressure"] = &state.pressure(); -// std::vector cell_velocity; -// Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); -// dm["velocity"] = &cell_velocity; -// Opm::writeVtkData(grid, dm, vtkfile); -// } +#ifdef DISABLE_OUTPUT + static void outputStateVtk(const UnstructuredGrid& grid, + const Opm::BlackoilState& state, + const int step, + const std::string& output_dir) + { + // Write data in VTK format. + std::ostringstream vtkfilename; + vtkfilename << output_dir << "/vtk_files"; + boost::filesystem::path fpath(vtkfilename.str()); + try { + create_directories(fpath); + } + catch (...) { + OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); + } + vtkfilename << "/output-" << std::setw(3) << std::setfill('0') << step << ".vtu"; + std::ofstream vtkfile(vtkfilename.str().c_str()); + if (!vtkfile) { + OPM_THROW(std::runtime_error, "Failed to open " << vtkfilename.str()); + } + Opm::DataMap dm; + dm["saturation"] = &state.saturation(); + dm["pressure"] = &state.pressure(); + std::vector cell_velocity; + Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); + dm["velocity"] = &cell_velocity; + Opm::writeVtkData(grid, dm, vtkfile); + } +#endif // 17.03.2016 Temporarily removed while moving functionality to opm-output -// static void outputStateMatlab(const UnstructuredGrid& grid, -// const Opm::BlackoilState& state, -// const int step, -// const std::string& output_dir) -// { -// Opm::DataMap dm; -// dm["saturation"] = &state.saturation(); -// dm["pressure"] = &state.pressure(); -// dm["surfvolume"] = &state.surfacevol(); -// std::vector cell_velocity; -// Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); -// dm["velocity"] = &cell_velocity; -// -// // Write data (not grid) in Matlab format -// for (Opm::DataMap::const_iterator it = dm.begin(); it != dm.end(); ++it) { -// std::ostringstream fname; -// fname << output_dir << "/" << it->first; -// boost::filesystem::path fpath = fname.str(); -// try { -// create_directories(fpath); -// } -// catch (...) { -// OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); -// } -// fname << "/" << std::setw(3) << std::setfill('0') << step << ".txt"; -// std::ofstream file(fname.str().c_str()); -// if (!file) { -// OPM_THROW(std::runtime_error, "Failed to open " << fname.str()); -// } -// file.precision(15); -// const std::vector& d = *(it->second); -// std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); -// } -// } +#ifdef DISABLE_OUTPUT + static void outputStateMatlab(const UnstructuredGrid& grid, + const Opm::BlackoilState& state, + const int step, + const std::string& output_dir) + { + Opm::DataMap dm; + dm["saturation"] = &state.saturation(); + dm["pressure"] = &state.pressure(); + dm["surfvolume"] = &state.surfacevol(); + std::vector cell_velocity; + Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); + dm["velocity"] = &cell_velocity; + + // Write data (not grid) in Matlab format + for (Opm::DataMap::const_iterator it = dm.begin(); it != dm.end(); ++it) { + std::ostringstream fname; + fname << output_dir << "/" << it->first; + boost::filesystem::path fpath = fname.str(); + try { + create_directories(fpath); + } + catch (...) { + OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); + } + fname << "/" << std::setw(3) << std::setfill('0') << step << ".txt"; + std::ofstream file(fname.str().c_str()); + if (!file) { + OPM_THROW(std::runtime_error, "Failed to open " << fname.str()); + } + file.precision(15); + const std::vector& d = *(it->second); + std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); + } + } +#endif static void outputWaterCut(const Opm::Watercut& watercut, @@ -351,9 +357,13 @@ namespace Opm timer.report(std::cout); if (output_ && (timer.currentStepNum() % output_interval_ == 0)) { if (output_vtk_) { -// outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +#ifdef DISABLE_OUTPUT + outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +#endif } -// outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +#ifdef DISABLE_OUTPUT + outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +#endif } SimulatorReport sreport; @@ -518,9 +528,13 @@ namespace Opm if (output_) { if (output_vtk_) { -// outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +#ifdef DISABLE_OUTPUT + outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +#endif } -// outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +#ifdef DISABLE_OUTPUT + outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +#endif outputWaterCut(watercut, output_dir_); if (wells_) { outputWellReport(wellreport, output_dir_); diff --git a/opm/core/simulator/SimulatorIncompTwophase.cpp b/opm/core/simulator/SimulatorIncompTwophase.cpp index a8efd60b..d55c4573 100644 --- a/opm/core/simulator/SimulatorIncompTwophase.cpp +++ b/opm/core/simulator/SimulatorIncompTwophase.cpp @@ -36,7 +36,9 @@ #include #include // 17.03.2016 Temporarily removed while moving functionality to opm-output -// #include +#ifdef DISABLE_OUTPUT +#include +#endif #include #include @@ -182,35 +184,36 @@ namespace Opm // 17.03.2016 Temporarily removed while moving functionality to opm-output -// static void outputStateVtk(const UnstructuredGrid& grid, -// const Opm::TwophaseState& state, -// const int step, -// const std::string& output_dir) -// { -// // Write data in VTK format. -// std::ostringstream vtkfilename; -// vtkfilename << output_dir << "/vtk_files"; -// boost::filesystem::path fpath(vtkfilename.str()); -// try { -// create_directories(fpath); -// } -// catch (...) { -// OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); -// } -// vtkfilename << "/output-" << std::setw(3) << std::setfill('0') << step << ".vtu"; -// std::ofstream vtkfile(vtkfilename.str().c_str()); -// if (!vtkfile) { -// OPM_THROW(std::runtime_error, "Failed to open " << vtkfilename.str()); -// } -// Opm::DataMap dm; -// dm["saturation"] = &state.saturation(); -// dm["pressure"] = &state.pressure(); -// std::vector cell_velocity; -// Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); -// dm["velocity"] = &cell_velocity; -// Opm::writeVtkData(grid, dm, vtkfile); -// } - +#ifdef DISABLE_OUTPUT + static void outputStateVtk(const UnstructuredGrid& grid, + const Opm::TwophaseState& state, + const int step, + const std::string& output_dir) + { + // Write data in VTK format. + std::ostringstream vtkfilename; + vtkfilename << output_dir << "/vtk_files"; + boost::filesystem::path fpath(vtkfilename.str()); + try { + create_directories(fpath); + } + catch (...) { + OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); + } + vtkfilename << "/output-" << std::setw(3) << std::setfill('0') << step << ".vtu"; + std::ofstream vtkfile(vtkfilename.str().c_str()); + if (!vtkfile) { + OPM_THROW(std::runtime_error, "Failed to open " << vtkfilename.str()); + } + Opm::DataMap dm; + dm["saturation"] = &state.saturation(); + dm["pressure"] = &state.pressure(); + std::vector cell_velocity; + Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); + dm["velocity"] = &cell_velocity; + Opm::writeVtkData(grid, dm, vtkfile); + } +#endif static void outputVectorMatlab(const std::string& name, const std::vector& vec, const int step, @@ -234,39 +237,41 @@ namespace Opm } // 17.03.2016 Temporarily removed while moving functionality to opm-output -// static void outputStateMatlab(const UnstructuredGrid& grid, -// const Opm::TwophaseState& state, -// const int step, -// const std::string& output_dir) -// { -// Opm::DataMap dm; -// dm["saturation"] = &state.saturation(); -// dm["pressure"] = &state.pressure(); -// std::vector cell_velocity; -// Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); -// dm["velocity"] = &cell_velocity; -// -// // Write data (not grid) in Matlab format -// for (Opm::DataMap::const_iterator it = dm.begin(); it != dm.end(); ++it) { -// std::ostringstream fname; -// fname << output_dir << "/" << it->first; -// boost::filesystem::path fpath = fname.str(); -// try { -// create_directories(fpath); -// } -// catch (...) { -// OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); -// } -// fname << "/" << std::setw(3) << std::setfill('0') << step << ".txt"; -// std::ofstream file(fname.str().c_str()); -// if (!file) { -// OPM_THROW(std::runtime_error, "Failed to open " << fname.str()); -// } -// file.precision(15); -// const std::vector& d = *(it->second); -// std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); -// } -// } +#ifdef DISABLE_OUTPUT + static void outputStateMatlab(const UnstructuredGrid& grid, + const Opm::TwophaseState& state, + const int step, + const std::string& output_dir) + { + Opm::DataMap dm; + dm["saturation"] = &state.saturation(); + dm["pressure"] = &state.pressure(); + std::vector cell_velocity; + Opm::estimateCellVelocity(grid, state.faceflux(), cell_velocity); + dm["velocity"] = &cell_velocity; + + // Write data (not grid) in Matlab format + for (Opm::DataMap::const_iterator it = dm.begin(); it != dm.end(); ++it) { + std::ostringstream fname; + fname << output_dir << "/" << it->first; + boost::filesystem::path fpath = fname.str(); + try { + create_directories(fpath); + } + catch (...) { + OPM_THROW(std::runtime_error, "Creating directories failed: " << fpath); + } + fname << "/" << std::setw(3) << std::setfill('0') << step << ".txt"; + std::ofstream file(fname.str().c_str()); + if (!file) { + OPM_THROW(std::runtime_error, "Failed to open " << fname.str()); + } + file.precision(15); + const std::vector& d = *(it->second); + std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); + } + } +#endif static void outputWaterCut(const Opm::Watercut& watercut, @@ -464,9 +469,13 @@ namespace Opm timer.report(*log_); if (output_ && (timer.currentStepNum() % output_interval_ == 0)) { if (output_vtk_) { -// outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +#ifdef DISABLE_OUTPUT + outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +#endif } -// outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +#ifdef DISABLE_OUTPUT + outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +#endif if (use_reorder_) { // This use of dynamic_cast is not ideal, but should be safe. outputVectorMatlab(std::string("reorder_it"), @@ -624,9 +633,13 @@ namespace Opm if (output_) { if (output_vtk_) { -// outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +#ifdef DISABLE_OUTPUT + outputStateVtk(grid_, state, timer.currentStepNum(), output_dir_); +#endif } -// outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +#ifdef DISABLE_OUTPUT + outputStateMatlab(grid_, state, timer.currentStepNum(), output_dir_); +#endif if (use_reorder_) { // This use of dynamic_cast is not ideal, but should be safe. outputVectorMatlab(std::string("reorder_it"), diff --git a/opm/core/simulator/SimulatorOutput.cpp b/opm/core/simulator/SimulatorOutput.cpp index 923c80b3..5c5bae83 100644 --- a/opm/core/simulator/SimulatorOutput.cpp +++ b/opm/core/simulator/SimulatorOutput.cpp @@ -22,7 +22,9 @@ // we need complete definitions for these types #include // 17.03.2016 Temporarily removed while moving functionality to opm-output -// #include +#ifdef DISABLE_OUTPUT +#include +#endif #include #include // partial_sum @@ -47,14 +49,17 @@ SimulatorOutputBase::SimulatorOutputBase ( // process parameters into a writer. we don't setup a new chain in // every timestep! // 17.03.2016 Temporarily removed while moving functionality to opm-output -// , writer_ (std::move (OutputWriter::create (params, eclipseState, phaseUsage, grid))) - +#ifdef DISABLE_OUTPUT + , writer_ (std::move (OutputWriter::create (params, eclipseState, phaseUsage, grid))) +#endif // always start from the first timestep , next_ (0) { // write the static initialization files, even before simulation starts // 17.03.2016 Temporarily removed while moving functionality to opm-output -// writer_->writeInit (*timer); +#ifdef DISABLE_OUTPUT + writer_->writeInit (*timer); +#endif } // default destructor is OK, just need to be defined @@ -92,7 +97,9 @@ SimulatorOutputBase::writeOutput () { // relay the request to the handlers (setup in the constructor // from parameters) // 17.03.2016 Temporarily removed while moving functionality to opm-output -// writer_->writeTimeStep (*timer_, *reservoirState_, *wellState_ , false); +#ifdef DISABLE_OUTPUT + writer_->writeTimeStep (*timer_, *reservoirState_, *wellState_ , false); +#endif // advance to the next reporting time ++next_; diff --git a/opm/core/simulator/SimulatorOutput.hpp b/opm/core/simulator/SimulatorOutput.hpp index c5833aac..ade767c1 100644 --- a/opm/core/simulator/SimulatorOutput.hpp +++ b/opm/core/simulator/SimulatorOutput.hpp @@ -35,7 +35,9 @@ namespace Opm { class Deck; class EclipseState; // 17.03.2016 Temporarily removed while moving functionality to opm-output -//class OutputWriter; +#ifdef DISABLE_OUTPUT +class OutputWriter; +#endif namespace parameter { class ParameterGroup; } class SimulationDataContainer; class SimulatorTimer; @@ -85,8 +87,9 @@ protected: std::shared_ptr wellState_; /// Created locally and destructed together with us -// 17.03.2016 Temporarily removed while moving functionality to opm-output -// std::unique_ptr writer_; +#ifdef DISABLE_OUTPUT + std::unique_ptr writer_; +#endif /// Call the writers that were created based on the parameters virtual void writeOutput (); diff --git a/tutorials/tutorial1.cpp b/tutorials/tutorial1.cpp index e4a57243..43216d17 100644 --- a/tutorials/tutorial1.cpp +++ b/tutorials/tutorial1.cpp @@ -38,7 +38,9 @@ #include #include // 17.03.2016 Temporarily removed while moving functionality to opm-output -//#include +#ifdef DISABLE_OUTPUT +#include +#endif #include #include #include @@ -96,7 +98,9 @@ try /// \snippet tutorial1.cpp data map /// \internal [data map] // 17.03.2016 Temporarily removed while moving functionality to opm-output -// Opm::DataMap dm; +#ifdef DISABLE_OUTPUT + Opm::DataMap dm; +#endif /// \internal [data map] /// \endinternal /// \page tutorial1 @@ -104,7 +108,9 @@ try /// \snippet tutorial1.cpp write vtk /// \internal [write vtk] // 17.03.2016 Temporarily removed while moving functionality to opm-output -// Opm::writeVtkData(*grid.c_grid(), dm, vtkfile); +#ifdef DISABLE_OUTPUT + Opm::writeVtkData(*grid.c_grid(), dm, vtkfile); +#endif /// \internal [write vtk] /// \endinternal } diff --git a/tutorials/tutorial2.cpp b/tutorials/tutorial2.cpp index 445e76a6..6ab81813 100644 --- a/tutorials/tutorial2.cpp +++ b/tutorials/tutorial2.cpp @@ -35,7 +35,9 @@ #include #include // 17.03.2016 Temporarily removed while moving functionality to opm-output -//#include +#ifdef DISABLE_OUTPUT +#include +#endif #include #include #include @@ -188,13 +190,15 @@ try /// \snippet tutorial2.cpp write output /// \internal [write output] // 17.03.2016 Temporarily removed while moving functionality to opm-output -// std::ofstream vtkfile("tutorial2.vtu"); -// Opm::DataMap dm; -// dm["pressure"] = &state.pressure(); -// std::vector cell_velocity; -// Opm::estimateCellVelocity(*grid.c_grid(), state.faceflux(), cell_velocity); -// dm["velocity"] = &cell_velocity; -// Opm::writeVtkData(*grid.c_grid(), dm, vtkfile); +#ifdef DISABLE_OUTPUT + std::ofstream vtkfile("tutorial2.vtu"); + Opm::DataMap dm; + dm["pressure"] = &state.pressure(); + std::vector cell_velocity; + Opm::estimateCellVelocity(*grid.c_grid(), state.faceflux(), cell_velocity); + dm["velocity"] = &cell_velocity; + Opm::writeVtkData(*grid.c_grid(), dm, vtkfile); +#endif /// \internal [write output] /// \endinternal } diff --git a/tutorials/tutorial3.cpp b/tutorials/tutorial3.cpp index 398dbc44..5bdb1958 100644 --- a/tutorials/tutorial3.cpp +++ b/tutorials/tutorial3.cpp @@ -30,7 +30,9 @@ #include #include // 17.03.2016 Temporarily removed while moving functionality to opm-output -//#include +#ifdef DISABLE_OUTPUT +#include +#endif #include #include #include @@ -315,11 +317,13 @@ try vtkfilename.str(""); vtkfilename << "tutorial3-" << std::setw(3) << std::setfill('0') << i << ".vtu"; // 17.03.2016 Temporarily removed while moving functionality to opm-output -// std::ofstream vtkfile(vtkfilename.str().c_str()); -// Opm::DataMap dm; -// dm["saturation"] = &state.saturation(); -// dm["pressure"] = &state.pressure(); -// Opm::writeVtkData(grid, dm, vtkfile); +#ifdef DISABLE_OUTPUT + std::ofstream vtkfile(vtkfilename.str().c_str()); + Opm::DataMap dm; + dm["saturation"] = &state.saturation(); + dm["pressure"] = &state.pressure(); + Opm::writeVtkData(grid, dm, vtkfile); +#endif } } catch (const std::exception &e) { diff --git a/tutorials/tutorial4.cpp b/tutorials/tutorial4.cpp index eee0de49..d1212a02 100644 --- a/tutorials/tutorial4.cpp +++ b/tutorials/tutorial4.cpp @@ -30,7 +30,9 @@ #include #include // 17.03.2016 Temporarily removed while moving functionality to opm-output -//#include +#ifdef DISABLE_OUTPUT +#include +#endif #include #include #include @@ -436,10 +438,12 @@ try vtkfilename << "tutorial4-" << std::setw(3) << std::setfill('0') << i << ".vtu"; std::ofstream vtkfile(vtkfilename.str().c_str()); // 17.03.2016 Temporarily removed while moving functionality to opm-output -// Opm::DataMap dm; -// dm["saturation"] = &state.saturation(); -// dm["pressure"] = &state.pressure(); -// Opm::writeVtkData(grid, dm, vtkfile); +#ifdef DISABLE_OUTPUT + Opm::DataMap dm; + dm["saturation"] = &state.saturation(); + dm["pressure"] = &state.pressure(); + Opm::writeVtkData(grid, dm, vtkfile); +#endif } destroy_wells(wells);