<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>
#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>
#include <lac/matrix_lib.h>
#include <lac/full_matrix.h>
#include <lac/vector.h>
-#include <lac/vector_memory.h>
double Adata[] =
#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>
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);
#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>
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,
#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>
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);
#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>
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);
#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>
{
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);
#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>
// 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
// 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,
#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>
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));
{
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),
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());
#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
// 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
#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>
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());
#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>
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
#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>
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);
#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>
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);
#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>
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);
#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>
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);
#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 */
{
public:
/**
- * Constructor. Assign a control
+ * Constructor. Takes a control
* object which evaluates the
* conditions for convergence,
* and an object to provide
* 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
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.
*/
/*@}*/
/*-------------------------------- 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 &
}
-template<class VECTOR>
-inline
-Solver<VECTOR>::Solver(SolverControl &cn, VectorMemory<VECTOR> &mem)
- : cntrl(cn),
- memory(mem)
-{}
#endif
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.
*/
+template<class VECTOR>
+SolverBicgstab<VECTOR>::SolverBicgstab (SolverControl &cn,
+ const AdditionalData &data)
+ :
+ Solver<VECTOR>(cn),
+ additional_data(data)
+{}
+
+
+
template<class VECTOR>
SolverBicgstab<VECTOR>::~SolverBicgstab ()
{}
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.
*/
{}
+
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)
+template <class VECTOR>
+SolverCG<VECTOR>::SolverCG (SolverControl &cn,
+ const AdditionalData &data)
+ :
+ Solver<VECTOR>(cn),
+ additional_data(data)
+{}
+
+
+
template <class VECTOR>
SolverCG<VECTOR>::~SolverCG ()
{}
*/
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());
/**
*/
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());
/**
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
}
+
template<class VECTOR>
template<class MATRIX, class PRECONDITIONER>
void
}
+
template<class VECTOR>
double
SolverGMRES<VECTOR>::criterion ()
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
*/
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());
/**
{}
+
+template<class VECTOR>
+SolverMinRes<VECTOR>::SolverMinRes (SolverControl &cn,
+ const AdditionalData &)
+ :
+ Solver<VECTOR>(cn)
+{}
+
+
template<class VECTOR>
SolverMinRes<VECTOR>::~SolverMinRes ()
{}
* 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$
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()
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.
*/
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()
{}
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
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
--- /dev/null
+//---------------------------------------------------------------------------
+// $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