From e218eafed977631d4950352749b6308da678e3a4 Mon Sep 17 00:00:00 2001 From: wolf Date: Thu, 25 Nov 1999 12:18:46 +0000 Subject: [PATCH] Use default template parameters extensively for Solver*, Precondition*, and MG* classes. Adapt many places in the library to use them. git-svn-id: https://svn.dealii.org/trunk@1946 0785d39b-7218-0410-832d-ea1e28bc413d --- .../common/scripts/forward_declarations.pl | 33 +++++++++++++++++-- deal.II/deal.II/include/multigrid/mg_base.h | 2 +- deal.II/deal.II/source/numerics/base.cc | 8 ++--- deal.II/deal.II/source/numerics/vectors.cc | 16 ++++----- deal.II/lac/include/lac/mgbase.h | 2 +- deal.II/lac/include/lac/precondition.h | 10 ++++-- deal.II/lac/include/lac/precondition_block.h | 6 ++-- .../lac/include/lac/precondition_selector.h | 3 +- deal.II/lac/include/lac/solver.h | 9 +++-- deal.II/lac/include/lac/solver_bicgstab.h | 3 +- deal.II/lac/include/lac/solver_cg.h | 3 +- deal.II/lac/include/lac/solver_gmres.h | 27 +++++++-------- deal.II/lac/include/lac/solver_qmrs.h | 3 +- deal.II/lac/include/lac/solver_richardson.h | 3 +- deal.II/lac/include/lac/solver_selector.h | 3 +- deal.II/lac/include/lac/vector_memory.h | 21 ++++++++---- tests/deal.II/mg.cc | 19 +++++------ tests/deal.II/mglocal.cc | 13 ++++---- tests/deal.II/wave-test-3.cc | 8 ++--- tests/lac/mg.cc | 27 ++++++++------- tests/lac/solver.cc | 22 ++++++------- tests/lac/solver.checked | 2 +- tests/lac/testmatrix.cc | 6 ++-- tests/lac/testmatrix.h | 6 ++-- 24 files changed, 150 insertions(+), 105 deletions(-) diff --git a/deal.II/common/scripts/forward_declarations.pl b/deal.II/common/scripts/forward_declarations.pl index 25703b19dc..3e9c055971 100644 --- a/deal.II/common/scripts/forward_declarations.pl +++ b/deal.II/common/scripts/forward_declarations.pl @@ -40,16 +40,39 @@ sub parse_class_declarations { open (FILE, $filename); while () { - # if the lines contains a "template" at the # beginning and no semicolon at the end: join it # with the next line. if ( /^\s*template/ && !/;\s*$/ ) { + # if this is a multiline template, then join following lines as well + # find out by looking at the last char of the line + while ( ! />\s*$/ ) { + # more concise check: find out whether the number of '<' is + # not equal to the number of '>' + my ( $tmp1, $tmp2 ) = ($_, $_); + $tmp1 =~ s/[^<]//g; + $tmp2 =~ s/[^>]//g; + if ( length ($tmp1) == length ($tmp2) ) { + last; + } + else + { + s/\n//; + $_ = $_ . " " . ; + s/ \s+/ /g; + } + } s/\n//; + s/ \s+/ /g; $_ = $_ . " " . ; } - if ( /^\s*((template\s*<(([-\w,_\s]|<([-\w_,+\s])+>)+)>\s*)?(class|struct))(.*)/ ) { + if ( /^\s*((template\s*< + ( + ([-\w,_=\s] | + <([-\w_,=\s])+> )+ + )>\s*)? + (class|struct))(.*)/x ) { # this is the declaration of a class, possibly a template $basepart = $1; $rest = $7; @@ -58,7 +81,7 @@ sub parse_class_declarations { # $rest contains the name of the class and what comes after that # # first extract the name of the class - $rest =~ /([\w_]+(\s*<(([-\w,_\s]|<([-\w,+\s])+>)+)>)?)(.*)/; + $rest =~ /([\w_]+(\s*<(([-\w,_\s]|<([-\w,\s])+>)+)>)?)(.*)/; $name = $1; $rest = $6; @@ -75,6 +98,10 @@ sub parse_class_declarations { $declaration =~ s/^\s+//g; $declaration =~ s/\s+$//g; + # strip default template parameters + while ( $declaration =~ s/<(.*)=\s[-\w,_\s]+(<[^.]*>)?(.*)>/<$1 $3 > / ) { + } + # impose a negativ-list of names we do not want to # have in this file. this is due to a compiler bug in # egcs 1.1.2, which does not gracefully handle local diff --git a/deal.II/deal.II/include/multigrid/mg_base.h b/deal.II/deal.II/include/multigrid/mg_base.h index 270eb9e5d6..77ab520edf 100644 --- a/deal.II/deal.II/include/multigrid/mg_base.h +++ b/deal.II/deal.II/include/multigrid/mg_base.h @@ -253,7 +253,7 @@ class MGVector : public Subscriptor * with a subscriptor for smart pointers. * @author Guido Kanschat, 1999 */ -template +template > class MGMatrix : public Subscriptor, public vector { diff --git a/deal.II/deal.II/source/numerics/base.cc b/deal.II/deal.II/source/numerics/base.cc index 4d1d23afee..53f7fe8fa2 100644 --- a/deal.II/deal.II/source/numerics/base.cc +++ b/deal.II/deal.II/source/numerics/base.cc @@ -144,11 +144,11 @@ void ProblemBase::solve () { Assert ((tria!=0) && (dof_handler!=0), ExcNoTriaSelected()); - SolverControl control(4000, 1e-16); - PrimitiveVectorMemory > memory; - SolverCG,Vector > cg(control,memory); + SolverControl control(4000, 1e-16); + PrimitiveVectorMemory<> memory; + SolverCG<> cg(control,memory); - PreconditionRelaxation, Vector > + PreconditionRelaxation<> prec(system_matrix, &SparseMatrix::template precondition_SSOR, 1.2); diff --git a/deal.II/deal.II/source/numerics/vectors.cc b/deal.II/deal.II/source/numerics/vectors.cc index 5332e93e73..041ab2d634 100644 --- a/deal.II/deal.II/source/numerics/vectors.cc +++ b/deal.II/deal.II/source/numerics/vectors.cc @@ -412,11 +412,11 @@ void VectorTools::project (const DoFHandler &dof, MatrixTools::apply_boundary_values (boundary_values, mass_matrix, vec, tmp); - SolverControl control(1000,1e-16); - PrimitiveVectorMemory > memory; - SolverCG,Vector > cg(control,memory); + SolverControl control(1000,1e-16); + PrimitiveVectorMemory<> memory; + SolverCG<> cg(control,memory); - PreconditionRelaxation, Vector > + PreconditionRelaxation<> prec(mass_matrix, &SparseMatrix::template precondition_SSOR, 1.2); @@ -684,11 +684,11 @@ VectorTools::project_boundary_values (const DoFHandler &dof, Vector boundary_projection (rhs.size()); - SolverControl control(1000, 1e-16); - PrimitiveVectorMemory > memory; - SolverCG,Vector > cg(control,memory); + SolverControl control(1000, 1e-16); + PrimitiveVectorMemory<> memory; + SolverCG<> cg(control,memory); - PreconditionRelaxation, Vector > + PreconditionRelaxation<> prec(mass_matrix, &SparseMatrix::template precondition_SSOR, 1.2); diff --git a/deal.II/lac/include/lac/mgbase.h b/deal.II/lac/include/lac/mgbase.h index 80ee32f174..3044152cb6 100644 --- a/deal.II/lac/include/lac/mgbase.h +++ b/deal.II/lac/include/lac/mgbase.h @@ -253,7 +253,7 @@ class MGVector : public Subscriptor * with a subscriptor for smart pointers. * @author Guido Kanschat, 1999 */ -template +template > class MGMatrix : public Subscriptor, public vector { diff --git a/deal.II/lac/include/lac/precondition.h b/deal.II/lac/include/lac/precondition.h index 8ad05245a3..6897d521cf 100644 --- a/deal.II/lac/include/lac/precondition.h +++ b/deal.II/lac/include/lac/precondition.h @@ -59,7 +59,8 @@ class PreconditionIdentity * * @author Guido Kanschat, Wolfgang Bangerth, 1999 */ -template +template, + class VECTOR = Vector > class PreconditionUseMatrix { public: @@ -131,7 +132,8 @@ class PreconditionUseMatrix * * @author Guido Kanschat, Wolfgang Bangerth, 1999 */ -template +template, + class VECTOR = Vector > class PreconditionRelaxation { public: @@ -194,7 +196,9 @@ class PreconditionRelaxation * * @author Guido Kanschat, 1999 */ -template +template, + class PRECONDITION = PreconditionIdentity> class PreconditionLACSolver { public: diff --git a/deal.II/lac/include/lac/precondition_block.h b/deal.II/lac/include/lac/precondition_block.h index 8eb6353cc0..ee42d8f416 100644 --- a/deal.II/lac/include/lac/precondition_block.h +++ b/deal.II/lac/include/lac/precondition_block.h @@ -43,7 +43,8 @@ * store the inverted blocks with less accuracy than the original matrix; * for example, #number==double, inverse_type=float# might be a viable choice. */ -template +template class PreconditionBlock: public Subscriptor { public: @@ -197,7 +198,8 @@ class PreconditionBlock: public Subscriptor * (of arbitray structure) below the diagonal blocks are used * in the #operator ()# function of this class. */ -template +template class PreconditionBlockSOR : public PreconditionBlock { public: diff --git a/deal.II/lac/include/lac/precondition_selector.h b/deal.II/lac/include/lac/precondition_selector.h index 6e113e90c1..72c08e56b3 100644 --- a/deal.II/lac/include/lac/precondition_selector.h +++ b/deal.II/lac/include/lac/precondition_selector.h @@ -70,7 +70,8 @@ * * @author Ralf Hartmann, 1999 */ -template +template , + class Vector = Vector > class PreconditionSelector { public: diff --git a/deal.II/lac/include/lac/solver.h b/deal.II/lac/include/lac/solver.h index 73d868cf27..35e8e73bf5 100644 --- a/deal.II/lac/include/lac/solver.h +++ b/deal.II/lac/include/lac/solver.h @@ -6,16 +6,14 @@ -// forward declaration -class SolverControl; -template class VectorMemory; +#include /** * Base class for iterative solvers. * - * HAS TO BE UPDATED! +//TODO: * HAS TO BE UPDATED! * * This class defines possible * return states of linear solvers and provides interfaces to a memory @@ -124,7 +122,8 @@ template class VectorMemory; * * @author Wolfgang Bangerth, Guido Kanschat, Ralf Hartmann, 1997, 1998, 1999 */ -template +template , + class Vector = Vector > class Solver { public: diff --git a/deal.II/lac/include/lac/solver_bicgstab.h b/deal.II/lac/include/lac/solver_bicgstab.h index a291ac19b6..55ea11569a 100644 --- a/deal.II/lac/include/lac/solver_bicgstab.h +++ b/deal.II/lac/include/lac/solver_bicgstab.h @@ -26,7 +26,8 @@ * has a default argument, so you may call it without the additional * parameter. */ -template +template , + class Vector = Vector > class SolverBicgstab : public Solver { public: diff --git a/deal.II/lac/include/lac/solver_cg.h b/deal.II/lac/include/lac/solver_cg.h index 7a3c383a3f..e71d502691 100644 --- a/deal.II/lac/include/lac/solver_cg.h +++ b/deal.II/lac/include/lac/solver_cg.h @@ -36,7 +36,8 @@ * * @author Original implementation by G. Kanschat, R. Becker and F.-T. Suttmeier, reworking and documentation by Wolfgang Bangerth */ -template +template , + class Vector = Vector > class SolverCG : public Solver { public: diff --git a/deal.II/lac/include/lac/solver_gmres.h b/deal.II/lac/include/lac/solver_gmres.h index c192016029..3d35082310 100644 --- a/deal.II/lac/include/lac/solver_gmres.h +++ b/deal.II/lac/include/lac/solver_gmres.h @@ -60,8 +60,9 @@ * * @author Wolfgang Bangerth */ -template -class SolverGMRES : public Solver +template , + class VECTOR = Vector > +class SolverGMRES : public Solver { public: /** @@ -98,7 +99,7 @@ class SolverGMRES : public Solver * Solver method. */ template - typename Solver::ReturnState solve (const Matrix &A, + typename Solver::ReturnState solve (const MATRIX &A, VECTOR &x, const VECTOR &b, const Preconditioner& precondition); @@ -139,11 +140,11 @@ class SolverGMRES : public Solver /* --------------------- Inline and template functions ------------------- */ -template -SolverGMRES::SolverGMRES (SolverControl &cn, +template +SolverGMRES::SolverGMRES (SolverControl &cn, VectorMemory &mem, const AdditionalData &data) : - Solver (cn,mem), + Solver (cn,mem), additional_data(data) { Assert (additional_data.max_n_tmp_vectors >= 10, @@ -152,10 +153,10 @@ SolverGMRES::SolverGMRES (SolverControl &cn, -template +template inline void -SolverGMRES::givens_rotation (Vector &h, +SolverGMRES::givens_rotation (Vector &h, Vector &b, Vector &ci, Vector &si, @@ -181,10 +182,10 @@ SolverGMRES::givens_rotation (Vector &h, -template +template template -typename Solver::ReturnState -SolverGMRES::solve (const Matrix& A, +typename Solver::ReturnState +SolverGMRES::solve (const MATRIX& A, VECTOR & x, const VECTOR& b, const Preconditioner& precondition) @@ -396,9 +397,9 @@ SolverGMRES::solve (const Matrix& A, -template +template double -SolverGMRES::criterion () +SolverGMRES::criterion () { // dummy implementation. this function is // not needed for the present implementation diff --git a/deal.II/lac/include/lac/solver_qmrs.h b/deal.II/lac/include/lac/solver_qmrs.h index 79cbaf0afd..1dcda2651e 100644 --- a/deal.II/lac/include/lac/solver_qmrs.h +++ b/deal.II/lac/include/lac/solver_qmrs.h @@ -41,7 +41,8 @@ * * @author Guido Kanschat, 1999 */ -template +template , + class Vector = Vector > class SolverQMRS : public Solver { public: diff --git a/deal.II/lac/include/lac/solver_richardson.h b/deal.II/lac/include/lac/solver_richardson.h index 8a529775d7..c1db38f4bf 100644 --- a/deal.II/lac/include/lac/solver_richardson.h +++ b/deal.II/lac/include/lac/solver_richardson.h @@ -31,7 +31,8 @@ * * @author Ralf Hartmann */ -template +template , + class Vector = Vector > class SolverRichardson : public Solver { public: diff --git a/deal.II/lac/include/lac/solver_selector.h b/deal.II/lac/include/lac/solver_selector.h index 000f250763..0f7855dc7e 100644 --- a/deal.II/lac/include/lac/solver_selector.h +++ b/deal.II/lac/include/lac/solver_selector.h @@ -77,7 +77,8 @@ * * @author Ralf Hartmann, 1999 */ -template +template , + class Vector = Vector > class SolverSelector { public: diff --git a/deal.II/lac/include/lac/vector_memory.h b/deal.II/lac/include/lac/vector_memory.h index 6a909b0720..aab9e728c5 100644 --- a/deal.II/lac/include/lac/vector_memory.h +++ b/deal.II/lac/include/lac/vector_memory.h @@ -24,7 +24,7 @@ * sophisticated management of vectors. One of these has to be * applied by the user according to his needs. */ -template +template > class VectorMemory : public Subscriptor { public: @@ -67,7 +67,7 @@ class VectorMemory : public Subscriptor * vectors as needed from the global heap, i.e. performs no * specially adapted actions to the purpose of this class. */ -template +template > class PrimitiveVectorMemory : public VectorMemory { public: @@ -103,8 +103,7 @@ class PrimitiveVectorMemory : public VectorMemory * * @author Guido Kanschat, 1999 */ - -template +template > class GrowingVectorMemory : public VectorMemory { public: @@ -115,6 +114,7 @@ class GrowingVectorMemory : public VectorMemory * of vectors. */ GrowingVectorMemory(unsigned int initial_size = 0); + /** * Destructor. * Release all vectors. @@ -127,6 +127,7 @@ class GrowingVectorMemory : public VectorMemory * are allocated vectors left. */ ~GrowingVectorMemory(); + /** * Return new vector from the pool. */ @@ -146,19 +147,25 @@ class GrowingVectorMemory : public VectorMemory * vector itself. */ typedef pair entry_type; + /** * Array of allocated vectors. */ vector pool; + /** * Overall number of allocations. */ unsigned int n_alloc; }; + + +/* --------------------- inline functions ---------------------- */ + + template -GrowingVectorMemory::GrowingVectorMemory(unsigned int - initial_size) +GrowingVectorMemory::GrowingVectorMemory(const unsigned int initial_size) : pool(initial_size) { for (vector::iterator i=pool.begin();i != pool.end() @@ -170,6 +177,8 @@ GrowingVectorMemory::GrowingVectorMemory(unsigned int n_alloc = 0; } + + template GrowingVectorMemory::~GrowingVectorMemory() { diff --git a/tests/deal.II/mg.cc b/tests/deal.II/mg.cc index b170f61545..79805a5ad9 100644 --- a/tests/deal.II/mg.cc +++ b/tests/deal.II/mg.cc @@ -100,9 +100,9 @@ int main() Vector u; u.reinit(f); - PrimitiveVectorMemory > mem; + PrimitiveVectorMemory<> mem; SolverControl control(100, 1.e-12); - SolverCG, Vector >solver(control, mem); + SolverCG<> solver(control, mem); PreconditionIdentity precondition; solver.solve(A,u,f,precondition); @@ -122,11 +122,10 @@ int main() equation.build_mgmatrix(mgA, mgdof, quadrature); SolverControl cgcontrol(10,0., false, false); - PrimitiveVectorMemory > cgmem; - SolverCG, Vector > cgsolver(cgcontrol, cgmem); + PrimitiveVectorMemory<> cgmem; + SolverCG<> cgsolver(cgcontrol, cgmem); PreconditionIdentity cgprec; - MGCoarseGridLACIteration, Vector >, - SparseMatrix, PreconditionIdentity> + MGCoarseGridLACIteration, SparseMatrix, PreconditionIdentity> coarse(cgsolver, mgA[0], cgprec); MGSmootherLAC smoother(mgA); @@ -135,7 +134,7 @@ int main() MG<2> multigrid(mgdof, hanging_nodes, mgA, transfer); - PreconditionMG, Vector > + PreconditionMG > mgprecondition(multigrid, smoother, smoother, coarse); solver.solve(A, u, f, mgprecondition); @@ -163,9 +162,9 @@ MGSmootherLAC::smooth (const unsigned int level, const Vector &rhs) const { SolverControl control(2,1.e-300,false,false); - PrimitiveVectorMemory > mem; - SolverRichardson , Vector > rich(control, mem); - PreconditionRelaxation , Vector > + PrimitiveVectorMemory<> mem; + SolverRichardson<> rich(control, mem); + PreconditionRelaxation<> prec((*matrices)[level], &SparseMatrix ::template precondition_SSOR, 1.); rich.solve((*matrices)[level], u, rhs, prec); diff --git a/tests/deal.II/mglocal.cc b/tests/deal.II/mglocal.cc index f5bfdda87d..7664abf5cd 100644 --- a/tests/deal.II/mglocal.cc +++ b/tests/deal.II/mglocal.cc @@ -116,9 +116,9 @@ int main() Vector u; u.reinit(f); - PrimitiveVectorMemory > mem; + PrimitiveVectorMemory<> mem; SolverControl control(20, 1.e-12, true); - SolverRichardson, Vector >solver(control, mem); + SolverRichardson<> solver(control, mem); vector mgstruct(tr.n_levels()); MGMatrix > mgA(0,tr.n_levels()-1); @@ -134,11 +134,10 @@ int main() equation.build_mgmatrix(mgA, mgdof, quadrature); SolverControl cgcontrol(20,0., false, false); - PrimitiveVectorMemory > cgmem; - SolverCG, Vector > cgsolver(cgcontrol, cgmem); + PrimitiveVectorMemory<> cgmem; + SolverCG<> cgsolver(cgcontrol, cgmem); PreconditionIdentity cgprec; - MGCoarseGridLACIteration, Vector >, - SparseMatrix, PreconditionIdentity> + MGCoarseGridLACIteration, SparseMatrix, PreconditionIdentity> coarse(cgsolver, mgA[tr.n_levels()-2], cgprec); MGSmootherRelaxation @@ -151,7 +150,7 @@ int main() MG<2> multigrid(mgdof, hanging_nodes, mgA, transfer, tr.n_levels()-2); - PreconditionMG, Vector > + PreconditionMG > mgprecondition(multigrid, smoother, smoother, coarse); u = 0.; diff --git a/tests/deal.II/wave-test-3.cc b/tests/deal.II/wave-test-3.cc index 0650686cce..b9f69f2525 100644 --- a/tests/deal.II/wave-test-3.cc +++ b/tests/deal.II/wave-test-3.cc @@ -5478,13 +5478,13 @@ template unsigned int TimeStep_Wave::solve (const UserMatrix &matrix, Vector &solution, const Vector &rhs) const { - SolverControl control(2000, 1.e-12); - PrimitiveVectorMemory > memory; - SolverCG > pcg(control,memory); + SolverControl control(2000, 1.e-12); + PrimitiveVectorMemory<> memory; + SolverCG pcg(control,memory); // solve pcg.solve (matrix, solution, rhs, - PreconditionUseMatrix > + PreconditionUseMatrix (matrix, &UserMatrix::precondition)); // distribute solution diff --git a/tests/lac/mg.cc b/tests/lac/mg.cc index 3fc185677c..23510ccf16 100644 --- a/tests/lac/mg.cc +++ b/tests/lac/mg.cc @@ -74,10 +74,9 @@ class MGSmootherLAC }; -typedef MGCoarseGridLACIteration , Vector >, -SparseMatrix, /*PreconditionRelaxation ,*/ - PreconditionIdentity > -Coarse; +typedef MGCoarseGridLACIteration, + SparseMatrix, + PreconditionIdentity > Coarse; int main() @@ -98,9 +97,9 @@ int main() FDMGTransfer transfer(maxsize, maxsize, maxlevel); // coarse grid solver - PrimitiveVectorMemory > cgmem; + PrimitiveVectorMemory<> cgmem; ReductionControl cgcontrol(100, 1.e-30, 1.e-2, false, false); - SolverCG , Vector > cgcg(cgcontrol,cgmem); + SolverCG<> cgcg(cgcontrol,cgmem); for (unsigned int level = 0; level <= maxlevel; ++level) @@ -132,7 +131,7 @@ int main() FDMG multigrid(level, A, transfer); -// PreconditionRelaxation , Vector > +// PreconditionRelaxation<> PreconditionIdentity cgprec;//(A[minlevel], &SparseMatrix ::template precondition_SSOR, 1.2); @@ -166,8 +165,8 @@ FDMG::FDMG(unsigned int maxlevel, MGMatrix >& matrices, { for (unsigned int level = minlevel; level<=maxlevel ; ++level) { - s[level].reinit(matrices[level].m()); - d[level].reinit(matrices[level].m()); + solution[level].reinit(matrices[level].m()); + defect[level].reinit(matrices[level].m()); } } @@ -184,13 +183,13 @@ FDMG::level_vmult (unsigned int level, void FDMG::copy_to_mg(const Vector& v) { - d[maxlevel] = v; + defect[maxlevel] = v; } void FDMG::copy_from_mg(Vector& v) { - v = s[maxlevel]; + v = solution[maxlevel]; } @@ -206,9 +205,9 @@ MGSmootherLAC::smooth (const unsigned int level, const Vector &rhs) const { SolverControl control(1,1.e-300,false,false); - PrimitiveVectorMemory > mem; - SolverRichardson , Vector > rich(control, mem); - PreconditionRelaxation , Vector > + PrimitiveVectorMemory<> mem; + SolverRichardson<> rich(control, mem); + PreconditionRelaxation<> prec((*matrices)[level], &SparseMatrix ::template precondition_SSOR, 1.); rich.solve((*matrices)[level], u, rhs, prec); diff --git a/tests/lac/solver.cc b/tests/lac/solver.cc index 8b61521f15..eb6e0d4f48 100644 --- a/tests/lac/solver.cc +++ b/tests/lac/solver.cc @@ -34,14 +34,14 @@ int main() deallog.attach(logfile); deallog.depth_console(0); - GrowingVectorMemory > mem; + GrowingVectorMemory<> mem; SolverControl control(100, 1.e-5); SolverControl verbose_control(100, 1.e-5, true); - SolverCG , Vector > cg(control, mem); - SolverGMRES , Vector > gmres(control, mem,20); - SolverBicgstab , Vector > bicgstab(control, mem); - SolverRichardson , Vector > rich(control, mem); - SolverQMRS , Vector > qmrs(control, mem); + SolverCG<> cg(control, mem); + SolverGMRES<> gmres(control, mem,20); + SolverBicgstab<> bicgstab(control, mem); + SolverRichardson<> rich(control, mem); + SolverQMRS<> qmrs(control, mem); for (unsigned int size=4; size <= 40; size *= 3) { @@ -54,15 +54,15 @@ int main() SparseMatrixStruct structure(dim, dim, 5); testproblem.build_structure(structure); structure.compress(); - SparseMatrix A(structure); + SparseMatrix A(structure); testproblem.laplacian(A); PreconditionIdentity prec_no; - PreconditionRelaxation , Vector > - prec_sor(A, &SparseMatrix ::template precondition_SOR, 1.2); - PreconditionRelaxation , Vector > - prec_ssor(A, &SparseMatrix ::template precondition_SSOR, 1.2); + PreconditionRelaxation<> + prec_sor(A, &SparseMatrix::template precondition_SOR, 1.2); + PreconditionRelaxation<> + prec_ssor(A, &SparseMatrix::template precondition_SSOR, 1.2); Vector f(dim); Vector u(dim); diff --git a/tests/lac/solver.checked b/tests/lac/solver.checked index 01317a3f9c..4def8b6bdf 100644 --- a/tests/lac/solver.checked +++ b/tests/lac/solver.checked @@ -77,7 +77,7 @@ DEAL:no:cg::Starting DEAL:no:cg::Starting DEAL:no:cg::Convergence step 58 DEAL:no:Bicgstab::Starting -DEAL:no:Bicgstab::Convergence step 45 +DEAL:no:Bicgstab::Convergence step 44 DEAL:no:GMRES::Starting DEAL:no:GMRES::Starting DEAL:no:GMRES::Failure step 100 diff --git a/tests/lac/testmatrix.cc b/tests/lac/testmatrix.cc index ac7840fb1a..9450261738 100644 --- a/tests/lac/testmatrix.cc +++ b/tests/lac/testmatrix.cc @@ -221,9 +221,9 @@ FDMGTransfer::prolongate (const unsigned int to_level, void -FDMGTransfer::restrict (const unsigned int from_level, - Vector &dst, - const Vector &src) const +FDMGTransfer::restrict_and_add (const unsigned int from_level, + Vector &dst, + const Vector &src) const { Assert((from_level>0) && (from_level<=matrices.size()), ExcIndexRange(from_level, 0, matrices.size()+1)); diff --git a/tests/lac/testmatrix.h b/tests/lac/testmatrix.h index 13e9ebf7be..885e6db18b 100644 --- a/tests/lac/testmatrix.h +++ b/tests/lac/testmatrix.h @@ -73,9 +73,9 @@ class FDMGTransfer * Implementation of abstract * function in #MGTranferBase#. */ - virtual void restrict (const unsigned int from_level, - Vector &dst, - const Vector &src) const; + virtual void restrict_and_add (const unsigned int from_level, + Vector &dst, + const Vector &src) const; /** * Exception. -- 2.39.5