From 1f39e6a9a966f662fbb84ba23fc07dbc8db1df26 Mon Sep 17 00:00:00 2001 From: Arne Morten Kvarving Date: Mon, 15 Apr 2024 22:38:04 +0200 Subject: [PATCH] BISAI: template Scalar type --- opm/simulators/linalg/bda/opencl/BISAI.cpp | 220 +++++++++++------- opm/simulators/linalg/bda/opencl/BISAI.hpp | 35 +-- .../linalg/bda/opencl/Preconditioner.cpp | 4 +- 3 files changed, 157 insertions(+), 102 deletions(-) diff --git a/opm/simulators/linalg/bda/opencl/BISAI.cpp b/opm/simulators/linalg/bda/opencl/BISAI.cpp index 9fc6df899..c3bc4402c 100644 --- a/opm/simulators/linalg/bda/opencl/BISAI.cpp +++ b/opm/simulators/linalg/bda/opencl/BISAI.cpp @@ -39,18 +39,20 @@ namespace Opm::Accelerator { using Opm::OpmLog; using Dune::Timer; -template -BISAI::BISAI(bool opencl_ilu_parallel_, int verbosity_) +template +BISAI::BISAI(bool opencl_ilu_parallel_, int verbosity_) : Base(verbosity_) { #if CHOW_PATEL OPM_THROW(std::logic_error, "Error --linear-solver=isai cannot be used if ChowPatelIlu is used, probably defined by CMake\n"); #endif - bilu0 = std::make_unique>(opencl_ilu_parallel_, verbosity_); + bilu0 = std::make_unique>(opencl_ilu_parallel_, verbosity_); } -template -void BISAI::setOpencl(std::shared_ptr& context_, std::shared_ptr& queue_) +template +void BISAI:: +setOpencl(std::shared_ptr& context_, + std::shared_ptr& queue_) { context = context_; queue = queue_; @@ -58,7 +60,9 @@ void BISAI::setOpencl(std::shared_ptr& context_, std::s bilu0->setOpencl(context, queue); } -std::vector buildCsrToCscOffsetMap(std::vector colPointers, std::vector rowIndices){ +std::vector +buildCsrToCscOffsetMap(std::vector colPointers, std::vector rowIndices) +{ std::vector aux(colPointers); // colPointers must be copied to this vector std::vector csrToCscOffsetMap(rowIndices.size()); // map must have the same size as the indices vector @@ -74,15 +78,15 @@ std::vector buildCsrToCscOffsetMap(std::vector colPointers, std::vecto return csrToCscOffsetMap; } -template -bool BISAI::analyze_matrix(BlockedMatrix* mat) +template +bool BISAI::analyze_matrix(BlockedMatrix* mat) { return analyze_matrix(mat, nullptr); } -template -bool BISAI::analyze_matrix(BlockedMatrix* mat, - BlockedMatrix* jacMat) +template +bool BISAI:: +analyze_matrix(BlockedMatrix* mat, BlockedMatrix* jacMat) { const unsigned int bs = block_size; auto *m = mat; @@ -103,21 +107,22 @@ bool BISAI::analyze_matrix(BlockedMatrix* mat, } } -template -void BISAI::buildLowerSubsystemsStructures(){ +template +void BISAI::buildLowerSubsystemsStructures() +{ lower.subsystemPointers.assign(Nb + 1, 0); Dune::Timer t_buildLowerSubsystemsStructures; - for(int tcol = 0; tcol < Nb; tcol++){ + for (int tcol = 0; tcol < Nb; tcol++) { int frow = diagIndex[tcol] + 1; int lrow = colPointers[tcol + 1]; int nx = lrow - frow; int nv = 0; - for(int sweep = 0; sweep < nx - 1; sweep++){ - for(int xid = sweep + 1; xid < nx; xid++){ - for(int ptr = diagIndex[rowIndices[frow + sweep]] + 1; ptr < colPointers[rowIndices[frow + sweep + 1]]; ptr++){ + for (int sweep = 0; sweep < nx - 1; sweep++) { + for (int xid = sweep + 1; xid < nx; xid++) { + for( int ptr = diagIndex[rowIndices[frow + sweep]] + 1; ptr < colPointers[rowIndices[frow + sweep + 1]]; ptr++) { if(rowIndices[ptr] == rowIndices[frow + xid]){ lower.nzIndices.push_back(csrToCscOffsetMap[ptr]); lower.knownRhsIndices.push_back(csrToCscOffsetMap[frow + sweep]); @@ -131,29 +136,31 @@ void BISAI::buildLowerSubsystemsStructures(){ lower.subsystemPointers[tcol + 1] = lower.subsystemPointers[tcol] + nv; } - if(verbosity >= 4){ + if (verbosity >= 4) { std::ostringstream out; - out << "BISAI buildLowerSubsystemsStructures time: " << t_buildLowerSubsystemsStructures.stop() << " s"; + out << "BISAI buildLowerSubsystemsStructures time: " + << t_buildLowerSubsystemsStructures.stop() << " s"; OpmLog::info(out.str()); } } -template -void BISAI::buildUpperSubsystemsStructures(){ +template +void BISAI::buildUpperSubsystemsStructures() +{ upper.subsystemPointers.assign(Nb + 1, 0); Dune::Timer t_buildUpperSubsystemsStructures; - for(int tcol = 0; tcol < Nb; tcol++){ + for (int tcol = 0; tcol < Nb; tcol++) { int frow = colPointers[tcol]; int lrow = diagIndex[tcol]; int nx = lrow - frow + 1; int nv = 0; - for(int sweep = 0; sweep < nx - 1; sweep++){ - for(int xid = 0; xid < nx; xid++){ - for(int ptr = colPointers[rowIndices[lrow - sweep]]; ptr < diagIndex[rowIndices[lrow - sweep]]; ptr++){ - if(rowIndices[ptr] == rowIndices[lrow - xid]){ + for (int sweep = 0; sweep < nx - 1; sweep++) { + for (int xid = 0; xid < nx; xid++) { + for (int ptr = colPointers[rowIndices[lrow - sweep]]; ptr < diagIndex[rowIndices[lrow - sweep]]; ptr++) { + if (rowIndices[ptr] == rowIndices[lrow - xid]) { upper.nzIndices.push_back(csrToCscOffsetMap[ptr]); upper.knownRhsIndices.push_back(csrToCscOffsetMap[lrow - sweep]); upper.unknownRhsIndices.push_back(csrToCscOffsetMap[lrow - xid]); @@ -166,17 +173,17 @@ void BISAI::buildUpperSubsystemsStructures(){ upper.subsystemPointers[tcol + 1] = upper.subsystemPointers[tcol] + nv; } - if(verbosity >= 4){ + if (verbosity >= 4) { std::ostringstream out; - out << "BISAI buildUpperSubsystemsStructures time: " << t_buildUpperSubsystemsStructures.stop() << " s"; + out << "BISAI buildUpperSubsystemsStructures time: " + << t_buildUpperSubsystemsStructures.stop() << " s"; OpmLog::info(out.str()); } } -template -bool BISAI:: -create_preconditioner(BlockedMatrix* mat, - BlockedMatrix* jacMat) +template +bool BISAI:: +create_preconditioner(BlockedMatrix* mat, BlockedMatrix* jacMat) { const unsigned int bs = block_size; @@ -199,48 +206,93 @@ create_preconditioner(BlockedMatrix* mat, buildLowerSubsystemsStructures(); buildUpperSubsystemsStructures(); - d_colPointers = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * colPointers.size()); - d_rowIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * rowIndices.size()); - d_csrToCscOffsetMap = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * csrToCscOffsetMap.size()); - d_diagIndex = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * diagIndex.size()); - d_invLvals = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(double) * nnzb * bs * bs); - d_invUvals = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(double) * nnzb * bs * bs); - d_invL_x = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(double) * Nb * bs); - d_lower.subsystemPointers = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * lower.subsystemPointers.size()); - d_upper.subsystemPointers = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * upper.subsystemPointers.size()); + d_colPointers = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * colPointers.size()); + d_rowIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * rowIndices.size()); + d_csrToCscOffsetMap = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * csrToCscOffsetMap.size()); + d_diagIndex = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * diagIndex.size()); + d_invLvals = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(Scalar) * nnzb * bs * bs); + d_invUvals = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(Scalar) * nnzb * bs * bs); + d_invL_x = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(Scalar) * Nb * bs); + d_lower.subsystemPointers = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * lower.subsystemPointers.size()); + d_upper.subsystemPointers = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * upper.subsystemPointers.size()); - if(!lower.nzIndices.empty()){ // knownRhsIndices and unknownRhsIndices will also be empty if nzIndices is empty - d_lower.nzIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * lower.nzIndices.size()); - d_lower.knownRhsIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * lower.knownRhsIndices.size()); - d_lower.unknownRhsIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * lower.unknownRhsIndices.size()); + if (!lower.nzIndices.empty()) { // knownRhsIndices and unknownRhsIndices will also be empty if nzIndices is empty + d_lower.nzIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * lower.nzIndices.size()); + d_lower.knownRhsIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * lower.knownRhsIndices.size()); + d_lower.unknownRhsIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * lower.unknownRhsIndices.size()); } - if(!upper.nzIndices.empty()){ // knownRhsIndices and unknownRhsIndices will also be empty if nzIndices is empty - d_upper.nzIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * upper.nzIndices.size()); - d_upper.knownRhsIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * upper.knownRhsIndices.size()); - d_upper.unknownRhsIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, sizeof(int) * upper.unknownRhsIndices.size()); + if (!upper.nzIndices.empty()) { // knownRhsIndices and unknownRhsIndices will also be empty if nzIndices is empty + d_upper.nzIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * upper.nzIndices.size()); + d_upper.knownRhsIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * upper.knownRhsIndices.size()); + d_upper.unknownRhsIndices = cl::Buffer(*context, CL_MEM_READ_WRITE, + sizeof(int) * upper.unknownRhsIndices.size()); } events.resize(6); - err = queue->enqueueWriteBuffer(d_colPointers, CL_FALSE, 0, colPointers.size() * sizeof(int), colPointers.data(), nullptr, &events[0]); - err |= queue->enqueueWriteBuffer(d_rowIndices, CL_FALSE, 0, rowIndices.size() * sizeof(int), rowIndices.data(), nullptr, &events[1]); - err |= queue->enqueueWriteBuffer(d_csrToCscOffsetMap, CL_FALSE, 0, csrToCscOffsetMap.size() * sizeof(int), csrToCscOffsetMap.data(), nullptr, &events[2]); - err |= queue->enqueueWriteBuffer(d_diagIndex, CL_FALSE, 0, diagIndex.size() * sizeof(int), diagIndex.data(), nullptr, &events[3]); - err |= queue->enqueueWriteBuffer(d_lower.subsystemPointers, CL_FALSE, 0, sizeof(int) * lower.subsystemPointers.size(), lower.subsystemPointers.data(), nullptr, &events[4]); - err |= queue->enqueueWriteBuffer(d_upper.subsystemPointers, CL_FALSE, 0, sizeof(int) * upper.subsystemPointers.size(), upper.subsystemPointers.data(), nullptr, &events[5]); + err = queue->enqueueWriteBuffer(d_colPointers, CL_FALSE, 0, + colPointers.size() * sizeof(int), + colPointers.data(), nullptr, &events[0]); + err |= queue->enqueueWriteBuffer(d_rowIndices, CL_FALSE, 0, + rowIndices.size() * sizeof(int), + rowIndices.data(), nullptr, &events[1]); + err |= queue->enqueueWriteBuffer(d_csrToCscOffsetMap, CL_FALSE, 0, + csrToCscOffsetMap.size() * sizeof(int), + csrToCscOffsetMap.data(), nullptr, &events[2]); + err |= queue->enqueueWriteBuffer(d_diagIndex, CL_FALSE, 0, + diagIndex.size() * sizeof(int), + diagIndex.data(), nullptr, &events[3]); + err |= queue->enqueueWriteBuffer(d_lower.subsystemPointers, CL_FALSE, 0, + sizeof(int) * lower.subsystemPointers.size(), + lower.subsystemPointers.data(), nullptr, &events[4]); + err |= queue->enqueueWriteBuffer(d_upper.subsystemPointers, CL_FALSE, 0, + sizeof(int) * upper.subsystemPointers.size(), + upper.subsystemPointers.data(), nullptr, &events[5]); - if(!lower.nzIndices.empty()){ + if (!lower.nzIndices.empty()) { events.resize(events.size() + 3); - err |= queue->enqueueWriteBuffer(d_lower.nzIndices, CL_FALSE, 0, sizeof(int) * lower.nzIndices.size(), lower.nzIndices.data(), nullptr, &events[events.size() - 3]); - err |= queue->enqueueWriteBuffer(d_lower.knownRhsIndices, CL_FALSE, 0, sizeof(int) * lower.knownRhsIndices.size(), lower.knownRhsIndices.data(), nullptr, &events[events.size() - 2]); - err |= queue->enqueueWriteBuffer(d_lower.unknownRhsIndices, CL_FALSE, 0, sizeof(int) * lower.unknownRhsIndices.size(), lower.unknownRhsIndices.data(), nullptr, &events[events.size() - 1]); + err |= queue->enqueueWriteBuffer(d_lower.nzIndices, CL_FALSE, 0, + sizeof(int) * lower.nzIndices.size(), + lower.nzIndices.data(), nullptr, + &events[events.size() - 3]); + err |= queue->enqueueWriteBuffer(d_lower.knownRhsIndices, CL_FALSE, 0, + sizeof(int) * lower.knownRhsIndices.size(), + lower.knownRhsIndices.data(), nullptr, + &events[events.size() - 2]); + err |= queue->enqueueWriteBuffer(d_lower.unknownRhsIndices, CL_FALSE, 0, + sizeof(int) * lower.unknownRhsIndices.size(), + lower.unknownRhsIndices.data(), nullptr, + &events[events.size() - 1]); } - if(!upper.nzIndices.empty()){ + if (!upper.nzIndices.empty()) { events.resize(events.size() + 3); - err |= queue->enqueueWriteBuffer(d_upper.nzIndices, CL_FALSE, 0, sizeof(int) * upper.nzIndices.size(), upper.nzIndices.data(), nullptr, &events[events.size() - 3]); - err |= queue->enqueueWriteBuffer(d_upper.knownRhsIndices, CL_FALSE, 0, sizeof(int) * upper.knownRhsIndices.size(), upper.knownRhsIndices.data(), nullptr, &events[events.size() - 2]); - err |= queue->enqueueWriteBuffer(d_upper.unknownRhsIndices, CL_FALSE, 0, sizeof(int) * upper.unknownRhsIndices.size(), upper.unknownRhsIndices.data(), nullptr, &events[events.size() - 1]); + err |= queue->enqueueWriteBuffer(d_upper.nzIndices, CL_FALSE, + 0, sizeof(int) * upper.nzIndices.size(), + upper.nzIndices.data(), nullptr, + &events[events.size() - 3]); + err |= queue->enqueueWriteBuffer(d_upper.knownRhsIndices, CL_FALSE, 0, + sizeof(int) * upper.knownRhsIndices.size(), + upper.knownRhsIndices.data(), nullptr, + &events[events.size() - 2]); + err |= queue->enqueueWriteBuffer(d_upper.unknownRhsIndices, CL_FALSE, 0, + sizeof(int) * upper.unknownRhsIndices.size(), + upper.unknownRhsIndices.data(), nullptr, + &events[events.size() - 1]); } cl::WaitForEvents(events); @@ -255,12 +307,14 @@ create_preconditioner(BlockedMatrix* mat, std::tie(d_LUvals, d_invDiagVals) = bilu0->get_preconditioner_data(); events.resize(2); - err = queue->enqueueFillBuffer(d_invLvals, 0, 0, sizeof(double) * nnzb * bs * bs, nullptr, &events[0]); - err |= queue->enqueueFillBuffer(d_invUvals, 0, 0, sizeof(double) * nnzb * bs * bs, nullptr, &events[1]); + err = queue->enqueueFillBuffer(d_invLvals, 0, 0, + sizeof(Scalar) * nnzb * bs * bs, nullptr, &events[0]); + err |= queue->enqueueFillBuffer(d_invUvals, 0, 0, + sizeof(Scalar) * nnzb * bs * bs, nullptr, &events[1]); cl::WaitForEvents(events); events.clear(); - OpenclKernels::isaiL(d_diagIndex, d_colPointers, d_csrToCscOffsetMap, + OpenclKernels::isaiL(d_diagIndex, d_colPointers, d_csrToCscOffsetMap, d_lower.subsystemPointers, d_lower.nzIndices, d_lower.unknownRhsIndices, d_lower.knownRhsIndices, d_LUvals, d_invLvals, Nb); @@ -270,7 +324,7 @@ create_preconditioner(BlockedMatrix* mat, d_upper.knownRhsIndices, d_LUvals, d_invDiagVals, d_invUvals, Nb); - if(verbosity >= 4){ + if (verbosity >= 4) { std::ostringstream out; out << "BISAI createPreconditioner time: " << t_preconditioner.stop() << " s"; OpmLog::info(out.str()); @@ -279,34 +333,34 @@ create_preconditioner(BlockedMatrix* mat, return true; } -template -bool BISAI::create_preconditioner(BlockedMatrix* mat) +template +bool BISAI:: +create_preconditioner(BlockedMatrix* mat) { return create_preconditioner(mat, nullptr); } -template -void BISAI::apply(const cl::Buffer& x, cl::Buffer& y){ +template +void BISAI::apply(const cl::Buffer& x, cl::Buffer& y) +{ const unsigned int bs = block_size; - OpenclKernels::spmv(d_invLvals, d_rowIndices, d_colPointers, + OpenclKernels::spmv(d_invLvals, d_rowIndices, d_colPointers, x, d_invL_x, Nb, bs, true, true); // application of isaiL is a simple spmv with addition // (to compensate for the unitary diagonal that is not // included in isaiL, for simplicity) - OpenclKernels::spmv(d_invUvals, d_rowIndices, d_colPointers, + OpenclKernels::spmv(d_invUvals, d_rowIndices, d_colPointers, d_invL_x, y, Nb, bs); // application of isaiU is a simple spmv } -#define INSTANTIATE_BDA_FUNCTIONS(n) \ -template class BISAI; +#define INSTANCE_TYPE(T) \ + template class BISAI; \ + template class BISAI; \ + template class BISAI; \ + template class BISAI; \ + template class BISAI; \ + template class BISAI; -INSTANTIATE_BDA_FUNCTIONS(1); -INSTANTIATE_BDA_FUNCTIONS(2); -INSTANTIATE_BDA_FUNCTIONS(3); -INSTANTIATE_BDA_FUNCTIONS(4); -INSTANTIATE_BDA_FUNCTIONS(5); -INSTANTIATE_BDA_FUNCTIONS(6); - -#undef INSTANTIATE_BDA_FUNCTIONS +INSTANCE_TYPE(double) } // namespace Opm::Accelerator diff --git a/opm/simulators/linalg/bda/opencl/BISAI.hpp b/opm/simulators/linalg/bda/opencl/BISAI.hpp index 9f251ada2..930ba1947 100644 --- a/opm/simulators/linalg/bda/opencl/BISAI.hpp +++ b/opm/simulators/linalg/bda/opencl/BISAI.hpp @@ -32,10 +32,10 @@ template class BlockedMatrix; /// This class implements a Blocked version of the Incomplete Sparse Approximate Inverse (ISAI) preconditioner. /// Inspired by the paper "Incomplete Sparse Approximate Inverses for Parallel Preconditioning" by Anzt et. al. -template -class BISAI : public Preconditioner +template +class BISAI : public Preconditioner { - using Base = Preconditioner; + using Base = Preconditioner; using Base::N; using Base::Nb; @@ -54,8 +54,8 @@ private: std::vector rowIndices; std::vector diagIndex; std::vector csrToCscOffsetMap; - std::vector invLvals; - std::vector invUvals; + std::vector invLvals; + std::vector invUvals; cl::Buffer d_colPointers; cl::Buffer d_rowIndices; @@ -68,10 +68,10 @@ private: cl::Buffer d_invL_x; bool opencl_ilu_parallel; - std::unique_ptr> bilu0; + std::unique_ptr> bilu0; /// Struct that holds the structure of the small subsystems for each column - typedef struct{ + struct subsystemStructure { /// This vector holds the cumulative sum for the number of non-zero blocks for each subsystem. /// Works similarly to row and column pointers for the CSR and CSC matrix representations. std::vector subsystemPointers; @@ -85,15 +85,15 @@ private: std::vector knownRhsIndices; /// This vector holds the indices of the unknown values of the right hand sides of the subsystems. std::vector unknownRhsIndices; - } subsystemStructure; + }; /// GPU version of subsystemStructure - typedef struct{ + struct subsystemStructureGPU { cl::Buffer subsystemPointers; cl::Buffer nzIndices; cl::Buffer knownRhsIndices; cl::Buffer unknownRhsIndices; - } subsystemStructureGPU; + } ; subsystemStructure lower, upper; subsystemStructureGPU d_lower, d_upper; @@ -110,17 +110,18 @@ public: BISAI(bool opencl_ilu_parallel, int verbosity); // set own Opencl variables, but also that of the bilu0 preconditioner - void setOpencl(std::shared_ptr& context, std::shared_ptr& queue) override; + void setOpencl(std::shared_ptr& context, + std::shared_ptr& queue) override; // analysis, extract parallelism - bool analyze_matrix(BlockedMatrix* mat) override; - bool analyze_matrix(BlockedMatrix* mat, - BlockedMatrix* jacMat) override; + bool analyze_matrix(BlockedMatrix* mat) override; + bool analyze_matrix(BlockedMatrix* mat, + BlockedMatrix* jacMat) override; // ilu_decomposition - bool create_preconditioner(BlockedMatrix* mat) override; - bool create_preconditioner(BlockedMatrix* mat, - BlockedMatrix* jacMat) override; + bool create_preconditioner(BlockedMatrix* mat) override; + bool create_preconditioner(BlockedMatrix* mat, + BlockedMatrix* jacMat) override; // apply preconditioner, x = prec(y) void apply(const cl::Buffer& y, cl::Buffer& x) override; diff --git a/opm/simulators/linalg/bda/opencl/Preconditioner.cpp b/opm/simulators/linalg/bda/opencl/Preconditioner.cpp index 8243a4957..8d16927d5 100644 --- a/opm/simulators/linalg/bda/opencl/Preconditioner.cpp +++ b/opm/simulators/linalg/bda/opencl/Preconditioner.cpp @@ -47,11 +47,11 @@ Preconditioner::create(Type type, bool opencl_ilu_parallel, i { switch (type ) { case Type::BILU0: - return std::make_unique >(opencl_ilu_parallel, verbosity); + return std::make_unique>(opencl_ilu_parallel, verbosity); case Type::CPR: return std::make_unique >(opencl_ilu_parallel, verbosity); case Type::BISAI: - return std::make_unique >(opencl_ilu_parallel, verbosity); + return std::make_unique>(opencl_ilu_parallel, verbosity); } OPM_THROW(std::logic_error,