|
using | matrix_type = typename std::remove_const<M>::type |
| The matrix type the preconditioner is for.
|
|
using | domain_type = X |
| The domain type of the preconditioner.
|
|
using | range_type = Y |
| The range type of the preconditioner.
|
|
using | field_type = typename X::field_type |
| The field type of the preconditioner.
|
|
using | domain_type_to = typename CudaPreconditionerType::domain_type |
|
using | range_type_to = typename CudaPreconditionerType::range_type |
| The range type of the preconditioner.
|
|
using | field_type_to = typename domain_type_to::field_type |
| The field type of the preconditioner.
|
|
using | block_type = typename domain_type::block_type |
|
using | XTo = Dune::BlockVector<Dune::FieldVector<field_type_to, block_type::dimension>> |
|
using | YTo = Dune::BlockVector<Dune::FieldVector<field_type_to, block_type::dimension>> |
|
using | matrix_type_to |
|
template<class CudaPreconditionerType, class M, class X, class Y, int l = 1>
class Opm::gpuistl::PreconditionerConvertFieldTypeAdapter< CudaPreconditionerType, M, X, Y, l >
Converts the field type (eg.
double to float) to benchmark single precision preconditioners
- Note
- This is not a fast conversion, it is simply meant to benchmark the potential of some preconditioners on consumer grade GPUs where the double precision performance is often artificially limited.
-
In theory this can do any field_type conversion that is meaningful, but it is only tested on double to float conversion.
-
Remember to set the underlying preconditioner with setUnderlyingPreconditioner (should use the matrix from getConvertedMatrix())
-
One could simply change the constructor design by accepting a creator function for the underlying preconditioner. For the current use cases this is however not needed.
To use this, use something like the following code:
#include <opm/simulators/linalg/gpuistl/PreconditionerConvertFieldTypeAdapter.hpp>
#include <opm/simulators/linalg/ParallelOverlappingILU0.hpp>
using XDouble = Dune::BlockVector<Dune::FieldVector<double, 2>>;
using MDouble = Dune::FieldMatrix<double, 2, 2>;
using SpMatrixDouble = Dune::BCRSMatrix<MDouble>;
using XFloat = Dune::BlockVector<Dune::FieldVector<float, 2>>;
using MFloat = Dune::FieldMatrix<float, 2, 2>;
using SpMatrixFloat = Dune::BCRSMatrix<MFloat>;
template<class ParallelInfo>
void applyILU0AsFloat(const MDouble& matrix, const XDouble& x, XDouble& y) {
XDouble, XDouble>;
auto doubleToFloatConverter = DoubleToFloatConverter(matrix);
const auto& convertedMatrix = doubleToFloatConverter.getConvertedMatrix();
auto floatILU0 = std::make_shared<FloatILU0>(convertedMatrix, 0, 1.0, 0);
doubleToFloatConverter.setUnderlyingPreconditioner(floatILU0);
doubleToFloatConverter.apply(x, y);
}
A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as.
Definition ParallelOverlappingILU0.hpp:131
Converts the field type (eg.
Definition PreconditionerConvertFieldTypeAdapter.hpp:86