From 87570637115daff6b0aca107feb82731a09597b2 Mon Sep 17 00:00:00 2001 From: Pratik Nayak Date: Fri, 26 Apr 2019 11:23:24 +0200 Subject: [PATCH] Add all solvers supported by Ginkgo. Format with clang-format. Update docs. --- include/deal.II/lac/ginkgo_solver.h | 272 +++++++++++++++++++++++++++- source/lac/ginkgo_solver.cc | 245 ++++++++++++++++++++++--- tests/ginkgo/solver.cc | 47 ++++- tests/ginkgo/solver.output | 5 + 4 files changed, 536 insertions(+), 33 deletions(-) diff --git a/include/deal.II/lac/ginkgo_solver.h b/include/deal.II/lac/ginkgo_solver.h index 861cb047f5..eaca09895e 100644 --- a/include/deal.II/lac/ginkgo_solver.h +++ b/include/deal.II/lac/ginkgo_solver.h @@ -53,7 +53,7 @@ namespace GinkgoWrappers /** * Constructor. * - * The @p executor defines the paradigm where the solution is computed. + * The @p exec_type defines the paradigm where the solution is computed. * Ginkgo currently supports three different executor types: * * + OmpExecutor specifies that the data should be stored and the @@ -226,7 +226,7 @@ namespace GinkgoWrappers * parameters and setup the CG solver from the CG factory which solves the * linear system. * - * @p executor The execution paradigm for the CG solver. + * @p exec_type The execution paradigm for the CG solver. */ SolverCG(SolverControl & solver_control, std::string exec_type, @@ -239,7 +239,9 @@ namespace GinkgoWrappers * parameters and setup the CG solver from the CG factory which solves the * linear system. * - * @p executor The execution paradigm for the CG solver. + * @p exec_type The execution paradigm for the CG solver. + * + * @p preconditioner The preconditioner for the solver. */ SolverCG(SolverControl & solver_control, std::string exec_type, @@ -253,6 +255,270 @@ namespace GinkgoWrappers }; + /** + * An implementation of the solver interface using the Ginkgo BICGSTAB solver. + * + * @ingroup GinkgoWrappers + */ + template + class SolverBICGSTAB : public SolverBase + { + public: + /** + * A standardized data struct to pipe additional data to the solver. + */ + struct AdditionalData + {}; + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the BICGSTAB solver from the BICGSTAB factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the BICGSTAB solver. + */ + SolverBICGSTAB(SolverControl & solver_control, + std::string exec_type, + const AdditionalData & data = AdditionalData()); + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the BICGSTAB solver from the BICGSTAB factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the BICGSTAB solver. + * + * @p preconditioner The preconditioner for the solver. + */ + SolverBICGSTAB(SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data = AdditionalData()); + protected: + /** + * Store a copy of the settings for this particular solver. + */ + const AdditionalData additional_data; + }; + + /** + * An implementation of the solver interface using the Ginkgo CGS solver. + * + * @ingroup GinkgoWrappers + */ + template + class SolverCGS: public SolverBase + { + public: + /** + * A standardized data struct to pipe additional data to the solver. + */ + struct AdditionalData + {}; + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the CGS solver from the CGS factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the CGS solver. + */ + SolverCGS(SolverControl & solver_control, + std::string exec_type, + const AdditionalData & data = AdditionalData()); + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the CGS solver from the CGS factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the CGS solver. + * + * @p preconditioner The preconditioner for the solver. + */ + SolverCGS(SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data = AdditionalData()); + protected: + /** + * Store a copy of the settings for this particular solver. + */ + const AdditionalData additional_data; + }; + + /** + * An implementation of the solver interface using the Ginkgo FCG solver. + * + * @ingroup GinkgoWrappers + */ + template + class SolverFCG: public SolverBase + { + public: + /** + * A standardized data struct to pipe additional data to the solver. + */ + struct AdditionalData + {}; + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the FCG solver from the FCG factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the FCG solver. + */ + SolverFCG(SolverControl & solver_control, + std::string exec_type, + const AdditionalData & data = AdditionalData()); + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the FCG solver from the FCG factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the FCG solver. + * + * @p preconditioner The preconditioner for the solver. + */ + SolverFCG(SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data = AdditionalData()); + protected: + /** + * Store a copy of the settings for this particular solver. + */ + const AdditionalData additional_data; + }; + + /** + * An implementation of the solver interface using the Ginkgo GMRES solver. + * + * @ingroup GinkgoWrappers + */ + template + class SolverGMRES: public SolverBase + { + public: + /** + * A standardized data struct to pipe additional data to the solver. + */ + struct AdditionalData + { + /** + * Constructor. By default, set the number of temporary vectors to 30, + * i.e. do a restart every 30 iterations. + */ + AdditionalData(const unsigned int restart_parameter = 30); + + /** + * Maximum number of tmp vectors. + */ + unsigned int restart_parameter; + + }; + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the GMRES solver from the GMRES factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the GMRES solver. + */ + SolverGMRES(SolverControl & solver_control, + std::string exec_type, + const AdditionalData & data = AdditionalData()); + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the GMRES solver from the GMRES factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the GMRES solver. + * + * @p preconditioner The preconditioner for the solver. + */ + SolverGMRES(SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data = AdditionalData()); + protected: + /** + * Store a copy of the settings for this particular solver. + */ + const AdditionalData additional_data; + }; + + /** + * An implementation of the solver interface using the Ginkgo IR solver. + * + * @ingroup GinkgoWrappers + */ + template + class SolverIR: public SolverBase + { + public: + /** + * A standardized data struct to pipe additional data to the solver. + */ + struct AdditionalData + {}; + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the IR solver from the IR factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the IR solver. + */ + SolverIR(SolverControl & solver_control, + std::string exec_type, + const AdditionalData & data = AdditionalData()); + + /** + * Constructor. + * + * @p solver_control The solver control object is then used to set the + * parameters and setup the IR solver from the IR factory which solves the + * linear system. + * + * @p exec_type The execution paradigm for the IR solver. + * + * @p inner_solver The Inner solver for the IR solver. + */ + SolverIR(SolverControl & solver_control, + std::string exec_type, + std::shared_ptr inner_solver, + const AdditionalData & data = AdditionalData()); + + protected: + /** + * Store a copy of the settings for this particular solver. + */ + const AdditionalData additional_data; + }; + + } // namespace GinkgoWrappers DEAL_II_NAMESPACE_CLOSE diff --git a/source/lac/ginkgo_solver.cc b/source/lac/ginkgo_solver.cc index 4dba46b52b..1e58be0670 100644 --- a/source/lac/ginkgo_solver.cc +++ b/source/lac/ginkgo_solver.cc @@ -31,24 +31,30 @@ DEAL_II_NAMESPACE_OPEN namespace GinkgoWrappers { template - SolverBase::SolverBase( - SolverControl & solver_control, - std::string exec_type) + SolverBase::SolverBase(SolverControl &solver_control, + std::string exec_type) : solver_control(solver_control) , exec_type(exec_type) { - if (exec_type == "reference") { - executor = gko::ReferenceExecutor::create(); - } else if (exec_type == "omp") { - executor = gko::OmpExecutor::create(); - } else if (exec_type == "cuda" && - gko::CudaExecutor::get_num_devices() > 0) { - executor = gko::CudaExecutor::create(0, gko::OmpExecutor::create()); - } else { - std::cerr << "exec_type needs to be one of the three strings: reference, cuda or omp, but provided with" - << exec_type<< std::endl; - std::exit(-1); - } + if (exec_type == "reference") + { + executor = gko::ReferenceExecutor::create(); + } + else if (exec_type == "omp") + { + executor = gko::OmpExecutor::create(); + } + else if (exec_type == "cuda" && gko::CudaExecutor::get_num_devices() > 0) + { + executor = gko::CudaExecutor::create(0, gko::OmpExecutor::create()); + } + else + { + std::cerr + << "exec_type needs to be one of the three strings: reference, cuda or omp, but provided with" + << exec_type << std::endl; + std::exit(-1); + } using ResidualCriterionFactory = gko::stop::ResidualNormReduction<>; residual_criterion = ResidualCriterionFactory::build() .with_reduction_factor(solver_control.tolerance()) @@ -288,10 +294,9 @@ namespace GinkgoWrappers /* ---------------------- SolverCG ------------------------ */ template - SolverCG::SolverCG( - SolverControl & solver_control, - std::string exec_type, - const AdditionalData & data) + SolverCG::SolverCG(SolverControl & solver_control, + std::string exec_type, + const AdditionalData &data) : SolverBase(solver_control, exec_type) , additional_data(data) { @@ -302,17 +307,182 @@ namespace GinkgoWrappers template SolverCG::SolverCG( - SolverControl & solver_control, - std::string exec_type, - std::shared_ptr preconditioner, - const AdditionalData & data) + SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data) : SolverBase(solver_control, exec_type) , additional_data(data) { - using cg = gko::solver::Cg; + using cg = gko::solver::Cg; + this->solver_gen = cg::build() + .with_criteria(this->combined_factory) + .with_preconditioner(preconditioner) + .on(this->executor); + } + + + /* ---------------------- SolverBICGSTAB ------------------------ */ + + template + SolverBICGSTAB::SolverBICGSTAB( + SolverControl & solver_control, + std::string exec_type, + const AdditionalData &data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using bicgstab = gko::solver::Bicgstab; + this->solver_gen = bicgstab::build() + .with_criteria(this->combined_factory) + .on(this->executor); + } + + template + SolverBICGSTAB::SolverBICGSTAB( + SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using bicgstab = gko::solver::Bicgstab; + this->solver_gen = bicgstab::build() + .with_criteria(this->combined_factory) + .with_preconditioner(preconditioner) + .on(this->executor); + } + + /* ---------------------- SolverCGS ------------------------ */ + + template + SolverCGS::SolverCGS(SolverControl &solver_control, + std::string exec_type, + const AdditionalData &data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using cgs = gko::solver::Cgs; + this->solver_gen = + cgs::build().with_criteria(this->combined_factory).on(this->executor); + } + + template + SolverCGS::SolverCGS( + SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using cgs = gko::solver::Cgs; + this->solver_gen = cgs::build() + .with_criteria(this->combined_factory) + .with_preconditioner(preconditioner) + .on(this->executor); + } + + /* ---------------------- SolverFCG ------------------------ */ + + template + SolverFCG::SolverFCG(SolverControl &solver_control, + std::string exec_type, + const AdditionalData &data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using fcg = gko::solver::Fcg; + this->solver_gen = + fcg::build().with_criteria(this->combined_factory).on(this->executor); + } + + template + SolverFCG::SolverFCG( + SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using fcg = gko::solver::Fcg; + this->solver_gen = fcg::build() + .with_criteria(this->combined_factory) + .with_preconditioner(preconditioner) + .on(this->executor); + } + + /* ---------------------- SolverGMRES ------------------------ */ + + template + SolverGMRES::AdditionalData::AdditionalData( + const unsigned int restart_parameter) + : restart_parameter(restart_parameter) + {} + + template + SolverGMRES::SolverGMRES(SolverControl &solver_control, + std::string exec_type, + const AdditionalData &data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using gmres = gko::solver::Gmres; + this->solver_gen = gmres::build() + .with_krylov_dim(additional_data.restart_parameter) + .with_criteria(this->combined_factory) + .on(this->executor); + } + + template + SolverGMRES::SolverGMRES( + SolverControl & solver_control, + std::string exec_type, + std::shared_ptr preconditioner, + const AdditionalData & data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using gmres = gko::solver::Gmres; + this->solver_gen = gmres::build() + .with_krylov_dim(additional_data.restart_parameter) + .with_criteria(this->combined_factory) + .with_preconditioner(preconditioner) + .on(this->executor); + } + + /* ---------------------- SolverIR ------------------------ */ + + template + SolverIR::SolverIR(SolverControl & solver_control, + std::string exec_type, + const AdditionalData &data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using ir = gko::solver::Ir; this->solver_gen = - cg::build().with_criteria(this->combined_factory).with_preconditioner(preconditioner).on(this->executor); + ir::build().with_criteria(this->combined_factory).on(this->executor); + } + + template + SolverIR::SolverIR( + SolverControl & solver_control, + std::string exec_type, + std::shared_ptr inner_solver, + const AdditionalData & data) + : SolverBase(solver_control, exec_type) + , additional_data(data) + { + using ir = gko::solver::Ir; + this->solver_gen = ir::build() + .with_criteria(this->combined_factory) + .with_solver(inner_solver) + .on(this->executor); } + // Explicit instantiations in GinkgoWrappers # define DEALII_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(_macro) \ template _macro(float, int32_t); \ @@ -330,6 +500,31 @@ namespace GinkgoWrappers DEALII_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(DECLARE_SOLVER_CG) # undef DECLARE_SOLVER_CG +# define DECLARE_SOLVER_BICGSTAB(ValueType, IndexType) \ + class SolverBICGSTAB + DEALII_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(DECLARE_SOLVER_BICGSTAB) +# undef DECLARE_SOLVER_BICGSTAB + +# define DECLARE_SOLVER_CGS(ValueType, IndexType) \ + class SolverCGS + DEALII_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(DECLARE_SOLVER_CGS) +# undef DECLARE_SOLVER_CGS + +# define DECLARE_SOLVER_FCG(ValueType, IndexType) \ + class SolverFCG + DEALII_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(DECLARE_SOLVER_FCG) +# undef DECLARE_SOLVER_FCG + +# define DECLARE_SOLVER_GMRES(ValueType, IndexType) \ + class SolverGMRES + DEALII_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(DECLARE_SOLVER_GMRES) +# undef DECLARE_SOLVER_GMRES + +# define DECLARE_SOLVER_IR(ValueType, IndexType) \ + class SolverIR + DEALII_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(DECLARE_SOLVER_IR) +# undef DECLARE_SOLVER_IR + } // namespace GinkgoWrappers diff --git a/tests/ginkgo/solver.cc b/tests/ginkgo/solver.cc index fc0059ce1d..b68c414d3b 100644 --- a/tests/ginkgo/solver.cc +++ b/tests/ginkgo/solver.cc @@ -33,7 +33,7 @@ main(int argc, char **argv) deallog << std::setprecision(4); { - SolverControl control(100, 1e-3); + SolverControl control(200, 1e-3); const unsigned int size = 32; unsigned int dim = (size - 1) * (size - 1); @@ -59,15 +59,52 @@ main(int argc, char **argv) auto executor = "reference"; std::shared_ptr exec = gko::ReferenceExecutor::create(); std::shared_ptr jacobi = gko::preconditioner::Jacobi<>::build().on(exec); - GinkgoWrappers::SolverCG<> solver(control, executor); - GinkgoWrappers::SolverCG<> solver_with_jacobi_precond(control, executor, jacobi); + std::shared_ptr inner_cg = gko::solver::Cg<>::build() + .with_criteria( + gko::stop::Iteration::build().with_max_iters(45u).on(exec), + gko::stop::ResidualNormReduction<>::build() + .with_reduction_factor(1e-5) + .on(exec)) + .on(exec); + GinkgoWrappers::SolverCG<> cg_solver(control, executor); + GinkgoWrappers::SolverBICGSTAB<> bicgstab_solver(control, executor); + GinkgoWrappers::SolverCGS<> cgs_solver(control, executor); + GinkgoWrappers::SolverFCG<> fcg_solver(control, executor); + GinkgoWrappers::SolverGMRES<> gmres_solver(control, executor); + GinkgoWrappers::SolverIR<> ir_solver_cg(control, executor, inner_cg); + GinkgoWrappers::SolverCG<> cg_solver_with_jacobi_precond(control, executor, jacobi); - check_solver_within_range(solver.solve(A, u, f), + check_solver_within_range(cg_solver.solve(A, u, f), control.last_step(), 35, 39); u = 0.; - check_solver_within_range(solver_with_jacobi_precond.solve(A, u, f), + check_solver_within_range(bicgstab_solver.solve(A, u, f), + control.last_step(), + 59, + 65); + u = 0.; + check_solver_within_range(cgs_solver.solve(A, u, f), + control.last_step(), + 72, + 79); + u = 0.; + check_solver_within_range(fcg_solver.solve(A, u, f), + control.last_step(), + 33, + 39); + u = 0.; + check_solver_within_range(gmres_solver.solve(A, u, f), + control.last_step(), + 23, + 29); + u = 0.; + check_solver_within_range(ir_solver_cg.solve(A, u, f), + control.last_step(), + 0, + 2); + u = 0.; + check_solver_within_range(cg_solver_with_jacobi_precond.solve(A, u, f), control.last_step(), 29, 33); diff --git a/tests/ginkgo/solver.output b/tests/ginkgo/solver.output index 06d6d30d58..bef1f4d706 100644 --- a/tests/ginkgo/solver.output +++ b/tests/ginkgo/solver.output @@ -1,4 +1,9 @@ DEAL::Size 32 Unknowns 961 DEAL::Solver stopped within 35 - 39 iterations +DEAL::Solver stopped within 59 - 65 iterations +DEAL::Solver stopped within 72 - 79 iterations +DEAL::Solver stopped within 33 - 39 iterations +DEAL::Solver stopped within 23 - 29 iterations +DEAL::Solver stopped within 0 - 2 iterations DEAL::Solver stopped within 29 - 33 iterations -- 2.39.5