]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add additional constructors to the solver classes where one does not have to give...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 28 Nov 2005 17:11:58 +0000 (17:11 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 28 Nov 2005 17:11:58 +0000 (17:11 +0000)
git-svn-id: https://svn.dealii.org/trunk@11798 0785d39b-7218-0410-832d-ea1e28bc413d

26 files changed:
deal.II/doc/news/changes.html
deal.II/examples/doxygen/block_matrix_array.cc
deal.II/examples/doxygen/product_matrix.cc
deal.II/examples/step-11/step-11.cc
deal.II/examples/step-12/step-12.cc
deal.II/examples/step-13/step-13.cc
deal.II/examples/step-14/step-14.cc
deal.II/examples/step-15/step-15.cc
deal.II/examples/step-16/step-16.cc
deal.II/examples/step-20/step-20.cc
deal.II/examples/step-3/step-3.cc
deal.II/examples/step-4/step-4.cc
deal.II/examples/step-5/step-5.cc
deal.II/examples/step-6/step-6.cc
deal.II/examples/step-7/step-7.cc
deal.II/examples/step-8/step-8.cc
deal.II/examples/step-9/step-9.cc
deal.II/lac/include/lac/solver.h
deal.II/lac/include/lac/solver_bicgstab.h
deal.II/lac/include/lac/solver_cg.h
deal.II/lac/include/lac/solver_gmres.h
deal.II/lac/include/lac/solver_minres.h
deal.II/lac/include/lac/solver_qmrs.h
deal.II/lac/include/lac/solver_richardson.h
deal.II/lac/include/lac/solver_selector.h
deal.II/lac/source/solver.cc [new file with mode: 0644]

index cb26090d199b1e8babe453fcf51735230a456128..5d5f7b13f0dd55176fd78434fbe79ebd0fddea93 100644 (file)
@@ -177,6 +177,24 @@ inconvenience this causes.
 <h3>lac</h3>
 
 <ol>
+  <li> <p>
+       New: All solver classes took an argument of type <code
+       class="class">VectorMemory</code> 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 <code class="class">PrimitiveVectorMemory</code>
+       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.
+       <br> 
+       (WB 2005/11/23)
+       </p>
+
   <li> <p>
        Fixed: The <code class="class">SparseMatrix</code> iterators had a
        problem when one wrote code like <code>iterator->value()=0</code>
index 03dc2ea18982aec6d207e7d2f55c37e957e8e649..4a0cec66d8e1baf29d2c2e1522bd80650957fe27 100644 (file)
@@ -18,7 +18,6 @@
 #include <lac/full_matrix.h>
 #include <lac/vector.h>
 #include <lac/block_vector.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <lac/solver_cg.h>
 #include <lac/solver_gmres.h>
index b6ffaa8da687952da0e341b0830a7e3c739c6a0d..7e92d864bd76e9e71504e7f6ef27f3f18078836d 100644 (file)
@@ -17,7 +17,6 @@
 #include <lac/matrix_lib.h>
 #include <lac/full_matrix.h>
 #include <lac/vector.h>
-#include <lac/vector_memory.h>
 
 
 double Adata[] =
index 02f1388810a4dfae31ab8fefbbd501136add14da..fb0ae9afcd0833fddaf11f7f196968361f84c7ce 100644 (file)
@@ -22,7 +22,6 @@
 #include <lac/vector.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
@@ -583,8 +582,7 @@ template <int dim>
 void LaplaceProblem<dim>::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);
index 761bf736c57cabbe7a0e2a33d12fcd51bb540106..14657c8fa465a374764d0a7476b9af85d9b856db 100644 (file)
@@ -19,7 +19,6 @@
 #include <base/function.h>
 #include <lac/vector.h>
 #include <lac/sparse_matrix.h>
-#include <lac/vector_memory.h>
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
 #include <grid/grid_out.h>
@@ -1431,8 +1430,7 @@ template <int dim>
 void DGMethod<dim>::solve (Vector<double> &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,
index 85a3c8d0df1c25173d157b28ddc63c42e1351de6..365415c99ae21a6fcb911613c30c2cab1dabd979 100644 (file)
@@ -30,7 +30,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
@@ -1513,8 +1512,7 @@ namespace LaplaceSolver
   Solver<dim>::LinearSystem::solve (Vector<double> &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);
index a06c14a1794b0f872e6d0f4656ec594e51957818..6126068763942c16dc91c79099d315dde8221320 100644 (file)
@@ -21,7 +21,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
@@ -761,8 +760,7 @@ namespace LaplaceSolver
   Solver<dim>::LinearSystem::solve (Vector<double> &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);
index 11ba91ff552844e80ed2949373230dd0657478cd..31b4c040820db88b4487b25d5d42b9a406a18be7 100644 (file)
@@ -21,7 +21,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
@@ -804,8 +803,7 @@ void MinimizationProblem<dim>::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);
index 06ca843a165f3700d76b01a01362ccd25058502c..dde5cf2cea2e4343fc52f69b4b04db54a8900305 100644 (file)
@@ -21,7 +21,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <grid/tria_accessor.h>
@@ -393,8 +392,8 @@ void LaplaceProblem<dim>::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<dim>::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,
index e23f8d0bd72e977f0dd80bc04d64005044e41259..be5230e3b4ac5454e6e7d8dd04938eea746ae047 100644 (file)
@@ -40,7 +40,6 @@ const unsigned int degree = 2;
 #include <lac/block_sparse_matrix.h>
 #include <lac/solver_cg.h>
 #include <lac/solver_gmres.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 
 #include <numerics/data_out.h>
@@ -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<dim>::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<dim>::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());
index aa83fc01dd995d2b81eff8a3ba6ae890df15c136..071dcc93508dc7c074827c06a6f51d524613f436 100644 (file)
@@ -79,7 +79,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 
                                 // 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
index 9e1be30a087c68ac4c944335e72c01bd3dd0a0fc..9e068f90250c6db0afcae11cb2a07afc960615d9 100644 (file)
@@ -33,7 +33,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 
 #include <numerics/data_out.h>
@@ -497,8 +496,7 @@ template <int dim>
 void LaplaceProblem<dim>::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());
 
index 93d2697133c48032f25b33f12e731fa5651a789c..749bdcad840a5da7153136af8161063b7a280037 100644 (file)
@@ -21,7 +21,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <dofs/dof_handler.h>
@@ -536,8 +535,7 @@ template <int dim>
 void LaplaceProblem<dim>::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
index dda25221d7e876ee04575762cc0ed555c4b40be6..d70b1495b4698da6feece633e03707e29de471a7 100644 (file)
@@ -22,7 +22,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <dofs/dof_handler.h>
@@ -596,8 +595,7 @@ template <int dim>
 void LaplaceProblem<dim>::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);
index 01f571c7b60a418dceaa4f9a3de93a6ff27db4fd..175a77ac7d131d27c1608b8382bdc88d5e2c3ec0 100644 (file)
@@ -22,7 +22,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
@@ -1057,8 +1056,7 @@ template <int dim>
 void LaplaceProblem<dim>::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);
index 3ae01d2ef1e90dd0d5a611afb4bfa4f732d2e577..cdd271daae13ace965f9f4fc5700dbd9250af5bc 100644 (file)
@@ -21,7 +21,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_cg.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
@@ -761,8 +760,7 @@ template <int dim>
 void ElasticProblem<dim>::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);
index 6b9e4c1c613bda37e4a8c105f8763ab8f624b33b..173d0fcc097ef87b6991901b67c8e2753ab113bd 100644 (file)
@@ -22,7 +22,6 @@
 #include <lac/full_matrix.h>
 #include <lac/sparse_matrix.h>
 #include <lac/solver_bicgstab.h>
-#include <lac/vector_memory.h>
 #include <lac/precondition.h>
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
@@ -1330,8 +1329,7 @@ template <int dim>
 void AdvectionProblem<dim>::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);
index c4977b6d1730071f86f8bd78ef07f9e0d2dc4b00..9a9e00e54e0a3fceae88074f3d8ed3f516b356f0 100644 (file)
 
 #include <base/config.h>
 #include <base/subscriptor.h>
+#include <lac/vector_memory.h>
+
 
 template <typename number> class Vector;
-template <class VECTOR>    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<VECTOR> &);
+    Solver (SolverControl        &solver_control,
+            VectorMemory<VECTOR> &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<VECTOR> static_vector_memory;
+    
                                     /**
                                      * Control structure.
                                      */
@@ -183,6 +211,27 @@ class Solver : public Subscriptor
 /*@}*/
 /*-------------------------------- Inline functions ------------------------*/
 
+template<class VECTOR>
+inline
+Solver<VECTOR>::Solver (SolverControl        &solver_control,
+                       VectorMemory<VECTOR> &vector_memory)
+               :
+               cntrl(solver_control),
+               memory(vector_memory)
+{}
+
+
+
+template<class VECTOR>
+inline
+Solver<VECTOR>::Solver (SolverControl        &solver_control)
+               :
+               cntrl(solver_control),
+               memory(static_vector_memory)
+{}
+
+
+
 template <class VECTOR>
 inline
 SolverControl &
@@ -192,12 +241,6 @@ Solver<VECTOR>::control() const
 }
 
 
-template<class VECTOR>
-inline
-Solver<VECTOR>::Solver(SolverControl &cn, VectorMemory<VECTOR> &mem)
-               : cntrl(cn),
-                 memory(mem)
-{}
 
 
 #endif
index 98627bc0559bf631fc7f830ed9da41c1a551c748..4198d90c3e7b004aeda0136c8e03a867bb7f949b 100644 (file)
@@ -102,6 +102,14 @@ class SolverBicgstab : public Solver<VECTOR>
                    VectorMemory<VECTOR> &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<VECTOR>::SolverBicgstab (SolverControl &cn,
 
 
 
+template<class VECTOR>
+SolverBicgstab<VECTOR>::SolverBicgstab (SolverControl &cn,
+                                       const AdditionalData &data)
+               :
+               Solver<VECTOR>(cn),
+               additional_data(data)
+{}
+
+
+
 template<class VECTOR>
 SolverBicgstab<VECTOR>::~SolverBicgstab ()
 {}
index 671f678d0d4ad2c6012fe54bdbf617e0e1621693..171731815b58e38737a1563946a85ec2a823975a 100644 (file)
@@ -98,6 +98,14 @@ class SolverCG : public Solver<VECTOR>
              VectorMemory<VECTOR> &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 <class VECTOR>
-SolverCG<VECTOR>::SolverCG(SolverControl        &cn,
-                          VectorMemory<VECTOR> &mem,
-                          const AdditionalData &data)
+SolverCG<VECTOR>::SolverCG (SolverControl        &cn,
+                           VectorMemory<VECTOR> &mem,
+                           const AdditionalData &data)
                :
                 Solver<VECTOR>(cn,mem),
                 additional_data(data)
@@ -195,6 +204,16 @@ SolverCG<VECTOR>::SolverCG(SolverControl        &cn,
 
 
 
+template <class VECTOR>
+SolverCG<VECTOR>::SolverCG (SolverControl        &cn,
+                           const AdditionalData &data)
+               :
+                Solver<VECTOR>(cn),
+                additional_data(data)
+{}
+
+
+
 template <class VECTOR>
 SolverCG<VECTOR>::~SolverCG ()
 {}
index ffe284b13c4f1354aa66d6f614a3c7221a6adfb8..bd6d5ef3eca7b36c85af5b21bcc46c726195697d 100644 (file)
@@ -231,6 +231,14 @@ class SolverGMRES : public Solver<VECTOR>
                                      */
     SolverGMRES (SolverControl        &cn,
                 VectorMemory<VECTOR> &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<VECTOR>
                                      */
     SolverFGMRES (SolverControl        &cn,
                  VectorMemory<VECTOR> &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 <class VECTOR>
 SolverGMRES<VECTOR>::SolverGMRES (SolverControl        &cn,
                                   VectorMemory<VECTOR> &mem,
-                                  const AdditionalData &data) :
+                                  const AdditionalData &data)
+               :
                Solver<VECTOR> (cn,mem),
                additional_data(data)
 {}
 
 
+
+template <class VECTOR>
+SolverGMRES<VECTOR>::SolverGMRES (SolverControl        &cn,
+                                  const AdditionalData &data) :
+               Solver<VECTOR> (cn),
+               additional_data(data)
+{}
+
+
+
 template <class VECTOR>
 inline
 void
@@ -484,6 +512,7 @@ SolverGMRES<VECTOR>::givens_rotation (Vector<double> &h,
 }
 
 
+
 template<class VECTOR>
 template<class MATRIX, class PRECONDITIONER>
 void
@@ -795,6 +824,7 @@ SolverGMRES<VECTOR>::solve (const MATRIX         &A,
 }
 
 
+
 template<class VECTOR>
 double
 SolverGMRES<VECTOR>::criterion () 
@@ -815,10 +845,21 @@ SolverFGMRES<VECTOR>::SolverFGMRES (SolverControl        &cn,
                                    const AdditionalData &data)
                :
                Solver<VECTOR> (cn, mem),
-  additional_data(data)
+               additional_data(data)
 {}
 
 
+
+template <class VECTOR>
+SolverFGMRES<VECTOR>::SolverFGMRES (SolverControl        &cn,
+                                   const AdditionalData &data)
+               :
+               Solver<VECTOR> (cn),
+               additional_data(data)
+{}
+
+
+
 template<class VECTOR>
 template<class MATRIX, class PRECONDITIONER>
 void
index 1385fbe1121a64cdc00f5bea5385aedc715d5275..1dc07b43cd61f77d7e9c3b56a82de3d7e9059615 100644 (file)
@@ -73,6 +73,14 @@ class SolverMinRes : public Solver<VECTOR>
                                      */
     SolverMinRes (SolverControl &cn,
                  VectorMemory<VECTOR> &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<VECTOR>::SolverMinRes (SolverControl &cn,
 {}
 
 
+
+template<class VECTOR>
+SolverMinRes<VECTOR>::SolverMinRes (SolverControl &cn,
+                                   const AdditionalData &)
+               :
+               Solver<VECTOR>(cn)
+{}
+
+
 template<class VECTOR>
 SolverMinRes<VECTOR>::~SolverMinRes ()
 {}
index f3df068bb8b87e03678a00c8b68cd0c2f8a6ccd9..532bee0dec6bfc04a10588510ef705ce43bb7ded 100644 (file)
@@ -112,8 +112,16 @@ class SolverQMRS : public Solver<VECTOR>
                                      * Constructor.
                                      */
     SolverQMRS (SolverControl &cn,
-             VectorMemory<VECTOR> &mem,
-             const AdditionalData &data=AdditionalData());
+               VectorMemory<VECTOR> &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<VECTOR>
 template<class VECTOR>
 SolverQMRS<VECTOR>::SolverQMRS(SolverControl &cn,
                               VectorMemory<VECTOR> &mem,
-                              const AdditionalData &data) :
+                              const AdditionalData &data)
+               :
                Solver<VECTOR>(cn,mem),
                additional_data(data)
 {}
 
 
+
+template<class VECTOR>
+SolverQMRS<VECTOR>::SolverQMRS(SolverControl &cn,
+                              const AdditionalData &data)
+               :
+               Solver<VECTOR>(cn),
+               additional_data(data)
+{}
+
+
+
 template<class VECTOR>
 double
 SolverQMRS<VECTOR>::criterion()
index f14a9456d230522bef5cba61038562ec2497d8d6..f9ccaa06c2d2334dcc46bd039bf45d2446ba99b8 100644 (file)
@@ -84,6 +84,14 @@ class SolverRichardson : public Solver<VECTOR>
                      VectorMemory<VECTOR> &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 <class VECTOR>
 SolverRichardson<VECTOR>::SolverRichardson(SolverControl &cn,
                                           VectorMemory<VECTOR> &mem,
-                                          const AdditionalData &data):
+                                          const AdditionalData &data)
+               :
                Solver<VECTOR> (cn,mem),
                additional_data(data)
 {}
 
 
 
+template <class VECTOR>
+SolverRichardson<VECTOR>::SolverRichardson(SolverControl &cn,
+                                          const AdditionalData &data)
+               :
+               Solver<VECTOR> (cn),
+               additional_data(data)
+{}
+
+
+
 template <class VECTOR>
 SolverRichardson<VECTOR>::~SolverRichardson()
 {}
index a24f7e4b6e97e58c4be8b876eeaf452207d75606..470fc8a651650861999da1982b6ca2f242b3802f 100644 (file)
@@ -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<VECTOR> &vectorm);
+                   VectorMemory<VECTOR> &vector_memory);
   
                                     /**
                                      * Destructor
@@ -246,31 +246,31 @@ SolverSelector<VECTOR>::solve(const Matrix &A,
   if (solver_name=="richardson")
     {
       SolverRichardson<VECTOR> solver(*control,*vector_memory,
-                                            richardson_data);
+                                     richardson_data);
       solver.solve(A,x,b,precond);
     }       
   else if (solver_name=="cg")
     {
       SolverCG<VECTOR> solver(*control,*vector_memory,
-                                    cg_data);
+                             cg_data);
       solver.solve(A,x,b,precond);
     }
   else if (solver_name=="bicgstab")
     {
       SolverBicgstab<VECTOR> solver(*control,*vector_memory,
-                                          bicgstab_data);
+                                   bicgstab_data);
       solver.solve(A,x,b,precond);
     }
   else if (solver_name=="gmres")
     {
       SolverGMRES<VECTOR> solver(*control,*vector_memory,
-                                       gmres_data);
+                                gmres_data);
       solver.solve(A,x,b,precond);
     }
   else if (solver_name=="fgmres")
     {
       SolverFGMRES<VECTOR> 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 (file)
index 0000000..7159c8c
--- /dev/null
@@ -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 <lac/solver.h>
+#include <lac/vector_memory.h>
+#include <lac/vector.h>
+#include <lac/block_vector.h>
+#include <lac/petsc_vector.h>
+#include <lac/petsc_block_vector.h>
+
+
+// a few instantiations of static members. Hope that we catch all that
+// are required
+
+template <class VECTOR>
+PrimitiveVectorMemory<VECTOR>
+Solver<VECTOR>::static_vector_memory;
+
+template class Solver<Vector<double> >;
+template class Solver<BlockVector<double> >;
+
+template class Solver<Vector<float> >;
+template class Solver<BlockVector<float> >;
+
+#ifdef DEAL_II_USE_PETSC
+template class Solver<PETScWrappers::Vector>;
+template class Solver<PETScWrappers::BlockVector>;
+#endif

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.