From: Wolfgang Bangerth Date: Mon, 28 Nov 2005 17:11:58 +0000 (+0000) Subject: Add additional constructors to the solver classes where one does not have to give... X-Git-Tag: v8.0.0~12837 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=872389ab4cb359a94808740d1745a2bcaaeadbdd;p=dealii.git Add additional constructors to the solver classes where one does not have to give a vector memory object, and a static one is chosen instead. git-svn-id: https://svn.dealii.org/trunk@11798 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index cb26090d19..5d5f7b13f0 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -177,6 +177,24 @@ inconvenience this causes.

lac

    +
  1. + New: All solver classes took an argument of type VectorMemory which they use to allocate + memory for temporary vectors. A clever implementation of this + class might allow reusing temporary vectors, and thus reduce + the overhead of repeatedly allocating and freeing + memory. However, almost all instances of using these classes + used the PrimitiveVectorMemory + class. The solver class have now been changed so that if one + omits the respective argument, they fall back to using such a + memory object. Also, all example programs that did not + specifically use a different memory allocation class, have been + changed to not specify anything at all, and thus fall back to + the default. +
    + (WB 2005/11/23) +

    +
  2. Fixed: The SparseMatrix iterators had a problem when one wrote code like iterator->value()=0 diff --git a/deal.II/examples/doxygen/block_matrix_array.cc b/deal.II/examples/doxygen/block_matrix_array.cc index 03dc2ea189..4a0cec66d8 100644 --- a/deal.II/examples/doxygen/block_matrix_array.cc +++ b/deal.II/examples/doxygen/block_matrix_array.cc @@ -18,7 +18,6 @@ #include #include #include -#include #include #include #include diff --git a/deal.II/examples/doxygen/product_matrix.cc b/deal.II/examples/doxygen/product_matrix.cc index b6ffaa8da6..7e92d864bd 100644 --- a/deal.II/examples/doxygen/product_matrix.cc +++ b/deal.II/examples/doxygen/product_matrix.cc @@ -17,7 +17,6 @@ #include #include #include -#include double Adata[] = diff --git a/deal.II/examples/step-11/step-11.cc b/deal.II/examples/step-11/step-11.cc index 02f1388810..fb0ae9afcd 100644 --- a/deal.II/examples/step-11/step-11.cc +++ b/deal.II/examples/step-11/step-11.cc @@ -22,7 +22,6 @@ #include #include #include -#include #include #include #include @@ -583,8 +582,7 @@ template void LaplaceProblem::solve () { SolverControl solver_control (1000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); PreconditionSSOR<> preconditioner; preconditioner.initialize(system_matrix, 1.2); diff --git a/deal.II/examples/step-12/step-12.cc b/deal.II/examples/step-12/step-12.cc index 761bf736c5..14657c8fa4 100644 --- a/deal.II/examples/step-12/step-12.cc +++ b/deal.II/examples/step-12/step-12.cc @@ -19,7 +19,6 @@ #include #include #include -#include #include #include #include @@ -1431,8 +1430,7 @@ template void DGMethod::solve (Vector &solution) { SolverControl solver_control (1000, 1e-12, false, false); - PrimitiveVectorMemory<> vector_memory; - SolverRichardson<> solver (solver_control, vector_memory); + SolverRichardson<> solver (solver_control); // Here we create the // preconditioner, diff --git a/deal.II/examples/step-13/step-13.cc b/deal.II/examples/step-13/step-13.cc index 85a3c8d0df..365415c99a 100644 --- a/deal.II/examples/step-13/step-13.cc +++ b/deal.II/examples/step-13/step-13.cc @@ -30,7 +30,6 @@ #include #include #include -#include #include #include #include @@ -1513,8 +1512,7 @@ namespace LaplaceSolver Solver::LinearSystem::solve (Vector &solution) const { SolverControl solver_control (1000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); PreconditionSSOR<> preconditioner; preconditioner.initialize(matrix, 1.2); diff --git a/deal.II/examples/step-14/step-14.cc b/deal.II/examples/step-14/step-14.cc index a06c14a179..6126068763 100644 --- a/deal.II/examples/step-14/step-14.cc +++ b/deal.II/examples/step-14/step-14.cc @@ -21,7 +21,6 @@ #include #include #include -#include #include #include #include @@ -761,8 +760,7 @@ namespace LaplaceSolver Solver::LinearSystem::solve (Vector &solution) const { SolverControl solver_control (5000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); PreconditionSSOR<> preconditioner; preconditioner.initialize(matrix, 1.2); diff --git a/deal.II/examples/step-15/step-15.cc b/deal.II/examples/step-15/step-15.cc index 11ba91ff55..31b4c04082 100644 --- a/deal.II/examples/step-15/step-15.cc +++ b/deal.II/examples/step-15/step-15.cc @@ -21,7 +21,6 @@ #include #include #include -#include #include #include #include @@ -804,8 +803,7 @@ void MinimizationProblem::do_step () { SolverControl solver_control (residual.size(), 1e-2*residual.l2_norm()); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> solver (solver_control, vector_memory); + SolverCG<> solver (solver_control); PreconditionSSOR<> preconditioner; preconditioner.initialize(matrix); diff --git a/deal.II/examples/step-16/step-16.cc b/deal.II/examples/step-16/step-16.cc index 06ca843a16..dde5cf2cea 100644 --- a/deal.II/examples/step-16/step-16.cc +++ b/deal.II/examples/step-16/step-16.cc @@ -21,7 +21,6 @@ #include #include #include -#include #include #include #include @@ -393,8 +392,8 @@ void LaplaceProblem::solve () // Create a memory handler for // regular vectors. Note, that // GrowingVectorMemory is more time - // efficient than - // PrimitiveVectorMemory. + // efficient than the + // PrimitiveVectorMemory class. GrowingVectorMemory<> vector_memory; // Now, create an object handling @@ -466,7 +465,7 @@ void LaplaceProblem::solve () // Finally, create the solver // object and solve the system SolverControl solver_control (1000, 1e-12); - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); cg.solve (system_matrix, solution, system_rhs, diff --git a/deal.II/examples/step-20/step-20.cc b/deal.II/examples/step-20/step-20.cc index e23f8d0bd7..be5230e3b4 100644 --- a/deal.II/examples/step-20/step-20.cc +++ b/deal.II/examples/step-20/step-20.cc @@ -40,7 +40,6 @@ const unsigned int degree = 2; #include #include #include -#include #include #include @@ -545,8 +544,7 @@ class SchurComplement SolverControl solver_control (tmp1.size(), 1e-8*tmp1.l2_norm()); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); PreconditionSSOR<> precondition; precondition.initialize(A.block(0,0)); @@ -582,8 +580,7 @@ void LaplaceProblem::solve () { SolverControl solver_control (system_matrix.block(0,0).m(), 1e-6*system_rhs.block(1).l2_norm()); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); cg.solve (SchurComplement(system_matrix), solution.block(1), system_rhs.block(1), @@ -604,8 +601,7 @@ void LaplaceProblem::solve () SolverControl solver_control (system_matrix.block(0,0).m(), 1e-6*tmp.l2_norm()); - PrimitiveVectorMemory<> vector_memory; - SolverGMRES<> cg (solver_control, vector_memory); + SolverGMRES<> cg (solver_control); cg.solve (system_matrix.block(0,0), solution.block(0), tmp, PreconditionIdentity()); diff --git a/deal.II/examples/step-3/step-3.cc b/deal.II/examples/step-3/step-3.cc index aa83fc01dd..071dcc9350 100644 --- a/deal.II/examples/step-3/step-3.cc +++ b/deal.II/examples/step-3/step-3.cc @@ -79,7 +79,6 @@ #include #include #include -#include #include // Finally, this is for output to a @@ -719,30 +718,14 @@ void LaplaceProblem::solve () // latter criterion will be the one // which stops the iteration. SolverControl solver_control (1000, 1e-12); - // Furthermore, the CG algorithm - // needs some space for temporary - // vectors. Rather than allocating - // it on the stack or heap itself, - // it relies on helper objects, - // which can sometimes do a better - // job at this. The - // PrimitiveVectorMemory class is - // such a helper class which the - // solver can ask for memory. The - // angle brackets ``<>'' indicate - // that this class really takes a - // template parameter (here the - // data type of the vectors we - // use), which however has a - // default value, which is - // appropriate here. - PrimitiveVectorMemory<> vector_memory; // Then we need the solver // itself. The template parameters // here are the matrix type and the - // type of the vectors. They - // default to the ones we use here. - SolverCG<> cg (solver_control, vector_memory); + // type of the vectors, but the + // empty angle brackets indicate + // that we simply take the default + // arguments. + SolverCG<> cg (solver_control); // Now solve the system of // equations. The CG solver takes a diff --git a/deal.II/examples/step-4/step-4.cc b/deal.II/examples/step-4/step-4.cc index 9e1be30a08..9e068f9025 100644 --- a/deal.II/examples/step-4/step-4.cc +++ b/deal.II/examples/step-4/step-4.cc @@ -33,7 +33,6 @@ #include #include #include -#include #include #include @@ -497,8 +496,7 @@ template void LaplaceProblem::solve () { SolverControl solver_control (1000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); cg.solve (system_matrix, solution, system_rhs, PreconditionIdentity()); diff --git a/deal.II/examples/step-5/step-5.cc b/deal.II/examples/step-5/step-5.cc index 93d2697133..749bdcad84 100644 --- a/deal.II/examples/step-5/step-5.cc +++ b/deal.II/examples/step-5/step-5.cc @@ -21,7 +21,6 @@ #include #include #include -#include #include #include #include @@ -536,8 +535,7 @@ template void LaplaceProblem::solve () { SolverControl solver_control (1000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); // The only thing we have to alter // is that we need an object which diff --git a/deal.II/examples/step-6/step-6.cc b/deal.II/examples/step-6/step-6.cc index dda25221d7..d70b1495b4 100644 --- a/deal.II/examples/step-6/step-6.cc +++ b/deal.II/examples/step-6/step-6.cc @@ -22,7 +22,6 @@ #include #include #include -#include #include #include #include @@ -596,8 +595,7 @@ template void LaplaceProblem::solve () { SolverControl solver_control (1000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); PreconditionSSOR<> preconditioner; preconditioner.initialize(system_matrix, 1.2); diff --git a/deal.II/examples/step-7/step-7.cc b/deal.II/examples/step-7/step-7.cc index 01f571c7b6..175a77ac7d 100644 --- a/deal.II/examples/step-7/step-7.cc +++ b/deal.II/examples/step-7/step-7.cc @@ -22,7 +22,6 @@ #include #include #include -#include #include #include #include @@ -1057,8 +1056,7 @@ template void LaplaceProblem::solve () { SolverControl solver_control (1000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); PreconditionSSOR<> preconditioner; preconditioner.initialize(system_matrix, 1.2); diff --git a/deal.II/examples/step-8/step-8.cc b/deal.II/examples/step-8/step-8.cc index 3ae01d2ef1..cdd271daae 100644 --- a/deal.II/examples/step-8/step-8.cc +++ b/deal.II/examples/step-8/step-8.cc @@ -21,7 +21,6 @@ #include #include #include -#include #include #include #include @@ -761,8 +760,7 @@ template void ElasticProblem::solve () { SolverControl solver_control (1000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverCG<> cg (solver_control, vector_memory); + SolverCG<> cg (solver_control); PreconditionSSOR<> preconditioner; preconditioner.initialize(system_matrix, 1.2); diff --git a/deal.II/examples/step-9/step-9.cc b/deal.II/examples/step-9/step-9.cc index 6b9e4c1c61..173d0fcc09 100644 --- a/deal.II/examples/step-9/step-9.cc +++ b/deal.II/examples/step-9/step-9.cc @@ -22,7 +22,6 @@ #include #include #include -#include #include #include #include @@ -1330,8 +1329,7 @@ template void AdvectionProblem::solve () { SolverControl solver_control (1000, 1e-12); - PrimitiveVectorMemory<> vector_memory; - SolverBicgstab<> bicgstab (solver_control, vector_memory); + SolverBicgstab<> bicgstab (solver_control); PreconditionJacobi<> preconditioner; preconditioner.initialize(system_matrix, 1.0); diff --git a/deal.II/lac/include/lac/solver.h b/deal.II/lac/include/lac/solver.h index c4977b6d17..9a9e00e54e 100644 --- a/deal.II/lac/include/lac/solver.h +++ b/deal.II/lac/include/lac/solver.h @@ -15,9 +15,10 @@ #include #include +#include + template class Vector; -template class VectorMemory; class SolverControl; /*!@addtogroup Solvers */ @@ -145,7 +146,7 @@ class Solver : public Subscriptor { public: /** - * Constructor. Assign a control + * Constructor. Takes a control * object which evaluates the * conditions for convergence, * and an object to provide @@ -158,8 +159,28 @@ class Solver : public Subscriptor * least as long as that of the solver * object. */ - Solver (SolverControl &, - VectorMemory &); + Solver (SolverControl &solver_control, + VectorMemory &vector_memory); + + /** + * Constructor. Takes a control + * object which evaluates the + * conditions for convergence. In + * contrast to the other + * constructor, this constructor + * denotes an internal object of + * type PrimitiveVectorMemory to + * allocate memory. + * + * A reference to the control + * object is stored, so it is the + * user's responsibility to + * guarantee that the lifetime of + * the two arguments is at least + * as long as that of the solver + * object. + */ + Solver (SolverControl &solver_control); /** * Access to object that controls @@ -168,7 +189,14 @@ class Solver : public Subscriptor SolverControl & control() const; protected: - + /** + * A static vector memory object + * to be used whenever no such + * object has been given to the + * constructor. + */ + static PrimitiveVectorMemory static_vector_memory; + /** * Control structure. */ @@ -183,6 +211,27 @@ class Solver : public Subscriptor /*@}*/ /*-------------------------------- Inline functions ------------------------*/ +template +inline +Solver::Solver (SolverControl &solver_control, + VectorMemory &vector_memory) + : + cntrl(solver_control), + memory(vector_memory) +{} + + + +template +inline +Solver::Solver (SolverControl &solver_control) + : + cntrl(solver_control), + memory(static_vector_memory) +{} + + + template inline SolverControl & @@ -192,12 +241,6 @@ Solver::control() const } -template -inline -Solver::Solver(SolverControl &cn, VectorMemory &mem) - : cntrl(cn), - memory(mem) -{} #endif diff --git a/deal.II/lac/include/lac/solver_bicgstab.h b/deal.II/lac/include/lac/solver_bicgstab.h index 98627bc055..4198d90c3e 100644 --- a/deal.II/lac/include/lac/solver_bicgstab.h +++ b/deal.II/lac/include/lac/solver_bicgstab.h @@ -102,6 +102,14 @@ class SolverBicgstab : public Solver VectorMemory &mem, const AdditionalData &data=AdditionalData()); + /** + * Constructor. Use an object of + * type PrimitiveVectorMemory as + * a default to allocate memory. + */ + SolverBicgstab (SolverControl &cn, + const AdditionalData &data=AdditionalData()); + /** * Virtual destructor. */ @@ -247,6 +255,16 @@ SolverBicgstab::SolverBicgstab (SolverControl &cn, +template +SolverBicgstab::SolverBicgstab (SolverControl &cn, + const AdditionalData &data) + : + Solver(cn), + additional_data(data) +{} + + + template SolverBicgstab::~SolverBicgstab () {} diff --git a/deal.II/lac/include/lac/solver_cg.h b/deal.II/lac/include/lac/solver_cg.h index 671f678d0d..171731815b 100644 --- a/deal.II/lac/include/lac/solver_cg.h +++ b/deal.II/lac/include/lac/solver_cg.h @@ -98,6 +98,14 @@ class SolverCG : public Solver VectorMemory &mem, const AdditionalData &data = AdditionalData()); + /** + * Constructor. Use an object of + * type PrimitiveVectorMemory as + * a default to allocate memory. + */ + SolverCG (SolverControl &cn, + const AdditionalData &data=AdditionalData()); + /** * Virtual destructor. */ @@ -184,10 +192,11 @@ AdditionalData (const bool log_coefficients) {} + template -SolverCG::SolverCG(SolverControl &cn, - VectorMemory &mem, - const AdditionalData &data) +SolverCG::SolverCG (SolverControl &cn, + VectorMemory &mem, + const AdditionalData &data) : Solver(cn,mem), additional_data(data) @@ -195,6 +204,16 @@ SolverCG::SolverCG(SolverControl &cn, +template +SolverCG::SolverCG (SolverControl &cn, + const AdditionalData &data) + : + Solver(cn), + additional_data(data) +{} + + + template SolverCG::~SolverCG () {} diff --git a/deal.II/lac/include/lac/solver_gmres.h b/deal.II/lac/include/lac/solver_gmres.h index ffe284b13c..bd6d5ef3ec 100644 --- a/deal.II/lac/include/lac/solver_gmres.h +++ b/deal.II/lac/include/lac/solver_gmres.h @@ -231,6 +231,14 @@ class SolverGMRES : public Solver */ SolverGMRES (SolverControl &cn, VectorMemory &mem, + const AdditionalData &data=AdditionalData()); + + /** + * Constructor. Use an object of + * type PrimitiveVectorMemory as + * a default to allocate memory. + */ + SolverGMRES (SolverControl &cn, const AdditionalData &data=AdditionalData()); /** @@ -345,6 +353,14 @@ class SolverFGMRES : public Solver */ SolverFGMRES (SolverControl &cn, VectorMemory &mem, + const AdditionalData &data=AdditionalData()); + + /** + * Constructor. Use an object of + * type PrimitiveVectorMemory as + * a default to allocate memory. + */ + SolverFGMRES (SolverControl &cn, const AdditionalData &data=AdditionalData()); /** @@ -448,15 +464,27 @@ AdditionalData (const unsigned int max_n_tmp_vectors, use_default_residual(use_default_residual) {} + template SolverGMRES::SolverGMRES (SolverControl &cn, VectorMemory &mem, - const AdditionalData &data) : + const AdditionalData &data) + : Solver (cn,mem), additional_data(data) {} + +template +SolverGMRES::SolverGMRES (SolverControl &cn, + const AdditionalData &data) : + Solver (cn), + additional_data(data) +{} + + + template inline void @@ -484,6 +512,7 @@ SolverGMRES::givens_rotation (Vector &h, } + template template void @@ -795,6 +824,7 @@ SolverGMRES::solve (const MATRIX &A, } + template double SolverGMRES::criterion () @@ -815,10 +845,21 @@ SolverFGMRES::SolverFGMRES (SolverControl &cn, const AdditionalData &data) : Solver (cn, mem), - additional_data(data) + additional_data(data) {} + +template +SolverFGMRES::SolverFGMRES (SolverControl &cn, + const AdditionalData &data) + : + Solver (cn), + additional_data(data) +{} + + + template template void diff --git a/deal.II/lac/include/lac/solver_minres.h b/deal.II/lac/include/lac/solver_minres.h index 1385fbe112..1dc07b43cd 100644 --- a/deal.II/lac/include/lac/solver_minres.h +++ b/deal.II/lac/include/lac/solver_minres.h @@ -73,6 +73,14 @@ class SolverMinRes : public Solver */ SolverMinRes (SolverControl &cn, VectorMemory &mem, + const AdditionalData &data=AdditionalData()); + + /** + * Constructor. Use an object of + * type PrimitiveVectorMemory as + * a default to allocate memory. + */ + SolverMinRes (SolverControl &cn, const AdditionalData &data=AdditionalData()); /** @@ -157,6 +165,15 @@ SolverMinRes::SolverMinRes (SolverControl &cn, {} + +template +SolverMinRes::SolverMinRes (SolverControl &cn, + const AdditionalData &) + : + Solver(cn) +{} + + template SolverMinRes::~SolverMinRes () {} diff --git a/deal.II/lac/include/lac/solver_qmrs.h b/deal.II/lac/include/lac/solver_qmrs.h index f3df068bb8..532bee0dec 100644 --- a/deal.II/lac/include/lac/solver_qmrs.h +++ b/deal.II/lac/include/lac/solver_qmrs.h @@ -112,8 +112,16 @@ class SolverQMRS : public Solver * Constructor. */ SolverQMRS (SolverControl &cn, - VectorMemory &mem, - const AdditionalData &data=AdditionalData()); + VectorMemory &mem, + const AdditionalData &data=AdditionalData()); + + /** + * Constructor. Use an object of + * type PrimitiveVectorMemory as + * a default to allocate memory. + */ + SolverQMRS (SolverControl &cn, + const AdditionalData &data=AdditionalData()); /** * Solve the linear system $Ax=b$ @@ -202,12 +210,24 @@ class SolverQMRS : public Solver template SolverQMRS::SolverQMRS(SolverControl &cn, VectorMemory &mem, - const AdditionalData &data) : + const AdditionalData &data) + : Solver(cn,mem), additional_data(data) {} + +template +SolverQMRS::SolverQMRS(SolverControl &cn, + const AdditionalData &data) + : + Solver(cn), + additional_data(data) +{} + + + template double SolverQMRS::criterion() diff --git a/deal.II/lac/include/lac/solver_richardson.h b/deal.II/lac/include/lac/solver_richardson.h index f14a9456d2..f9ccaa06c2 100644 --- a/deal.II/lac/include/lac/solver_richardson.h +++ b/deal.II/lac/include/lac/solver_richardson.h @@ -84,6 +84,14 @@ class SolverRichardson : public Solver VectorMemory &mem, const AdditionalData &data=AdditionalData()); + /** + * Constructor. Use an object of + * type PrimitiveVectorMemory as + * a default to allocate memory. + */ + SolverRichardson (SolverControl &cn, + const AdditionalData &data=AdditionalData()); + /** * Virtual destructor. */ @@ -188,16 +196,28 @@ AdditionalData (const double omega, use_preconditioned_residual(use_preconditioned_residual) {} + template SolverRichardson::SolverRichardson(SolverControl &cn, VectorMemory &mem, - const AdditionalData &data): + const AdditionalData &data) + : Solver (cn,mem), additional_data(data) {} +template +SolverRichardson::SolverRichardson(SolverControl &cn, + const AdditionalData &data) + : + Solver (cn), + additional_data(data) +{} + + + template SolverRichardson::~SolverRichardson() {} diff --git a/deal.II/lac/include/lac/solver_selector.h b/deal.II/lac/include/lac/solver_selector.h index a24f7e4b6e..470fc8a651 100644 --- a/deal.II/lac/include/lac/solver_selector.h +++ b/deal.II/lac/include/lac/solver_selector.h @@ -95,13 +95,13 @@ class SolverSelector : public Subscriptor public: /** - * Constructor. Takes the @p SolverName, - * the @p SolverControl and the - * @p VectorMemory as argument. + * Constructor. Use the arguments + * to initialize actual solver + * objects. */ SolverSelector (const std::string &solvername, SolverControl &control, - VectorMemory &vectorm); + VectorMemory &vector_memory); /** * Destructor @@ -246,31 +246,31 @@ SolverSelector::solve(const Matrix &A, if (solver_name=="richardson") { SolverRichardson solver(*control,*vector_memory, - richardson_data); + richardson_data); solver.solve(A,x,b,precond); } else if (solver_name=="cg") { SolverCG solver(*control,*vector_memory, - cg_data); + cg_data); solver.solve(A,x,b,precond); } else if (solver_name=="bicgstab") { SolverBicgstab solver(*control,*vector_memory, - bicgstab_data); + bicgstab_data); solver.solve(A,x,b,precond); } else if (solver_name=="gmres") { SolverGMRES solver(*control,*vector_memory, - gmres_data); + gmres_data); solver.solve(A,x,b,precond); } else if (solver_name=="fgmres") { SolverFGMRES solver(*control,*vector_memory, - fgmres_data); + fgmres_data); solver.solve(A,x,b,precond); } else diff --git a/deal.II/lac/source/solver.cc b/deal.II/lac/source/solver.cc new file mode 100644 index 0000000000..7159c8cece --- /dev/null +++ b/deal.II/lac/source/solver.cc @@ -0,0 +1,38 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2005 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- + +#include +#include +#include +#include +#include +#include + + +// a few instantiations of static members. Hope that we catch all that +// are required + +template +PrimitiveVectorMemory +Solver::static_vector_memory; + +template class Solver >; +template class Solver >; + +template class Solver >; +template class Solver >; + +#ifdef DEAL_II_USE_PETSC +template class Solver; +template class Solver; +#endif