From 1bc24136a02d3592783ab76649c4590656497c6d Mon Sep 17 00:00:00 2001 From: turcksin Date: Thu, 19 Dec 2013 22:54:11 +0000 Subject: [PATCH] Add more functions to paralution_vector, paralution_solver, and paralution_solver. Modify others functions in deal.II to use ParalutionWrappers. git-svn-id: https://svn.dealii.org/branches/branch_paralution@32066 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/deal.II/lac/paralution_solver.h | 20 +- .../deal.II/lac/paralution_sparse_matrix.h | 53 +- .../include/deal.II/lac/paralution_vector.h | 492 +++++++++++++++++- deal.II/include/deal.II/lac/vector.h | 2 +- .../deal.II/numerics/vector_tools.templates.h | 1 + deal.II/source/fe/fe_tools_interpolate.cc | 1 + deal.II/source/fe/fe_values.cc | 1 + deal.II/source/lac/paralution_solver.cc | 39 +- deal.II/source/lac/solver.cc | 2 +- deal.II/source/lac/vector_memory.cc | 1 + deal.II/source/multigrid/multigrid.cc | 1 + .../numerics/derivative_approximation.cc | 1 + .../source/numerics/dof_output_operator.cc | 1 + deal.II/source/numerics/error_estimator.cc | 1 + deal.II/source/numerics/error_estimator_1d.cc | 1 + deal.II/source/numerics/fe_field_function.cc | 1 + .../source/numerics/point_value_history.cc | 1 + deal.II/source/numerics/solution_transfer.cc | 1 + 18 files changed, 582 insertions(+), 38 deletions(-) diff --git a/deal.II/include/deal.II/lac/paralution_solver.h b/deal.II/include/deal.II/lac/paralution_solver.h index ff52c026ac..1e559b2b6b 100644 --- a/deal.II/include/deal.II/lac/paralution_solver.h +++ b/deal.II/include/deal.II/lac/paralution_solver.h @@ -78,13 +78,15 @@ namespace ParalutionWrappers SolverCG (SolverControl &cn); /** - * Solve the linear system Ax=b using the CG solver of Paralution. + * Solve the linear system Ax=b using the CG solver of + * Paralution. If the flag @p move_to_accelerator is set to true, the + * solver is built on the accelerator. */ - //TODO: add a flag to move to the accelerator template void solve (const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator=false); }; @@ -103,12 +105,14 @@ namespace ParalutionWrappers /** * Solve the linear system Ax=b using the BiCGStab solver of - * Paralution. + * Paralution. If the flag @p move_to_accelerator is set to true, the + * solver is built on the accelerator. */ template void solve (const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator=false); }; @@ -146,12 +150,14 @@ namespace ParalutionWrappers /** * Solve the linear system Ax=b using the GMRES solver of - * Paralution. + * Paralution. If the flag @p move_to_accelerator is set to true, the + * solver is built on the accelerator. */ template void solve (const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator=false); private: /** diff --git a/deal.II/include/deal.II/lac/paralution_sparse_matrix.h b/deal.II/include/deal.II/lac/paralution_sparse_matrix.h index 217d0a90a5..ebc1e12308 100644 --- a/deal.II/include/deal.II/lac/paralution_sparse_matrix.h +++ b/deal.II/include/deal.II/lac/paralution_sparse_matrix.h @@ -54,9 +54,10 @@ namespace ParalutionWrappers { public : /** - * Declare type for container size. + * Declare some of the standard types used. */ - typedef types::global_dof_index size_type; + typedef dealii::types::global_dof_index size_type; + typedef Number value_type; /** * @name 1: Constructors and initialization @@ -304,8 +305,21 @@ namespace ParalutionWrappers //@} /** - * @name 4: Access to underlying Paralution data + * @name 4: Access to underlying Paralution data and move data to + * accelerator/host */ + /** + * Move the SparseMatrix to the accelerator. This function should only be + * called after convert_to_paralution_csr. + */ + void move_to_accelerator(); + + /** + * Move the SparseMatrix to the host. This function should only be called + * after convert_to_paralution_csr. + */ + void move_to_host(); + /** * Return a constant reference to the underlying Paralution LocalMatrix * data. @@ -319,6 +333,12 @@ namespace ParalutionWrappers //@} private : + /** + * This flag is true if local_matrix is used. It becomes true after + * convert_to_paralution_csr has been called. + */ + bool is_local_matrix; + /** * Underlying Paralution LocalMatrix. */ @@ -335,7 +355,10 @@ namespace ParalutionWrappers // ------------------- inline functions -------------- template - inline SparseMatrix::SparseMatrix() {} + inline SparseMatrix::SparseMatrix() + : + is_local_matrix(false) + {} @@ -360,6 +383,7 @@ namespace ParalutionWrappers { local_matrix.Clear(); sparse_matrix.reinit(sparsity_pattern); + is_local_matrix = false; } @@ -368,6 +392,7 @@ namespace ParalutionWrappers { local_matrix.clear(); sparse_matrix.clear(); + is_local_matrix = false; } @@ -375,7 +400,7 @@ namespace ParalutionWrappers template inline typename SparseMatrix::size_type SparseMatrix::m() const { - return local_matrix.get_nrows(); + return (is_local_matrix ? local_matrix.get_nrow() : sparse_matrix.m()); } @@ -383,7 +408,7 @@ namespace ParalutionWrappers template inline typename SparseMatrix::size_type SparseMatrix::n() const { - return local_matrix.get_ncols(); + return (is_local_matrix ? local_matrix.get_ncol() : sparse_matrix.n()); } @@ -479,6 +504,22 @@ namespace ParalutionWrappers + template + inline void SparseMatrix::move_to_accelerator() + { + local_matrix.MoveToAccelerator(); + } + + + + template + inline void SparseMatrix::move_to_host() + { + local_matrix.MoveToHost(); + } + + + template inline paralution::LocalMatrix const &SparseMatrix::paralution_matrix() const { diff --git a/deal.II/include/deal.II/lac/paralution_vector.h b/deal.II/include/deal.II/lac/paralution_vector.h index e743d8929c..366c1ee724 100644 --- a/deal.II/include/deal.II/lac/paralution_vector.h +++ b/deal.II/include/deal.II/lac/paralution_vector.h @@ -22,7 +22,10 @@ #ifdef DEAL_II_WITH_PARALUTION #include +#include #include +#include +#include #include "paralution.hpp" @@ -57,19 +60,30 @@ namespace ParalutionWrappers { public : /** - * Declare some the standard types used in all containers. These types + * Declare some of the standard types used in all containers. These types * parallel those in the C standard libraries * vector<...> class. */ typedef dealii::types::global_dof_index size_type; + typedef Number value_type; typedef Number *iterator; typedef const Number *const_iterator; + /** + * A variable that indicates whether this vector supports distributed data + * storage. If true, then this vector also needs an appropriate compress() + * function that allows communicating recent set or add operations to + * individual elements to be communicated to other processors. + * + * For the current class, the variable equals false, since it does not + * support parallel data storage. + */ + static const bool supports_distributed_data = false; + /** * @name 1: Basic Object-handling */ //@{ - /** * Default constructor that generates and empy (zero size) vector. The * function reinit() will have to give the vector the correct @@ -114,6 +128,23 @@ namespace ParalutionWrappers */ ~Vector(); + /** + * This function does nothing but is there for compatibility with the + * @p PETScWrappers::Vector class. + * + * For the PETSC vector wrapper class, this function compresses the + * underlying representation of the PETSc object, i.e. flusehes the + * buffers of the vector obkect if it has any. This function is necessary + * after writing into a vector element-by-element and before anything else + * can be done on it. + * + * However, for the implementation of this class, it is immaterial and + * thus an empty function. + */ + void compress (::dealii::VectorOperation::values operation + =::dealii::VectorOperation::unknown) const; + + /** * Change the dimension of the vector to @p N. The vector is filled with * zeros. @@ -160,11 +191,68 @@ namespace ParalutionWrappers */ const_iterator end() const; + /** + * Set all components of the vector to the given number s. + */ + Vector& operator= (const Number s); + + /** + * Copy the given vector in the present one. Resize if necessary. + */ + Vector& operator= (const Vector &v); + + /** + * Mean value of the elements of this vector. + */ + Number mean_value() const; + /** * $l_2$-norm of the vector. The ssquare root of the sum of the squares of * the elements. */ Number l2_norm() const; + + /** + * Return true if the vector contains ghost elements. Since this not a + * distributed vector the method always returns false. + */ + bool has_ghost_elements() const; + + /** + * Returns true if the given global index is in the local range of this + * processor. Since this is not a distributed vector the method always + * returns true. + */ + bool in_local_range (const size_type global_index) const; + + /** + * Return an index set that describes wich elements of this vector are + * owned by the current processor. Note that this index set does not + * include elements this vector may store locally as ghost elements but + * that are in fact owned by another processor. As a consequence, the + * index sets returned returned on different processors if this is a + * distributed vector will form disjoint sets that add up to the the + * complete index set. Obviously, if a vector is created on only one + * processor, then the result would satisfy + * @code + * vec.locally_owned_elements() == complete_index_set(vec.size()) + * @endcode + * + * Since the current data type does not support parallel data storage + * across different processors, the returned index set is the complete + * index set. + */ + IndexSet locally_owned_elements() const; + + /** + * Move the Vector to the accelerator. + */ + void move_to_accelerator(); + + /** + * Move the Vector to the host. + */ + void move_to_host(); //@} /** @@ -196,6 +284,24 @@ namespace ParalutionWrappers */ Number &operator[] (const size_type i); + /** + * A collective get operation: instead of getting individual elements of a + * vector, this function allows to get a whole set of elements at once. + * The indices of the elements to be read are stated in the first + * argument, the corresponding values are returned in the second. + */ + void extract_subvector_to (const std::vector &indices, + std::vector &values) const; + + /** + * Just as the above, but with pointers. Useful in minimizing copying of + * the data around. + */ + template + void extract_subvector_to (ForwardIterator indices_begin, + const ForwardIterator indices_end, + OutputIterator values_begin) const; + /** * Return a constant reference to the underlying Paralution LocalVector * class. @@ -223,17 +329,78 @@ namespace ParalutionWrappers */ Vector& operator-= (const Vector &v); + /** + * A collective add operation: This function adds a whole set of values + * stored in @p values to the vector components specified by @p indices. + */ + template + void add (const std::vector &indices, + const std::vector &values); + + /** + * This is a second collective add operation. As a difference, this + * function takes a deal.II vector of values. + */ + template + void add (const std::vector &indices, + const ::dealii::Vector &values); + + /** + * Take an address where n_eleements are stored contiguously and + * add them into the vector. Handles all case which are not covered by the + * other two add() functions above. + */ + template + void add (const size_type n_elements, + const size_type *indices, + const Number2 *values); + /** * Addition of @p s to all components. Note that @p s is a scalar and * not a vector. */ - void add(const Number s); + void add (const Number s); + + /** + * Simple vector addition, equal to the operator +=. + */ + void add (const Vector &v); /** * Simple addition of a multiple of a vector, i.e. *this += a*V. */ - void add(const Number a, - const Vector &v); + void add (const Number a, const Vector &V); + + /** + * Scaling and simple vector addition, i.e. *this += a*V + v*W. + */ + void add (const Number a, const Vector &V, + const Number b, const Vector &W); + + /** + * Scaling and simple vector addition, i.e. *this = s*(*this)+V. + */ + void sadd (const Number s, const Vector &V); + + /** + * Scaling and simple addition, i.e. *this = s*(*this)+a*V. + */ + void sadd (const Number s, const Number a, const Vector &V); + + /** + * Scaling and multiple addition. + * *this = s*(*this) + a*v + b*W. + */ + void sadd (const Number s, const Number a, const Vector &V, + const Number b, const Vector &W); + + /** + * Scaling and multiple addition. + * *this = s*(*this) + a*v + b*W + c*X. + */ + void sadd (const Number s, const Number a, const Vector &V, + const Number v, const Vector &W, const Number c, + const Vector &X); /** * Scale each element of the vector by a constant value. @@ -244,6 +411,35 @@ namespace ParalutionWrappers * Scale each element of the vector by the inverse of the given value. */ Vector& operator/= (const Number factor); + + /** + * Assignement *this = a*u. + */ + void equ (const Number a, const Vector &u); + + /** + * Assignement *this = a*u + b*v. + */ + void equ (const Number a, const Vector &u, + const Number b, const Vector &v); + + /** + * Assignement *this = a*u + b*v + c*w. + */ + void equ (const Number a, const Vector &u, + const Number b, const Vector &v, + const Number c, const Vector &w); + //@} + + /** + * @name 4: Mixed stuff + */ + //@{ + /** + * Determine an estimate for the memory consumption (in bytes) of this + * object. + */ + std::size_t memory_consumption () const; //@} private : @@ -285,7 +481,9 @@ namespace ParalutionWrappers local_vector.Clear(); } - + template + inline void Vector:: compress(::dealii::VectorOperation::values operation) const + {} template void Vector::reinit(const size_type n) @@ -345,6 +543,42 @@ namespace ParalutionWrappers + template + inline Vector& Vector::operator= (const Number s) + { + Assert(numbers::is_finite(s),ExcNumberNotFinite()); + + local_vector.SetValues(s); + + return *this; + } + + + + template + inline Vector& Vector::operator= (const Vector &v) + { + local_vector.CopyFrom(v.paralution_vector()); + + return *this; + } + + + + template + inline Number Vector::mean_value() const + { + Number mean(0.); + unsigned int i_max(size()); + for (unsigned int i=0; i inline Number Vector::l2_norm() const { @@ -353,6 +587,44 @@ namespace ParalutionWrappers + template + inline bool Vector::has_ghost_elements() const + { + return false; + } + + template + inline bool Vector::in_local_range(const size_type global_index) const + { + return true; + } + + + + template + inline IndexSet Vector::locally_owned_elements() const + { + return complete_index_set(size()); + } + + + + template + inline void Vector::move_to_accelerator() + { + local_vector.MoveToAccelerator(); + } + + + + template + inline void Vector::move_to_host() + { + local_vector.MoveToHost(); + } + + + template inline Number Vector::operator() (const size_type i) const { @@ -393,6 +665,32 @@ namespace ParalutionWrappers + template + inline void Vector::extract_subvector_to (const std::vector &indices, + std::vector &values) const + { + for (size_type i=0; i + template + inline void Vector::extract_subvector_to (ForwardIterator indices_begin, + const ForwardIterator indices_end, + OutputIterator values_begin) const + { + while (indices_begin != indices_end) + { + *values_begin = operator()(*indices_begin); + ++indices_begin; + ++values_begin; + } + } + + + template inline paralution::LocalVector const &Vector::paralution_vector() const { @@ -414,7 +712,7 @@ namespace ParalutionWrappers { Assert(size()==v.size(),ExcDimensionMismatch(size(),v.size())); - local_vector.ScaleAdd(1.,v.paralution_vector()); + local_vector.AddScale(v.paralution_vector(),1.); return *this; } @@ -426,13 +724,53 @@ namespace ParalutionWrappers { Assert(size()==v.size(),ExcDimensionMismatch(size(),v.size())); - local_vector.ScaleAddScale(1.,v.paralution_vector(),-1.); + local_vector.AddScale(v.paralution_vector(),-1.); return *this; } + template + template + inline void Vector::add (const std::vector &indices, + const std::vector &values) + { + Assert(indices.size()==values.size(), + ExcDimensionMismatch(indices.size(),values.size())); + add(indices.size(),&indices[0],&values[0]); + } + + + + template + template + inline void Vector::add (const std::vector &indices, + const ::dealii::Vector &values) + { + Assert(indices.size()==values.size(), + ExcDimensionMismatch(indices.size(),values.size())); + add(indices.size(),&indices[0],values.begin()); + } + + + + template + template + inline void Vector::add (const size_type n_indices,const size_type *indices, + const Number2 *values) + { + for (size_type i=0; i inline void Vector::add(const Number s) { @@ -445,12 +783,92 @@ namespace ParalutionWrappers template - inline void Vector::add(const Number a,const Vector &v) + inline void Vector::add(const Vector &V) + { + Assert(size()==V.size(),ExcDimensionMismatch(size(),V.size())); + local_vector.AddScale(V.paralution_vector(),1.); + } + + + + template + inline void Vector::add(const Number a,const Vector &V) { Assert(numbers::is_finite(a),ExcNumberNotFinite()); - Assert(size()==v.size(),ExcDimensionMismatch(size(),v.size())); + Assert(size()==V.size(),ExcDimensionMismatch(size(),V.size())); - local_vector.ScaleAddScale(1.,v.paralution_vector(),a); + local_vector.ScaleAddScale(1.,V.paralution_vector(),a); + } + + + + template + inline void Vector::add(const Number a,const Vector &V, + const Number b,const Vector &W) + { + Assert(numbers::is_finite(a),ExcNumberNotFinite()); + Assert(numbers::is_finite(b),ExcNumberNotFinite()); + Assert(size()==V.size(),ExcDimensionMismatch(size(),V.size())); + Assert(size()==W.size(),ExcDimensionMismatch(size(),W.size())); + + local_vector.ScaleAdd2(1.,V.paralution_vector(),a,W.paralution_vector(),b); + } + + + + template + inline void Vector::sadd(const Number s, const Vector &V) + { + Assert(numbers::is_finite(s),ExcNumberNotFinite()); + Assert(size()==V.size(),ExcDimensionMismatch(size(),V.size())); + + local_vector.ScaleAdd(s,V.paralution_vector()); + } + + + + template + inline void Vector::sadd(const Number s, const Number a, const Vector &V) + { + Assert(numbers::is_finite(s),ExcNumberNotFinite()); + Assert(numbers::is_finite(a),ExcNumberNotFinite()); + Assert(size()==V.size(),ExcDimensionMismatch(size(),V.size())); + + local_vector.ScaleAddScale(s,V.paralution_vector(),a); + } + + + + template + inline void Vector::sadd(const Number s, const Number a, const Vector &V, + const Number b, const Vector &W) + { + Assert(numbers::is_finite(s),ExcNumberNotFinite()); + Assert(numbers::is_finite(a),ExcNumberNotFinite()); + Assert(numbers::is_finite(b),ExcNumberNotFinite()); + Assert(size()==V.size(),ExcDimensionMismatch(size(),V.size())); + Assert(size()==W.size(),ExcDimensionMismatch(size(),W.size())); + + local_vector.ScaleAdd2(s,V.paralution_vector(),a,W.paralution_vector(),b); + } + + + + template + inline void Vector::sadd(const Number s, const Number a, const Vector &V, + const Number b, const Vector &W, const Number c, + const Vector &X) + { + Assert(numbers::is_finite(s),ExcNumberNotFinite()); + Assert(numbers::is_finite(a),ExcNumberNotFinite()); + Assert(numbers::is_finite(b),ExcNumberNotFinite()); + Assert(numbers::is_finite(c),ExcNumberNotFinite()); + Assert(size()==V.size(),ExcDimensionMismatch(size(),V.size())); + Assert(size()==W.size(),ExcDimensionMismatch(size(),W.size())); + Assert(size()==X.size(),ExcDimensionMismatch(size(),X.size())); + + local_vector.ScaleAdd2(s,V.paralution_vector(),a,W.paralution_vector(),b); + local_vector.AddScale(X.paralution_vector(),c); } @@ -480,6 +898,58 @@ namespace ParalutionWrappers return *this; } + + + + template + inline void Vector::equ (const Number a, const Vector &u) + { + Assert(numbers::is_finite(a), ExcNumberNotFinite()); + Assert(size()==u.size(),ExcDimensionMismatch(size(),u.size())); + + local_vector.ScaleAddScale(0.,u.paralution_vector(),a); + } + + + + template + inline void Vector::equ (const Number a, const Vector &u, + const Number b, const Vector &v) + { + Assert(numbers::is_finite(a), ExcNumberNotFinite()); + Assert(numbers::is_finite(b), ExcNumberNotFinite()); + Assert(size()==u.size(),ExcDimensionMismatch(size(),u.size())); + Assert(size()==v.size(),ExcDimensionMismatch(size(),v.size())); + + local_vector.ScaleAdd2(0.,u.paralution_vector(),a,v.paralution_vector(),b); + } + + + + template + inline void Vector::equ (const Number a, const Vector &u, + const Number b, const Vector &v, + const Number c, const Vector &w) + { + Assert(numbers::is_finite(a), ExcNumberNotFinite()); + Assert(numbers::is_finite(b), ExcNumberNotFinite()); + Assert(numbers::is_finite(c), ExcNumberNotFinite()); + Assert(size()==u.size(),ExcDimensionMismatch(size(),u.size())); + Assert(size()==v.size(),ExcDimensionMismatch(size(),v.size())); + Assert(size()==w.size(),ExcDimensionMismatch(size(),w.size())); + + local_vector.ScaleAdd2(0.,u.paralution_vector(),a,v.paralution_vector(),b); + local_vector.AddScale(w.paralution_vector(),c); + } + + + + template + std::size_t Vector::memory_consumption () const + { + //TODO: This is a very poor approzimation. + return sizeof(*this) + (size() * sizeof(Number)); + } } /*@}*/ diff --git a/deal.II/include/deal.II/lac/vector.h b/deal.II/include/deal.II/lac/vector.h index 8dae825247..d6271d4615 100644 --- a/deal.II/include/deal.II/lac/vector.h +++ b/deal.II/include/deal.II/lac/vector.h @@ -762,7 +762,7 @@ public: /** * A collective add operation: - * This funnction adds a whole + * This function adds a whole * set of values stored in @p * values to the vector * components specified by @p diff --git a/deal.II/include/deal.II/numerics/vector_tools.templates.h b/deal.II/include/deal.II/numerics/vector_tools.templates.h index 757c2a35fe..a6c7aa08d6 100644 --- a/deal.II/include/deal.II/numerics/vector_tools.templates.h +++ b/deal.II/include/deal.II/numerics/vector_tools.templates.h @@ -29,6 +29,7 @@ #include #include #include +#include #include #include #include diff --git a/deal.II/source/fe/fe_tools_interpolate.cc b/deal.II/source/fe/fe_tools_interpolate.cc index 13a0dace97..3eb3271180 100644 --- a/deal.II/source/fe/fe_tools_interpolate.cc +++ b/deal.II/source/fe/fe_tools_interpolate.cc @@ -28,6 +28,7 @@ #include #include #include +#include #include #include #include diff --git a/deal.II/source/fe/fe_values.cc b/deal.II/source/fe/fe_values.cc index 639d937fd9..d650f823e1 100644 --- a/deal.II/source/fe/fe_values.cc +++ b/deal.II/source/fe/fe_values.cc @@ -25,6 +25,7 @@ #include #include #include +#include #include #include #include diff --git a/deal.II/source/lac/paralution_solver.cc b/deal.II/source/lac/paralution_solver.cc index 0a7e6564e8..3bf9e3668b 100644 --- a/deal.II/source/lac/paralution_solver.cc +++ b/deal.II/source/lac/paralution_solver.cc @@ -50,15 +50,18 @@ namespace ParalutionWrappers template void SolverCG::solve(const SparseMatrix &A, Vector &x, - const Vector &b) + const Vector &b, + bool move_to_accelerator) { paralution::CG, paralution::LocalVector,Number> solver; solver.SetOperator(A.paralution_matrix()); // Set absolute tolerance, relative tolerance, divergence tolerance, // maximum number of iterations. - solver.Init(solver_control.tolerance(),1e10,1e10, + solver.Init(solver_control.tolerance(),0.,1.e100, solver_control.max_steps()); + if (move_to_accelerator==true) + solver.MoveToAccelerator(); solver.Build(); solver.Solve(b.paralution_vector(),&(x.paralution_vector())); } @@ -77,14 +80,17 @@ namespace ParalutionWrappers template void SolverBicgstab::solve(const SparseMatrix &A, Vector &x, - const Vector &b) + const Vector &b, + bool move_to_accelerator) { paralution::BiCGStab, paralution::LocalVector,Number> solver; // Set absolute tolerance, relative tolerance, divergence tolerance, // maximum number of iterations. - solver.Init(solver_control.tolerance(),1e10,1e10, + solver.Init(solver_control.tolerance(),0.,1.e100, solver_control.max_steps()); + if (move_to_accelerator==true) + solver.MoveToAccelerator(); solver.Build(); solver.Solve(b.paralution_vector(),&(x.paralution_vector())); } @@ -105,15 +111,18 @@ namespace ParalutionWrappers template void SolverGMRES::solve(const SparseMatrix &A, Vector &x, - const Vector &b) + const Vector &b, + bool move_to_accelerator) { paralution::GMRES, paralution::LocalVector,Number> solver; // Set absolute tolerance, relative tolerance, divergence tolerance, // maximum number of iterations. - solver.Init(solver_control.tolerance(),1e10,1e10, + solver.Init(solver_control.tolerance(),0.,1.e100, solver_control.max_steps()); solver.SetBasisSize(additional_data.restart_parameter); + if (move_to_accelerator==true) + solver.MoveToAccelerator(); solver.Build(); solver.Solve(b.paralution_vector(),&(x.paralution_vector())); } @@ -128,27 +137,33 @@ namespace ParalutionWrappers { template void SolverCG::solve(const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator); template void SolverCG::solve(const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator); template void SolverBicgstab::solve(const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator); template void SolverBicgstab::solve(const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator); template void SolverGMRES::solve(const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator); template void SolverGMRES::solve(const SparseMatrix &A, Vector &x, - const Vector &b); + const Vector &b, + bool move_to_accelerator); } DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/source/lac/solver.cc b/deal.II/source/lac/solver.cc index 145347c6cf..1630185b90 100644 --- a/deal.II/source/lac/solver.cc +++ b/deal.II/source/lac/solver.cc @@ -23,7 +23,7 @@ #include #include #include - +#include #include #include diff --git a/deal.II/source/lac/vector_memory.cc b/deal.II/source/lac/vector_memory.cc index 9ed021a997..711dc53fee 100644 --- a/deal.II/source/lac/vector_memory.cc +++ b/deal.II/source/lac/vector_memory.cc @@ -25,6 +25,7 @@ #include #include #include +#include DEAL_II_NAMESPACE_OPEN diff --git a/deal.II/source/multigrid/multigrid.cc b/deal.II/source/multigrid/multigrid.cc index 2001252900..e151cce1aa 100644 --- a/deal.II/source/multigrid/multigrid.cc +++ b/deal.II/source/multigrid/multigrid.cc @@ -22,6 +22,7 @@ #include #include #include +#include #include #include #include diff --git a/deal.II/source/numerics/derivative_approximation.cc b/deal.II/source/numerics/derivative_approximation.cc index dd2a9ca41d..ece310a4f7 100644 --- a/deal.II/source/numerics/derivative_approximation.cc +++ b/deal.II/source/numerics/derivative_approximation.cc @@ -25,6 +25,7 @@ #include #include #include +#include #include #include #include diff --git a/deal.II/source/numerics/dof_output_operator.cc b/deal.II/source/numerics/dof_output_operator.cc index 740c8e89e7..22829fcfa4 100644 --- a/deal.II/source/numerics/dof_output_operator.cc +++ b/deal.II/source/numerics/dof_output_operator.cc @@ -29,6 +29,7 @@ #include #include #include +#include DEAL_II_NAMESPACE_OPEN diff --git a/deal.II/source/numerics/error_estimator.cc b/deal.II/source/numerics/error_estimator.cc index 83e89a9afa..12a2861a4a 100644 --- a/deal.II/source/numerics/error_estimator.cc +++ b/deal.II/source/numerics/error_estimator.cc @@ -26,6 +26,7 @@ #include #include #include +#include #include #include #include diff --git a/deal.II/source/numerics/error_estimator_1d.cc b/deal.II/source/numerics/error_estimator_1d.cc index 15b65d0bc8..d54175d12b 100644 --- a/deal.II/source/numerics/error_estimator_1d.cc +++ b/deal.II/source/numerics/error_estimator_1d.cc @@ -26,6 +26,7 @@ #include #include #include +#include #include #include #include diff --git a/deal.II/source/numerics/fe_field_function.cc b/deal.II/source/numerics/fe_field_function.cc index d6f176aac0..08df1b38bf 100644 --- a/deal.II/source/numerics/fe_field_function.cc +++ b/deal.II/source/numerics/fe_field_function.cc @@ -29,6 +29,7 @@ #include #include #include +#include DEAL_II_NAMESPACE_OPEN diff --git a/deal.II/source/numerics/point_value_history.cc b/deal.II/source/numerics/point_value_history.cc index c7b849693b..60841e6915 100644 --- a/deal.II/source/numerics/point_value_history.cc +++ b/deal.II/source/numerics/point_value_history.cc @@ -23,6 +23,7 @@ #include #include #include +#include #include diff --git a/deal.II/source/numerics/solution_transfer.cc b/deal.II/source/numerics/solution_transfer.cc index 5e0ac217ed..fb7ef2a30e 100644 --- a/deal.II/source/numerics/solution_transfer.cc +++ b/deal.II/source/numerics/solution_transfer.cc @@ -29,6 +29,7 @@ #include #include #include +#include #include DEAL_II_NAMESPACE_OPEN -- 2.39.5