59 using typename Dune::IterativeSolver<X, X>::domain_type;
60 using typename Dune::IterativeSolver<X, X>::range_type;
61 using typename Dune::IterativeSolver<X, X>::field_type;
62 using typename Dune::IterativeSolver<X, X>::real_type;
63 using typename Dune::IterativeSolver<X, X>::scalar_real_type;
64 static constexpr auto block_size = domain_type::block_type::dimension;
65 using XGPU = Opm::gpuistl::GpuVector<real_type>;
80 Dune::ScalarProduct<X>&
sp,
81 std::shared_ptr<Dune::Preconditioner<X, X>> prec,
86 : Dune::IterativeSolver<X, X>(
op,
sp, *prec, reduction, maxit,
verbose)
87 , m_opOnCPUWithMatrix(
op)
89 , m_underlyingSolver(constructSolver(prec, reduction, maxit,
verbose, comm))
93 "Currently we only support operators of type MatrixAdapter in the CUDA/HIP solver. "
94 "Use --matrix-add-well-contributions=true. "
95 "Using WellModelMatrixAdapter with SolverAdapter is not well-defined so it will not work well, or at all.");
98 virtual void apply(X& x, X&
b,
double reduction, Dune::InverseOperatorResult&
res)
override
104 if (!m_inputBuffer) {
105 m_inputBuffer.reset(
new XGPU(
b.dim()));
106 m_outputBuffer.reset(
new XGPU(x.dim()));
109 m_inputBuffer->copyFromHost(
b);
111 m_outputBuffer->copyFromHost(x);
113 m_underlyingSolver.apply(*m_outputBuffer, *m_inputBuffer, reduction,
res);
116 m_inputBuffer->copyToHost(
b);
117 m_outputBuffer->copyToHost(x);
119 virtual void apply(X& x, X&
b, Dune::InverseOperatorResult&
res)
override
124 if (!m_inputBuffer) {
125 m_inputBuffer.reset(
new XGPU(
b.dim()));
126 m_outputBuffer.reset(
new XGPU(x.dim()));
129 m_inputBuffer->copyFromHost(
b);
131 m_outputBuffer->copyFromHost(x);
133 m_underlyingSolver.apply(*m_outputBuffer, *m_inputBuffer,
res);
136 m_inputBuffer->copyToHost(
b);
137 m_outputBuffer->copyToHost(x);
141 Operator& m_opOnCPUWithMatrix;
142 GpuSparseMatrix<real_type> m_matrix;
150 template <
class Comm>
155 const Comm& communication)
161 auto precAsHolder = std::dynamic_pointer_cast<PreconditionerHolder<X, X>>(prec);
164 "The preconditioner needs to be a CUDA preconditioner (eg. GPUDILU) wrapped in a "
165 "Opm::gpuistl::PreconditionerAdapter wrapped in a "
166 "Opm::gpuistl::GpuBlockPreconditioner. If you are unsure what this means, set "
167 "preconditioner to 'GPUDILU'");
175 "The preconditioner needs to be a CUDA preconditioner (eg. GPUDILU) wrapped in a "
176 "Opm::gpuistl::PreconditionerAdapter wrapped in a "
177 "Opm::gpuistl::GpuBlockPreconditioner. If you are unsure what this means, set "
178 "preconditioner to 'GPUDILU'");
187#if defined(MPIX_CUDA_AWARE_SUPPORT) && MPIX_CUDA_AWARE_SUPPORT
191#if defined(MPIX_CUDA_AWARE_SUPPORT)
199 std::shared_ptr<Opm::gpuistl::GPUSender<real_type, Comm>>
gpuComm;
212 = Dune::OverlappingSchwarzOperator<GpuSparseMatrix<real_type>, XGPU, XGPU,
CudaCommunication>;
215 auto mpiPreconditioner = std::make_shared<GpuBlockPreconditioner<XGPU, XGPU, CudaCommunication>>(
218 auto scalarProduct = std::make_shared<Dune::ParallelScalarProduct<XGPU, CudaCommunication>>(
239 [[
maybe_unused]]
const Dune::Amg::SequentialInformation& communication)
242 return constructSolver(prec, reduction, maxit,
verbose);
255 auto precAsHolder = std::dynamic_pointer_cast<PreconditionerHolder<XGPU, XGPU>>(prec);
258 "The preconditioner needs to be a CUDA preconditioner wrapped in a "
259 "Opm::gpuistl::PreconditionerHolder (eg. GPUDILU).");
263 auto matrixOperator = std::make_shared<Dune::MatrixAdapter<GpuSparseMatrix<real_type>, XGPU, XGPU>>(m_matrix);
264 auto scalarProduct = std::make_shared<Dune::SeqScalarProduct<XGPU>>();
268 std::unique_ptr<XGPU> m_inputBuffer;
269 std::unique_ptr<XGPU> m_outputBuffer;
SolverAdapter(Operator &op, Dune::ScalarProduct< X > &sp, std::shared_ptr< Dune::Preconditioner< X, X > > prec, scalar_real_type reduction, int maxit, int verbose, const Comm &comm)
constructor
Definition SolverAdapter.hpp:79