20#ifndef OPM_WRITESYSTEMMATRIXHELPER_HEADER_INCLUDED
21#define OPM_WRITESYSTEMMATRIXHELPER_HEADER_INCLUDED
23#include <dune/istl/matrixmarket.hh>
25#include <opm/simulators/linalg/MatrixMarketSpecializations.hpp>
32namespace Opm::Helper {
62 template <
class SimulatorType,
class LinalgObjectType,
class Communicator>
68 namespace fs = std::filesystem;
70 const auto output_dir = fs::path {simulator.problem().outputDir()} /
"reports";
76 std::ostringstream
oss;
77 oss <<
"prob_" << simulator.episodeIndex()
78 <<
"_time_" << std::setprecision(15) << std::setw(12) << std::setfill(
'0') << simulator.time()
79 <<
"_nit_" << simulator.model().newtonMethod().numIterations()
85 if (comm !=
nullptr) {
120 template <
class SimulatorType,
class VectorType,
class Communicator>
122 const VectorType& rhs,
124 const Communicator* comm)
126 detail::writeMatrixMarketObject(simulator, rhs,
sysName +
"_vector", comm);
153 template <
class SimulatorType,
class MatrixType,
class Communicator>
155 const MatrixType& matrix,
157 const Communicator* comm)
159 detail::writeMatrixMarketObject(simulator, matrix,
sysName +
"_matrix", comm);
195 template <
class SimulatorType,
class MatrixType,
class VectorType,
class Communicator>
197 const MatrixType& matrix,
198 const VectorType& rhs,
200 const Communicator* comm)
202 writeMatrix(simulator, matrix,
sysName, comm);
203 writeVector(simulator, rhs,
sysName, comm);
237 template <
class SimulatorType,
class MatrixType,
class VectorType,
class Communicator>
239 const MatrixType& matrix,
240 const VectorType& rhs,
241 const Communicator* comm)
243 writeSystem(simulator, matrix, rhs,
"flow_", comm);
277 template <
class SimulatorType,
class MatrixType,
class VectorType,
class Communicator>
279 const MatrixType& matrix,
280 const VectorType& rhs,
281 const Communicator* comm)
283 writeSystem(simulator, matrix, rhs,
"mech_", comm);
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242