#include <numerics/base.h>
#include <numerics/assembler.h>
#include <numerics/vectors.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <map>
#include <fstream>
Equation<dim>(1),
right_hand_side (rhs) {};
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
protected:
template <>
-void PoissonEquation<2>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<2>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<2> &fe_values,
const DoFHandler<2>::cell_iterator &) const {
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void PoissonEquation<dim>::assemble (dVector &,
+void PoissonEquation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
solve ();
Solution<dim> sol;
- dVector l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell;
- dVector h1_seminorm_error_per_cell, h1_error_per_cell;
+ Vector<float> l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell;
+ Vector<float> h1_seminorm_error_per_cell, h1_error_per_cell;
cout << " Calculating L1 error... ";
VectorTools<dim>::integrate_difference (*dof_handler,
if (dof->n_dofs()<=5000)
{
- dVector l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof;
- dVector h1_seminorm_error_per_dof, h1_error_per_dof;
+ Vector<double> l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof;
+ Vector<double> h1_seminorm_error_per_dof, h1_error_per_dof;
dof->distribute_cell_to_dof_vector (l1_error_per_cell, l1_error_per_dof);
dof->distribute_cell_to_dof_vector (l2_error_per_cell, l2_error_per_dof);
dof->distribute_cell_to_dof_vector (linfty_error_per_cell,
h1_seminorm_error_per_dof);
dof->distribute_cell_to_dof_vector (h1_error_per_cell, h1_error_per_dof);
-// dVector projected_solution;
+// Vector<double> projected_solution;
// ConstraintMatrix constraints;
// constraints.close ();
// VectorTools<dim>::project (*dof, constraints, *fe,
#include <grid/tria_boundary.h>
#include <grid/tria_iterator.h>
#include <grid/tria_accessor.h>
-#include <lac/dsmatrix.h>
+#include <lac/sparsematrix.h>
#include <base/parameter_handler.h>
#include <grid/dof_constraints.h>
cout << " Renumbering degrees of freedom..." << endl;
dof->renumber_dofs (Cuthill_McKee, false);
- dSMatrixStruct sparsity (dof->n_dofs(),
- dof->max_couplings_between_dofs());
+ SparseMatrixStruct sparsity (dof->n_dofs(),
+ dof->max_couplings_between_dofs());
dof->make_sparsity_pattern (sparsity);
#include <numerics/assembler.h>
#include <numerics/vectors.h>
#include <numerics/error_estimator.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <map.h>
#include <fstream.h>
right_hand_side (rhs),
coefficient (coefficient) {};
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
protected:
template <>
-void PoissonEquation<2>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<2>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<2> &fe_values,
const DoFHandler<2>::cell_iterator &) const {
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void PoissonEquation<dim>::assemble (dVector &,
+void PoissonEquation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
cout << " Solving..." << endl;
solve ();
- dVector l2_error_per_cell, linfty_error_per_cell, h1_error_per_cell;
- dVector estimated_error_per_cell;
+ Vector<float> l2_error_per_cell, linfty_error_per_cell, h1_error_per_cell;
+ Vector<float> estimated_error_per_cell;
QGauss3<dim> q;
cout << " Calculating L2 error... ";
cout << estimated_error_per_cell.l2_norm() << endl;
estimated_error.push_back (estimated_error_per_cell.l2_norm());
- dVector l2_error_per_dof, linfty_error_per_dof;
- dVector h1_error_per_dof, estimated_error_per_dof;
- dVector error_ratio;
+ Vector<double> l2_error_per_dof, linfty_error_per_dof;
+ Vector<double> h1_error_per_dof, estimated_error_per_dof;
+ Vector<double> error_ratio;
dof->distribute_cell_to_dof_vector (l2_error_per_cell, l2_error_per_dof);
dof->distribute_cell_to_dof_vector (linfty_error_per_cell,
linfty_error_per_dof);
#include <numerics/vectors.h>
#include <numerics/matrices.h>
-#include <lac/dvector.h>
-#include <lac/dsmatrix.h>
+#include <lac/vector.h>
+#include <lac/sparsematrix.h>
#include <lac/solver_cg.h>
#include <lac/vector_memory.h>
Equation<dim>(1),
right_hand_side (rhs) {};
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
protected:
/**
* Sparsity pattern of the system matrix.
*/
- dSMatrixStruct system_sparsity;
+ SparseMatrixStruct system_sparsity;
/**
* System matrix.
*/
- dSMatrix system_matrix;
+ SparseMatrix<double> system_matrix;
/**
* Vector storing the right hand side.
*/
- dVector right_hand_side;
+ Vector<double> right_hand_side;
/**
* Solution vector.
*/
- dVector solution;
+ Vector<double> solution;
/**
* List of constraints introduced by
template <>
-void PoissonEquation<2>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<2>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<2> &fe_values,
const DoFHandler<2>::cell_iterator &) const {
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void PoissonEquation<dim>::assemble (dVector &,
+void PoissonEquation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
Assert ((tria!=0) && (dof!=0), ExcNoTriaSelected());
SolverControl control(4000, 1e-16);
- PrimitiveVectorMemory<dVector> memory;
- SolverCG<dSMatrix,dVector> cg(control,memory);
+ PrimitiveVectorMemory<Vector<double> > memory;
+ SolverCG<SparseMatrix<double>,Vector<double> > cg(control,memory);
// solve
cg.solve (system_matrix, solution, right_hand_side);
solve ();
Solution<dim> sol;
- dVector l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell;
- dVector h1_seminorm_error_per_cell, h1_error_per_cell;
+ Vector<float> l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell;
+ Vector<float> h1_seminorm_error_per_cell, h1_error_per_cell;
cout << " Calculating L1 error... ";
VectorTools<dim>::integrate_difference (*dof,
if (dof->DoFHandler<dim>::n_dofs()<=5000)
{
- dVector l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof;
- dVector h1_seminorm_error_per_dof, h1_error_per_dof;
+ Vector<double> l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof;
+ Vector<double> h1_seminorm_error_per_dof, h1_error_per_dof;
dof->distribute_cell_to_dof_vector (l1_error_per_cell, l1_error_per_dof);
dof->distribute_cell_to_dof_vector (l2_error_per_cell, l2_error_per_dof);
dof->distribute_cell_to_dof_vector (linfty_error_per_cell,
h1_seminorm_error_per_dof);
dof->distribute_cell_to_dof_vector (h1_error_per_cell, h1_error_per_dof);
-// dVector projected_solution;
+// Vector<double> projected_solution;
// ConstraintMatrix constraints;
// constraints.close ();
// VectorTools<dim>::project (*dof, constraints, *fe,
#include "poisson.h"
-#include <lac/dvector.h>
+#include <lac/vector.h>
#if deal_II_dimension == 1
template <>
-void PoissonEquation<1>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<1>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<1> &fe_values,
const DoFHandler<1>::cell_iterator &) const {
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
//#if deal_II_dimension >= 2
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &) const {
const vector<vector<Tensor<1,dim> > >&gradients = fe_values.get_shape_grads ();
- const dFMatrix &values = fe_values.get_shape_values ();
+ const FullMatrix<double> &values = fe_values.get_shape_values ();
vector<double> rhs_values (fe_values.n_quadrature_points);
const vector<double> &weights = fe_values.get_JxW_values ();
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void PoissonEquation<dim>::assemble (dVector &,
+void PoissonEquation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
#include "poisson.h"
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <base/quadrature_lib.h>
#include <numerics/base.h>
#include <numerics/assembler.h>
-#include <lac/dsmatrix.h>
+#include <lac/sparsematrix.h>
#include <map.h>
#include <fstream.h>
#include <cmath>
-extern "C" {
-# include <stdlib.h>
-}
+#include <cstdlib>
+
Equation<dim>(1),
right_hand_side (rhs) {};
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
protected:
#include "poisson.h"
-#include <lac/dvector.h>
+#include <lac/vector.h>
/*---------------------------- dof_iterator.h ---------------------------*/
-template <int dim> class Triangulation;
-template <int dim> class DoFHandler;
-class dVector;
-class dFMatrix;
-class dSMatrix;
-
-
#include <grid/tria_accessor.h>
+#include <lac/forward-declarations.h>
#include <vector>
* rather than just set, since this is
* usually what one wants.
*/
- void distribute_local_to_global (const dVector &local_source,
- dVector &global_destination) const;
+ void distribute_local_to_global (const Vector<double> &local_source,
+ Vector<double> &global_destination) const;
/**
* This function does much the same as the
* is supposed to have non-zero entry
* slots where required.
*/
- void distribute_local_to_global (const dFMatrix &local_source,
- dSMatrix &global_destination) const;
+ void distribute_local_to_global (const FullMatrix<double> &local_source,
+ SparseMatrix<double> &global_destination) const;
/**
* Implement the copy operator needed
* rather than just set, since this is
* usually what one wants.
*/
- void distribute_local_to_global (const dVector &local_source,
- dVector &global_destination) const;
+ void distribute_local_to_global (const Vector<double> &local_source,
+ Vector<double> &global_destination) const;
/**
* This function does much the same as the
* is supposed to have non-zero entry
* slots where required.
*/
- void distribute_local_to_global (const dFMatrix &local_source,
- dSMatrix &global_destination) const;
+ void distribute_local_to_global (const FullMatrix<double> &local_source,
+ SparseMatrix<double> &global_destination) const;
/**
* Implement the copy operator needed
* rather than just set, since this is
* usually what one wants.
*/
- void distribute_local_to_global (const dVector &local_source,
- dVector &global_destination) const;
+ void distribute_local_to_global (const Vector<double> &local_source,
+ Vector<double> &global_destination) const;
/**
* This function does much the same as the
* is supposed to have non-zero entry
* slots where required.
*/
- void distribute_local_to_global (const dFMatrix &local_source,
- dSMatrix &global_destination) const;
+ void distribute_local_to_global (const FullMatrix<double> &local_source,
+ SparseMatrix<double> &global_destination) const;
/**
* Implement the copy operator needed
* function is only callable for active
* cells.
*/
- void get_dof_values (const dVector &values,
- dVector &local_values) const;
+ void get_dof_values (const Vector<double> &values,
+ Vector<double> &local_values) const;
/**
* Return the interpolation of the given
* what the this function does in these
* cases.
*/
- void get_interpolated_dof_values (const dVector &values,
- dVector &interpolated_values) const;
+ void get_interpolated_dof_values (const Vector<double> &values,
+ Vector<double> &interpolated_values) const;
/**
* This function is the counterpart to
* It is assumed that both vectors already
* have the right size beforehand.
*/
- void set_dof_values (const dVector &local_values,
- dVector &values) const;
+ void set_dof_values (const Vector<double> &local_values,
+ Vector<double> &values) const;
/**
* This, again, is the counterpart to
* the prolongation matrices represent in
* this case.
*/
- void set_dof_values_by_interpolation (const dVector &local_values,
- dVector &values) const;
+ void set_dof_values_by_interpolation (const Vector<double> &local_values,
+ Vector<double> &values) const;
/**
* Exception
#include <vector>
#include <utility>
#include <base/exceptions.h>
+#include <lac/forward-declarations.h>
-class ostream;
-class dSMatrix;
-class dSMatrixStruct;
-class dVector;
-
-
/**
* This class represents the matrix denoting the distribution of the degrees
* of freedom of hanging nodes.
* The constraint matrix object must be
* closed to call this function.
*/
- void condense (const dSMatrixStruct &uncondensed,
- dSMatrixStruct &condensed) const;
+ void condense (const SparseMatrixStruct &uncondensed,
+ SparseMatrixStruct &condensed) const;
/**
* The matrix struct is compressed at the
* end of the function.
*/
- void condense (dSMatrixStruct &sparsity) const;
+ void condense (SparseMatrixStruct &sparsity) const;
/**
* Condense a given matrix. The associated
* The constraint matrix object must be
* closed to call this function.
*/
- void condense (const dSMatrix &uncondensed,
- dSMatrix &condensed) const;
+ void condense (const SparseMatrix<double> &uncondensed,
+ SparseMatrix<double> &condensed) const;
/**
* This function does much the same as
* documentation of this class for more
* detailed information.
*/
- void condense (dSMatrix &matrix) const;
+ void condense (SparseMatrix<double> &matrix) const;
/**
* Condense the given vector #uncondensed#
* responsibility to guarantee that all
* entries of #condensed# be zero!
*/
- void condense (const dVector &uncondensed,
- dVector &condensed) const;
+ void condense (const Vector<double> &uncondensed,
+ Vector<double> &condensed) const;
/**
* Condense the given vector in-place.
*/
- void condense (dVector &vec) const;
+ void condense (Vector<double> &vec) const;
/**
* Re-distribute the elements of the vector
* #condense# somehow, but it should be noted
* that it is not the inverse of #condense#-
*/
- void distribute (const dVector &condensed,
- dVector &uncondensed) const;
+ void distribute (const Vector<double> &condensed,
+ Vector<double> &uncondensed) const;
/**
* Re-distribute the elements of the vector
* in-place.
*/
- void distribute (dVector &vec) const;
+ void distribute (Vector<double> &vec) const;
/**
#include <vector>
#include <map>
#include <base/exceptions.h>
+#include <base/forward-declarations.h>
#include <base/smartpointer.h>
-
-template <int dim> class TriaAccessor;
-template <int dim> class LineAccessor;
-template <int dim> class QuadAccessor;
-template <int dim> class HexAccessor;
-template <int dim> class CellAccessor;
-
-template <int dim, class BaseClass> class DoFLineAccessor;
-template <int dim, class BaseClass> class DoFQuadAccessor;
-template <int dim, class BaseClass> class DoFHexAccessor;
-template <int dim> class DoFCellAccessor;
-
-template <int dim> class DoFLevel;
-
-
-template <int dim, class Accessor> class TriaRawIterator;
-template <int dim, class Accessor> class TriaIterator;
-template <int dim, class Accessor> class TriaActiveIterator;
-
-template <int dim> class Triangulation;
-template <int dim> class FiniteElement;
-template <int dim> class Function;
-
-class dVector;
-class dSMatrix;
-class dSMatrixStruct;
-class ConstraintMatrix;
-
-
-
+#include <basic/forward-declarations.h>
+#include <lac/forward-declarations.h>
*
* The renumbering algorithms need quite a lot of memory, since they have
* to store for each dof with which other dofs it couples. This is done
- * using a #dSMatrixStruct# object used to store the sparsity pattern of
+ * using a #SparseMatrixStruct# object used to store the sparsity pattern of
* matrices. It
* is not useful for the user to do anything between distributing the dofs
* and renumbering, i.e. the calls to #DoFHandler::distribute_dofs# and
*
* To build the matrix, you first have to call
* #make_transfer_matrix (old_dof_object, sparsity_pattern);#, then create a
- * sparse matrix out of this pattern, e.g. by #dSMatrix m(sparsity_pattern);#
+ * sparse matrix out of this pattern, e.g. by #SparseMatrix<double> m(sparsity_pattern);#
* and finally give this to the second run:
* #make_transfer_matrix (old_dof_object, m);#. The spasity pattern created
* by the first run is automatically compressed.
*
- * When creating the #dSMatrixStruct# sparsity pattern, you have to give the
+ * When creating the #SparseMatrixStruct# sparsity pattern, you have to give the
* dimension and the maximum number of entries per row. Obviously the image
* dimension is the number of dofs on the new grid (you can get this using the
* #n_dofs()# function), while the range dimension is the number of dofs on the
* added. However, if you don't want to call
* #ConstraintMatrix::condense(1)#, you
* have to compress the matrix yourself,
- * using #dSMatrixStruct::compress()#.
+ * using #SparseMatrixStruct::compress()#.
*/
- void make_sparsity_pattern (dSMatrixStruct &) const;
+ void make_sparsity_pattern (SparseMatrixStruct &) const;
/**
* Write the sparsity structure of the
* added. However, if you don't want to call
* #ConstraintMatrix::condense(1)#, you
* have to compress the matrix yourself,
- * using #dSMatrixStruct::compress()#.
+ * using #SparseMatrixStruct::compress()#.
*
* Since this function is obviously useless
* in one spatial dimension, it is not
* implemented.
*/
void make_boundary_sparsity_pattern (const vector<int> &dof_to_boundary_mapping,
- dSMatrixStruct &) const;
+ SparseMatrixStruct &) const;
/**
* Write the sparsity structure of the
* in one spatial dimension, it is not
* implemented.
*/
- void make_boundary_sparsity_pattern (const FunctionMap &boundary_indicators,
- const vector<int> &dof_to_boundary_mapping,
- dSMatrixStruct &sparsity) const;
+ void make_boundary_sparsity_pattern (const FunctionMap &boundary_indicators,
+ const vector<int> &dof_to_boundary_mapping,
+ SparseMatrixStruct &sparsity) const;
/**
* documentation for this class.
*/
void make_transfer_matrix (const DoFHandler<dim> &transfer_from,
- dSMatrixStruct &transfer_pattern) const;
+ SparseMatrixStruct &transfer_pattern) const;
/**
* Make up the transfer matrix which
* documentation for this class.
*/
void make_transfer_matrix (const DoFHandler<dim> &transfer_from,
- dSMatrix &transfer_matrix) const;
+ SparseMatrix<double> &transfer_matrix) const;
/**
* Return the maximum number of
* This is the maximum number of entries
* per line in the system matrix; this
* information can therefore be used upon
- * construction of the #dSMatrixStruct#
+ * construction of the #SparseMatrixStruct#
* object.
*
* The returned number is not really the
* number of active cells. The size of
* #dof_data# is adjusted to the right
* size.
+ *
+ * Note that the input vector may be
+ * a vector of any data type as long
+ * as it is convertible to #double#.
+ * The output vector, being a data
+ * vector on the grid, always consists
+ * of elements of type #double#.
*/
- void distribute_cell_to_dof_vector (const dVector &cell_data,
- dVector &dof_data) const;
+ template <typename Number>
+ void distribute_cell_to_dof_vector (const Vector<Number> &cell_data,
+ Vector<double> &dof_data) const;
/**
* Create a mapping from degree of freedom
*/
void transfer_cell (const cell_iterator &old_cell,
const cell_iterator &new_cell,
- dSMatrixStruct &transfer_pattern) const;
+ SparseMatrixStruct &transfer_pattern) const;
/**
* Make up part of the transfer matrix by
* looking at the two cells given.
*/
- void transfer_cell (const cell_iterator &old_cell,
- const cell_iterator &new_cell,
- dSMatrix &transfer_matrix) const;
+ void transfer_cell (const cell_iterator &old_cell,
+ const cell_iterator &new_cell,
+ SparseMatrix<double> &transfer_matrix) const;
/**
* Space to store the DoF numbers for the
* other types than those explicitely instantiated.
*/
template <int N>
-class DoFLevel;
+class DoFLevel {
+ private:
+ /**
+ * Make the constructor private to avoid
+ * that someone uses this class.
+ */
+ DoFLevel ();
+};
#include <grid/dof_accessor.h>
-// forward declaration
-template <int dim> class MGDoFHandler;
-
-
/**
* entries as there are degrees of
* freedom on this level.
*/
- void get_mg_dof_values (const dVector &values,
- dVector &dof_values) const;
+ void get_mg_dof_values (const Vector<double> &values,
+ Vector<double> &dof_values) const;
/**
-// forward declarations
-template <int dim, typename BaseClass> class MGDoFLineAccessor;
-template <int dim, typename BaseClass> class MGDoFQuadAccessor;
-template <int dim> class MGDoFCellAccessor;
-
-
-
-
/**
* Define some types which differ between the dimensions. This class
* is analogous to the \Ref{TriaDimensionInfo} class hierarchy.
* added. However, if you don't want to call
* #ConstraintMatrix::condense(1)#, you
* have to compress the matrix yourself,
- * using #dSMatrixStruct::compress()#.
+ * using #SparseMatrixStruct::compress()#.
*/
- void make_sparsity_pattern (const unsigned int level,
- dSMatrixStruct &sparsity) const;
+ void make_sparsity_pattern (const unsigned int level,
+ SparseMatrixStruct &sparsity) const;
/*--------------------------------------*/
#include <base/tensor.h>
#include <grid/dof.h>
#include <grid/geometry_info.h>
-#include <lac/dfmatrix.h>
+#include <lac/fullmatrix.h>
-template <int dim> class Boundary;
-template <int dim> class FiniteElementData;
-
-
/**
* Dimension independent data for finite elements. See the #FiniteElementBase#
* mother cell. See the #restriction# array
* for more information.
*/
- const dFMatrix & restrict (const unsigned int child) const;
+ const FullMatrix<double> & restrict (const unsigned int child) const;
/**
* Return a readonly reference to the
* mother cell to the child with the given
* number.
*/
- const dFMatrix & prolongate (const unsigned int child) const;
+ const FullMatrix<double> & prolongate (const unsigned int child) const;
/**
* Return a readonly reference to the
* one space dimension, since there are no
* constraints then.
*/
- const dFMatrix & constraints () const;
+ const FullMatrix<double> & constraints () const;
/**
* Comparison operator. We also check for
* matrix are discarded and will not fill
* up the transfer matrix.
*/
- dFMatrix restriction[GeometryInfo<dim>::children_per_cell];
+ FullMatrix<double> restriction[GeometryInfo<dim>::children_per_cell];
/**
* Have #N=2^dim# matrices keeping the
* matrix are discarded and will not fill
* up the transfer matrix.
*/
- dFMatrix prolongation[GeometryInfo<dim>::children_per_cell];
+ FullMatrix<double> prolongation[GeometryInfo<dim>::children_per_cell];
/**
* Specify the constraints which the
* This field is obviously useless in one
* space dimension.
*/
- dFMatrix interface_constraints;
+ FullMatrix<double> interface_constraints;
/**
* Map between linear dofs and component dofs.
const bool compute_support_points,
vector<Point<dim> > &q_points,
const bool compute_q_points,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_grads_transform,
const Boundary<dim> &boundary) const;
const bool compute_face_jacobians,
vector<Point<dim> > &normal_vectors,
const bool compute_normal_vectors,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_grads_transform,
const Boundary<dim> &boundary) const;
const bool compute_face_jacobians,
vector<Point<dim> > &normal_vectors,
const bool compute_normal_vectors,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_grads_transform,
const Boundary<dim> &boundary) const;
* documentation.
*/
virtual void get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
- const Boundary<dim> &boundary,
- dFMatrix &local_mass_matrix) const =0;
+ const Boundary<dim> &boundary,
+ FullMatrix<double> &local_mass_matrix) const =0;
/**
* Exception
#include <fe/fe.h>
#include <base/quadrature.h>
-#include <math.h>
+#include <cmath>
+
/**
* This class implements a rather unusual macro element, the so-called
* triangular element, which uses a linear mapping. The missing linearity
* makes assemblage of matrices a bit more complicated, since the gradient
* is not constant and we need more than one quadrature point, as well
- * as some other dubtle difficulties. This problem can, however, be cured
+ * as some other subtle difficulties. This problem can, however, be cured
* using a trick: the usual transformation from unit coordinates $\vec\xi$
* to real coordinates $\vec x(\vec\xi)$ looks like
* $$
*/
virtual void get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
const Boundary<dim> &boundary,
- dFMatrix &local_mass_matrix) const;
+ FullMatrix<double> &local_mass_matrix) const;
/**
* Return the value of the #i#th shape
const bool compute_support_points,
vector<Point<dim> > &q_points,
const bool compute_q_points,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_grad_transform,
const Boundary<dim> &boundary) const;
*/
virtual void get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
const Boundary<dim> &boundary,
- dFMatrix &local_mass_matrix) const;
+ FullMatrix<double> &local_mass_matrix) const;
/**
* Return a readonly reference to the
* matrices is not yet finally decided
* about.
*/
- const dFMatrix & restrict (const unsigned int) const;
+ const FullMatrix<double> & restrict (const unsigned int) const;
/**
* Exception
* matrices is not yet finally decided
* about.
*/
- const dFMatrix & restrict (const unsigned int) const;
+ const FullMatrix<double> & restrict (const unsigned int) const;
};
* matrices is not yet finally decided
* about.
*/
- const dFMatrix & restrict (const unsigned int) const;
+ const FullMatrix<double> & restrict (const unsigned int) const;
};
* matrices is not yet finally decided
* about.
*/
- const dFMatrix & restrict (const unsigned int) const;
+ const FullMatrix<double> & restrict (const unsigned int) const;
};
* matrices is not yet finally decided
* about.
*/
- const dFMatrix & restrict (const unsigned int) const;
+ const FullMatrix<double> & restrict (const unsigned int) const;
};
*/
virtual void get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
const Boundary<dim> &boundary,
- dFMatrix &local_mass_matrix) const;
+ FullMatrix<double> &local_mass_matrix) const;
private:
*/
virtual void get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
const Boundary<dim> &boundary,
- dFMatrix &local_mass_matrix) const;
+ FullMatrix<double> &local_mass_matrix) const;
private:
*/
virtual void get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
const Boundary<dim> &boundary,
- dFMatrix &local_mass_matrix) const;
+ FullMatrix<double> &local_mass_matrix) const;
private:
*/
virtual void get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
const Boundary<dim> &boundary,
- dFMatrix &local_mass_matrix) const;
+ FullMatrix<double> &local_mass_matrix) const;
private:
/*---------------------------- fe_linear_mapping.h ---------------------------*/
-#include <math.h>
+#include <cmath>
#include <fe/fe.h>
const bool compute_support_points,
vector<Point<dim> > &q_points,
const bool compute_q_points,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_grad_transform,
const Boundary<dim> &boundary) const;
#include <vector>
#include <pair.h>
+
+
/**
* This class provides an interface to group several equal elements together
* into one. To the outside world, the resulting object looks just like
*/
virtual void get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
const Boundary<dim> &boundary,
- dFMatrix &local_mass_matrix) const;
+ FullMatrix<double> &local_mass_matrix) const;
/**
* Return the value of the #i#th shape
const bool compute_support_points,
vector<Point<dim> > &q_points,
const bool compute_q_points,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_grad_transform,
const Boundary<dim> &boundary) const;
#include <base/exceptions.h>
#include <base/subscriptor.h>
-#include <lac/dfmatrix.h>
+#include <lac/fullmatrix.h>
#include <base/point.h>
#include <grid/tria.h>
#include <fe/fe_update_flags.h>
-// forward declarations
-template <int dim> class Boundary;
-template <int dim> class FiniteElement;
-template <int dim> class Quadrature;
-
/**
* For the format of this matrix, see the
* documentation for the matrix itself.
*/
- const dFMatrix & get_shape_values () const;
+ const FullMatrix<double> & get_shape_values () const;
/**
* Return the values of the finite
* #values# object already has the
* right size.
*/
- void get_function_values (const dVector &fe_function,
- vector<double> &values) const;
+ void get_function_values (const Vector<double> &fe_function,
+ vector<double> &values) const;
/**
* Access to vector valued finite element functions.
* but applied to multi-component
* elements.
*/
- void get_function_values (const dVector &fe_function,
- vector<vector<double> > &values) const;
+ void get_function_values (const Vector<double> &fe_function,
+ vector<vector<double> > &values) const;
/**
* Return the gradient of the #i#th shape
* #gradients# object already has the
* right size.
*/
- void get_function_grads (const dVector &fe_function,
+ void get_function_grads (const Vector<double> &fe_function,
vector<Tensor<1,dim> > &gradients) const;
/**
* #second_derivatives# object already has
* the right size.
*/
- void get_function_2nd_derivatives (const dVector &fe_function,
+ void get_function_2nd_derivatives (const Vector<double> &fe_function,
vector<Tensor<2,dim> > &second_derivatives) const;
/**
* is determined by the #selected_dataset#
* variable.
*/
- vector<dFMatrix> shape_values;
+ vector<FullMatrix<double> > shape_values;
/**
* Store the gradients of the shape
* of the #i#th transfer basis function
* at the #j#th quadrature point.
*/
- vector<dFMatrix> shape_values_transform;
+ vector<FullMatrix<double> > shape_values_transform;
/**
* Store which of the data sets in the
template <int dim>
inline
-const dFMatrix & FEValuesBase<dim>::get_shape_values () const {
+const FullMatrix<double> & FEValuesBase<dim>::get_shape_values () const {
Assert (selected_dataset<shape_values.size(),
ExcInvalidIndex (selected_dataset, shape_values.size()));
return shape_values[selected_dataset];
/*---------------------------- tria.h ---------------------------*/
#include <vector>
+#include <basic/forward-declarations.h>
+#include <lac/forward-declarations.h>
+#include <base/forward-declarations.h>
#include <base/point.h>
#include <grid/geometry_info.h>
#include <base/subscriptor.h>
-//forward declaration needed
-template <int dim> class Boundary;
-
-template <int dim> class TriaAccessor;
-template <int dim> class LineAccessor;
-template <int dim> class QuadAccessor;
-template <int dim> class HexAccessor;
-template <int dim> class CellAccessor;
-
-template <int dim> class TriangulationLevel;
-
-template <int dim, class Accessor> class TriaRawIterator;
-template <int dim, class Accessor> class TriaIterator;
-template <int dim, class Accessor> class TriaActiveIterator;
-
-template <int dim> class DoFHandler;
-template <int dim> class MGDoFHandler;
-
-class dVector;
-
-
-class istream;
-class ostream;
-
/**
* elements.
*
* The central function to this is
- * #refine (const dVector &criterion, const double threshold)#: it takes a
+ * #refine (const Vector<float> &criterion, const double threshold)#: it takes a
* vector of values, one per active cell, which denote the criterion according
* to which the triangulation is to be refined. It marks all cells for which
* the criterion is greater than the threshold being given as the second
* argument. Analogously,
- * #coarsen (const dVector &criterion, const double threshold)# flags those
+ * #coarsen (const Vector<float> &criterion, const double threshold)# flags those
* cells for coarsening for which the criterion is less than the treshold.
*
* There are two variations of these functions, which rely on #refine# and
* for mesh refinement; refer to the
* following functions and to the general
* doc for this class for more information.
+ *
+ * Note that this function takes a vector
+ * of #float#s, rather than the usual
+ * #double#s, since accuracy is not so
+ * much needed here and saving memory may
+ * be a good goal when using many cells.
*/
- void refine (const dVector &criteria,
- const double threshold);
+ void refine (const Vector<float> &criteria,
+ const double threshold);
/**
* Analogue to the #refine# function:
* which the absolute value of the
* criterion is less than the
* given threshold.
+ *
+ * Note that this function takes a vector
+ * of #float#s, rather than the usual
+ * #double#s, since accuracy is not so
+ * much needed here and saving memory may
+ * be a good goal when using many cells.
*/
- void coarsen (const dVector &criteria,
- const double threshold);
+ void coarsen (const Vector<float> &criteria,
+ const double threshold);
/**
* Refine the triangulation by refining
*
* Refer to the general doc of this class
* for more information.
+ *
+ * Note that this function takes a vector
+ * of #float#s, rather than the usual
+ * #double#s, since accuracy is not so
+ * much needed here and saving memory may
+ * be a good goal when using many cells.
*/
- void refine_and_coarsen_fixed_number (const dVector &criteria,
- const double top_fraction_of_cells,
- const double bottom_fraction_of_cells);
+ void refine_and_coarsen_fixed_number (const Vector<float> &criteria,
+ const double top_fraction_of_cells,
+ const double bottom_fraction_of_cells);
/**
* Refine the triangulation by flagging
*
* Refer to the general doc of this class
* for more information.
+ *
+ * Note that this function takes a vector
+ * of #float#s, rather than the usual
+ * #double#s, since accuracy is not so
+ * much needed here and saving memory may
+ * be a good goal when using many cells.
*/
- void refine_and_coarsen_fixed_fraction (const dVector &criteria,
- const double top_fraction,
- const double bottom_fraction);
+ void refine_and_coarsen_fixed_fraction (const Vector<float> &criteria,
+ const double top_fraction,
+ const double bottom_fraction);
/**
* Refine all cells on all levels which
#include <base/exceptions.h>
-
-
-// forward declaration needed
-class Line;
-class Quad;
-class Hexahedron;
-
-template <int dim> class Point;
-
-
-template <int dim> class TriaAccessor;
-template <int dim> class LineAccessor;
-template <int dim> class QuadAccessor;
-template <int dim> class CellAccessor;
-
-template <int dim> class DoFCellAccessor;
-template <int dim, class BaseClass> class DoFLineAccessor;
-
-template <int dim, class Accessor> class TriaRawIterator;
-template <int dim, class Accessor> class TriaIterator;
-template <int dim, class Accessor> class TriaActiveIterator;
-
-template <int dim> class Triangulation;
-
+#include <base/forward-declarations.h>
+#include <basic/forward-declarations.h>
#include <base/point.h>
#include <base/subscriptor.h>
#include <grid/geometry_info.h>
-
-
-template <int dim> class Triangulation;
+#include <basic/forward-declarations.h>
-// forward declaration needed
-template <int dim> class Triangulation;
-
-
-
-
-
-
/**
#include <base/exceptions.h>
-template <int dim> class Triangulation;
/**
#include <grid/dof_accessor.h>
-// forward declaration
-template <int dim> class MGDoFHandler;
-
-
/**
* entries as there are degrees of
* freedom on this level.
*/
- void get_mg_dof_values (const dVector &values,
- dVector &dof_values) const;
+ void get_mg_dof_values (const Vector<double> &values,
+ Vector<double> &dof_values) const;
/**
-// forward declarations
-template <int dim, typename BaseClass> class MGDoFLineAccessor;
-template <int dim, typename BaseClass> class MGDoFQuadAccessor;
-template <int dim> class MGDoFCellAccessor;
-
-
-
-
/**
* Define some types which differ between the dimensions. This class
* is analogous to the \Ref{TriaDimensionInfo} class hierarchy.
* added. However, if you don't want to call
* #ConstraintMatrix::condense(1)#, you
* have to compress the matrix yourself,
- * using #dSMatrixStruct::compress()#.
+ * using #SparseMatrixStruct::compress()#.
*/
- void make_sparsity_pattern (const unsigned int level,
- dSMatrixStruct &sparsity) const;
+ void make_sparsity_pattern (const unsigned int level,
+ SparseMatrixStruct &sparsity) const;
/*--------------------------------------*/
#define __multigrid_H
/*---------------------------- multigrid.h ---------------------------*/
-class dVector;
+
+#include <lac/forward-declarations.h>
+
+
/**
* Vector with data for each level.
*/
class MGVector
{
- MGVector(const MGVector&);
-public:
- /**
- * Constructor using the number of
- * levels.
- */
- MGVector(unsigned l);
- /**
- * Access to data on a specific
- level.
- */
- dVector& operator[](unsigned l);
-
- /**
- * Access to data on a specific
- level.
- */
- const dVector& operator[](unsigned l) const;
-
+ MGVector(const MGVector&);
+ public:
+ /**
+ * Constructor using the number of
+ * levels.
+ */
+ MGVector(unsigned l);
+ /**
+ * Access to data on a specific
+ * level.
+ */
+ Vector<double>& operator[](unsigned l);
+
+ /**
+ * Access to data on a specific
+ level.
+ */
+ const Vector<double>& operator[](unsigned l) const;
};
* corresponding levels of an
* MGVector.
*/
- void copy_to_mg(MGVector& dst, const dVector& src) const;
+ void copy_to_mg(MGVector& dst, const Vector<double>& src) const;
/**
* Transfer from MGVector to
- * dVector.
+ * Vector<double>.
*
* Copies data from active portions
* of an MGVector into the
* respective positions of a
- * dVector. All other entries of
+ * Vector<double>. All other entries of
* #src# are zero.
*/
- void copy_from_mg(dVector& dst, const MGVector& src) const;
+ void copy_from_mg(Vector<double>& dst, const MGVector& src) const;
/**
* The actual v-cycle multigrid method.
*
*/
virtual void smooth(unsigned level,
- dVector& x,
- const dVector& b,
+ Vector<double>& x,
+ const Vector<double>& b,
unsigned steps) const;
/**
* Defaults to #smooth#.
*/
virtual void post_smooth(unsigned level,
- dVector& dst,
- const dVector& src,
+ Vector<double>& dst,
+ const Vector<double>& src,
unsigned steps) const;
/**
*
*/
virtual void level_vmult(unsigned level,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
/**
* Apply operator on non-refined
* cells.
* normal fine-grid matrix.
*/
virtual void level_active_vmult(unsigned level,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
/**
* Restriction from #level#.
* and #dst# is on #level#-1.
*/
virtual void restriction(unsigned level,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
/**
* #level# and #src# on #level#-1.
*/
virtual void prolongation(unsigned level,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
/**
* Solve exactly on coarsest grid.
*/
virtual void coarse_grid_solution(unsigned l,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
* the outer iteration.
*
*/
- void vmult(dVector& dst, const dVector& src) const;
+ void vmult(Vector<double>& dst, const Vector<double>& src) const;
/** Multigrid preconditioning.
*
* This is the function, where the
* multigrid method is implemented.
*
*/
- void precondition(dVector& dst, const dVector& src) const;
+ void precondition(Vector<double>& dst, const Vector<double>& src) const;
};
#include <vector>
-// forward declarations
-template <int dim> class FiniteElement;
-template <int dim> class Quadrature;
-
-class dFMatrix;
-class dVector;
-
* #ExcWrongSize# and #ExcObjectNotEmpty#
* are declared.
*/
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
* #ExcWrongSize# and #ExcObjectNotEmpty#
* are declared.
*/
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
* #ExcWrongSize# and #ExcObjectNotEmpty#
* are declared.
*/
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
AssemblerData (const DoFHandler<dim> &dof,
const bool assemble_matrix,
const bool assemble_rhs,
- dSMatrix &matrix,
- dVector &rhs_vector,
+ SparseMatrix<double> &matrix,
+ Vector<double> &rhs_vector,
const Quadrature<dim> &quadrature,
const FiniteElement<dim> &fe,
- const UpdateFlags &update_flags,
+ const UpdateFlags &update_flags,
const Boundary<dim> &boundary);
/**
* to clear this object (set all entries
* to zero) before use.
*/
- dSMatrix &matrix;
+ SparseMatrix<double> &matrix;
/**
* Pointer to the vector to be assembled
* to clear this object (set all entries
* to zero) before use.
*/
- dVector &rhs_vector;
+ Vector<double> &rhs_vector;
/**
* Pointer to a quadrature object to be
/**
* Store a local cell matrix.
*/
- dFMatrix cell_matrix;
+ FullMatrix<double> cell_matrix;
/**
* Right hand side local to cell.
*/
- dVector cell_vector;
+ Vector<double> cell_vector;
/**
* Store whether to assemble the
* Pointer to the matrix to be assembled
* by this object.
*/
- dSMatrix &matrix;
+ SparseMatrix<double> &matrix;
/**
* Pointer to the vector to be assembled
* by this object.
*/
- dVector &rhs_vector;
+ Vector<double> &rhs_vector;
/**
* Pointer to the finite element used for
/*---------------------------- problem_base.h ---------------------------*/
-#include <lac/dsmatrix.h>
+#include <lac/sparsematrix.h>
#include <base/exceptions.h>
#include <grid/dof_constraints.h>
#include <fe/fe_update_flags.h>
#include <map>
-// forward declaration
-template <int dim> class Triangulation;
-template <int dim> class DoFHandler;
-template <int dim> class FiniteElement;
-template <int dim> class Quadrature;
-template <int dim> class DataOut;
-template <int dim> class Equation;
-template <int dim> class Assembler;
-template <int dim> class Boundary;
-template <int dim> class StraightBoundary;
-template <int dim> class Function;
-
/**
* Sparsity pattern of the system matrix.
*/
- dSMatrixStruct system_sparsity;
+ SparseMatrixStruct system_sparsity;
/**
* System matrix.
*/
- dSMatrix system_matrix;
+ SparseMatrix<double> system_matrix;
/**
* Vector storing the right hand side.
*/
- dVector right_hand_side;
+ Vector<double> right_hand_side;
/**
* Solution vector.
*/
- dVector solution;
+ Vector<double> solution;
/**
* List of constraints introduced by
/*---------------------------- data_io.h ---------------------------*/
#include <base/exceptions.h>
+#include <basic/forward-declarations.h>
+#include <lac/forward-declarations.h>
#include <vector>
#include <string>
-template <int dim> class Triangulation;
-template <int dim> class DoFHandler;
-
-class dVector;
-
-
* see which characters are valid and which
* are not.
*/
- void add_data_vector (const dVector &data,
- const string &name,
- const string &units="<dimensionless>");
+ void add_data_vector (const Vector<double> &data,
+ const string &name,
+ const string &units="<dimensionless>");
/**
* Release the pointers to the data
/**
* Constructor
*/
- DataEntry (const dVector *data, const string name, const string units);
+ DataEntry (const Vector<double> *data, const string name, const string units);
/**
* Pointer to the data vector.
*/
- const dVector *data;
+ const Vector<double> *data;
+
/**
* Name of this component.
*/
string name;
+
/**
* Physical unit name of this
* component.
#include <map>
-// forward declarations
-template <int dim> class DoFHandler;
-template <int dim> class Quadrature;
-template <int dim> class FiniteElement;
-template <int dim> class FEFaceValues;
-template <int dim> class FESubfaceValues;
-template <int dim> class Boundary;
-template <int dim> class Function;
-class dVector;
-
-
/**
* the conormal derivative $a\frac{du}{dn} = g$.
*
* The error estimator returns a vector of estimated errors per cell which
- * can be used to feed the #Triangulation<dim>::refine_*# functions.
+ * can be used to feed the #Triangulation<dim>::refine_*# functions. This
+ * vector contains elements of data type #float#, rather than #double#,
+ * since accuracy is not so important here, and since this can save rather
+ * a lot of memory, when using many cells.
*
*
* \subsection{Implementation}
const FiniteElement<dim> &fe,
const Boundary<dim> &boundary,
const FunctionMap &neumann_bc,
- const dVector &solution,
- dVector &error,
+ const Vector<double> &solution,
+ Vector<float> &error,
const Function<dim> *coefficient = 0);
/**
FEFaceValues<dim> &fe_face_values_cell,
FEFaceValues<dim> &fe_face_values_neighbor,
FaceIntegrals &face_integrals,
- const dVector &solution,
+ const Vector<double>&solution,
const Function<dim> *coefficient);
/**
FEFaceValues<dim> &fe_face_values,
FESubfaceValues<dim> &fe_subface_values,
FaceIntegrals &face_integrals,
- const dVector &solution,
+ const Vector<double> &solution,
const Function<dim> *coefficient);
};
#include <base/exceptions.h>
#include <map>
-template <int dim> class Triangulation;
-template <int dim> class DoFHandler;
-template <int dim> class FiniteElement;
-template <int dim> class FEValues;
-template <int dim> class Quadrature;
-template <int dim> class Function;
-template <int dim> class Boundary;
-template <int dim> class Equation;
-
-class dVector;
-class dFMatrix;
-class dSMatrix;
*
* All functions take a sparse matrix object to hold the matrix to be
* created. The functions assume that the matrix is initialized with a
- * sparsity pattern (#dSMatrixStruct#) corresponding to the given degree
+ * sparsity pattern (#SparseMatrixStruct#) corresponding to the given degree
* of freedom handler, i.e. the sparsity structure is already as needed.
* You can do this by calling the #DoFHandler<dim>::make_sparsity_pattern#
* function.
static void create_mass_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const Function<dim> *a = 0);
/**
static void create_mass_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const Function<dim> &rhs,
- dVector &rhs_vector,
+ Vector<double> &rhs_vector,
const Function<dim> *a = 0);
/**
*/
static void create_mass_matrix (const DoFHandler<dim> &dof,
const Boundary<dim> &boundary,
- dSMatrix &matrix);
+ SparseMatrix<double> &matrix);
/**
* Assemble the mass matrix and a right
static void create_boundary_mass_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim-1> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const FunctionMap &rhs,
- dVector &rhs_vector,
+ Vector<double> &rhs_vector,
vector<int> &vec_to_dof_mapping,
const Function<dim> *a = 0);
static void create_laplace_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const Function<dim> *a = 0);
/**
static void create_laplace_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const Function<dim> &rhs,
- dVector &rhs_vector,
+ Vector<double> &rhs_vector,
const Function<dim> *a = 0);
/**
*/
static void create_interpolation_matrix(const FiniteElement<dim> &high,
const FiniteElement<dim> &low,
- dFMatrix& result);
+ FullMatrix<double>& result);
/**
* documentation.
*/
static void apply_boundary_values (const map<int,double> &boundary_values,
- dSMatrix &matrix,
- dVector &solution,
- dVector &right_hand_side);
+ SparseMatrix<double> &matrix,
+ Vector<double> &solution,
+ Vector<double> &right_hand_side);
/**
* Exception
* and right hand side to have the right
* size and to be empty.
*/
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const typename DoFHandler<dim>::cell_iterator &) const;
* If a coefficient was given to the
* constructor, it is used.
*/
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const typename DoFHandler<dim>::cell_iterator &) const;
* constructor in order to call this
* function.
*/
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const typename DoFHandler<dim>::cell_iterator &) const;
* a coefficient was given to the
* constructor, it is used.
*/
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const typename DoFHandler<dim>::cell_iterator &) const;
* If a coefficient was given to the
* constructor, it is used.
*/
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const typename DoFHandler<dim>::cell_iterator &) const;
* constructor in order to call this
* function.
*/
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const typename DoFHandler<dim>::cell_iterator &) const;
#define __multigrid_H
/*---------------------------- multigrid.h ---------------------------*/
-class dVector;
+
+#include <lac/forward-declarations.h>
+
+
/**
* Vector with data for each level.
*/
class MGVector
{
- MGVector(const MGVector&);
-public:
- /**
- * Constructor using the number of
- * levels.
- */
- MGVector(unsigned l);
- /**
- * Access to data on a specific
- level.
- */
- dVector& operator[](unsigned l);
-
- /**
- * Access to data on a specific
- level.
- */
- const dVector& operator[](unsigned l) const;
-
+ MGVector(const MGVector&);
+ public:
+ /**
+ * Constructor using the number of
+ * levels.
+ */
+ MGVector(unsigned l);
+ /**
+ * Access to data on a specific
+ * level.
+ */
+ Vector<double>& operator[](unsigned l);
+
+ /**
+ * Access to data on a specific
+ level.
+ */
+ const Vector<double>& operator[](unsigned l) const;
};
* corresponding levels of an
* MGVector.
*/
- void copy_to_mg(MGVector& dst, const dVector& src) const;
+ void copy_to_mg(MGVector& dst, const Vector<double>& src) const;
/**
* Transfer from MGVector to
- * dVector.
+ * Vector<double>.
*
* Copies data from active portions
* of an MGVector into the
* respective positions of a
- * dVector. All other entries of
+ * Vector<double>. All other entries of
* #src# are zero.
*/
- void copy_from_mg(dVector& dst, const MGVector& src) const;
+ void copy_from_mg(Vector<double>& dst, const MGVector& src) const;
/**
* The actual v-cycle multigrid method.
*
*/
virtual void smooth(unsigned level,
- dVector& x,
- const dVector& b,
+ Vector<double>& x,
+ const Vector<double>& b,
unsigned steps) const;
/**
* Defaults to #smooth#.
*/
virtual void post_smooth(unsigned level,
- dVector& dst,
- const dVector& src,
+ Vector<double>& dst,
+ const Vector<double>& src,
unsigned steps) const;
/**
*
*/
virtual void level_vmult(unsigned level,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
/**
* Apply operator on non-refined
* cells.
* normal fine-grid matrix.
*/
virtual void level_active_vmult(unsigned level,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
/**
* Restriction from #level#.
* and #dst# is on #level#-1.
*/
virtual void restriction(unsigned level,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
/**
* #level# and #src# on #level#-1.
*/
virtual void prolongation(unsigned level,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
/**
* Solve exactly on coarsest grid.
*/
virtual void coarse_grid_solution(unsigned l,
- dVector& dst,
- const dVector& src) const;
+ Vector<double>& dst,
+ const Vector<double>& src) const;
* the outer iteration.
*
*/
- void vmult(dVector& dst, const dVector& src) const;
+ void vmult(Vector<double>& dst, const Vector<double>& src) const;
/** Multigrid preconditioning.
*
* This is the function, where the
* multigrid method is implemented.
*
*/
- void precondition(dVector& dst, const dVector& src) const;
+ void precondition(Vector<double>& dst, const Vector<double>& src) const;
};
/*---------------------------- vectors.h ---------------------------*/
-
+#include <base/forward-declarations.h>
#include <base/exceptions.h>
#include <map>
-template <int dim> class DoFHandler;
-template <int dim> class Function;
-template <int dim> class Quadrature;
-template <int dim> class QGauss2;
-template <int dim> class FiniteElement;
-template <int dim> class Boundary;
-template <int dim> class StraightBoundary;
-class ConstraintMatrix;
-class dVector;
-class VectorFunction;
/**
* full $H_1$ norm is the sum of the seminorm and the $L_2$ norm.
*
* To get the {\it global} L_1 error, you have to sum up the entries in
- * #difference#, e.g. using #dVector::l1_norm# function.
+ * #difference#, e.g. using #Vector<double>::l1_norm# function.
* For the global L_2 difference, you have to sum up the squares of the
- * entries and take the root of the sum, e.g. using #dVector::l2_norm#.
+ * entries and take the root of the sum, e.g. using #Vector<double>::l2_norm#.
* These two operations represent the
* l_1 and l_2 norms of the vectors, but you need not take the absolute
* value of each entry, since the cellwise norms are already positive.
*
* To get the global mean difference, simply sum up the elements as above.
* To get the L_\infty norm, take the maximum of the vector elements, e.g.
- * using the #dVector::linfty_norm# function.
+ * using the #Vector<double>::linfty_norm# function.
*
* For the global $H_1$ norm and seminorm, the same rule applies as for the
* $L_2$ norm: compute the $l_2$ norm of the cell error vector.
static void interpolate (const DoFHandler<dim> &dof,
const Boundary<dim> &boundary,
const Function<dim> &function,
- dVector &vec);
+ Vector<double> &vec);
/**
* Interpolate different finite
*/
static void interpolate(const DoFHandler<dim> &high_dof,
const DoFHandler<dim> &low_dof,
- const dFMatrix &transfer,
- const dVector &high,
- dVector &low);
+ const FullMatrix<double> &transfer,
+ const Vector<double> &high,
+ Vector<double> &low);
/**
* Compute the projection of
const Boundary<dim> &boundary,
const Quadrature<dim> &q,
const Function<dim> &function,
- dVector &vec,
+ Vector<double> &vec,
const bool enforce_zero_boundary = false,
const Quadrature<dim-1> &q_boundary = QGauss2<dim-1>(),
const bool project_to_boundary_first = false);
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
const Function<dim> &rhs,
- dVector &rhs_vector);
+ Vector<double> &rhs_vector);
/**
* Make up the list of node subject
* is given as a continuous function
* object.
*
+ * Note that this function returns
+ * its results in a vector of #float#s,
+ * rather than in a vector of #double#s,
+ * since accuracy is not that important
+ * here and to save memory. During
+ * computation of the results, the full
+ * accuracy of the #double# data type is
+ * used.
+ *
* See the general documentation of this
* class for more information.
*/
static void integrate_difference (const DoFHandler<dim> &dof,
- const dVector &fe_function,
+ const Vector<double> &fe_function,
const Function<dim> &exact_solution,
- dVector &difference,
+ Vector<float> &difference,
const Quadrature<dim> &q,
const NormType &norm,
- const Boundary<dim> &boundary=StraightBoundary<dim>());
+ const Boundary<dim> &boundary=StraightBoundary<dim>());
/**
* Compute the error for the solution of a system.
* See the other #integrate_difference#.
*/
- static void integrate_difference (const DoFHandler<dim> &dof,
- const dVector &fe_function,
- const VectorFunction &exact_solution,
- dVector &difference,
- const Quadrature<dim> &q,
- const FiniteElement<dim> &fe,
- const NormType &norm,
- const Boundary<dim> &boundary);
+// static void integrate_difference (const DoFHandler<dim> &dof,
+// const Vector<double> &fe_function,
+// const TensorFunction<1,dim>&exact_solution,
+// Vector<float> &difference,
+// const Quadrature<dim> &q,
+// const FiniteElement<dim> &fe,
+// const NormType &norm,
+// const Boundary<dim> &boundary);
/**
#include <grid/tria_iterator.templates.h>
#include <fe/fe.h>
-#include <lac/dvector.h>
-#include <lac/dfmatrix.h>
-#include <lac/dsmatrix.h>
+#include <lac/vector.h>
+#include <lac/fullmatrix.h>
+#include <lac/sparsematrix.h>
#include <vector>
template <int dim, typename BaseClass>
void DoFLineAccessor<dim,BaseClass>::
-distribute_local_to_global (const dVector &local_source,
- dVector &global_destination) const {
+distribute_local_to_global (const Vector<double> &local_source,
+ Vector<double> &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
Assert (local_source.size() == (2*dof_handler->get_fe().dofs_per_vertex +
template <int dim, typename BaseClass>
void DoFLineAccessor<dim,BaseClass>::
-distribute_local_to_global (const dFMatrix &local_source,
- dSMatrix &global_destination) const {
+distribute_local_to_global (const FullMatrix<double> &local_source,
+ SparseMatrix<double> &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
Assert (local_source.m() == (2*dof_handler->get_fe().dofs_per_vertex +
template <int dim, typename BaseClass>
void DoFQuadAccessor<dim,BaseClass>::
-distribute_local_to_global (const dVector &local_source,
- dVector &global_destination) const {
+distribute_local_to_global (const Vector<double> &local_source,
+ Vector<double> &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
Assert (local_source.size() == (4*dof_handler->get_fe().dofs_per_vertex +
template <int dim, typename BaseClass>
void DoFQuadAccessor<dim,BaseClass>::
-distribute_local_to_global (const dFMatrix &local_source,
- dSMatrix &global_destination) const {
+distribute_local_to_global (const FullMatrix<double> &local_source,
+ SparseMatrix<double> &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
Assert (local_source.m() == (4*dof_handler->get_fe().dofs_per_vertex +
template <int dim, typename BaseClass>
void DoFHexAccessor<dim,BaseClass>::
-distribute_local_to_global (const dVector &local_source,
- dVector &global_destination) const {
+distribute_local_to_global (const Vector<double> &local_source,
+ Vector<double> &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
Assert (local_source.size() == (8*dof_handler->get_fe().dofs_per_vertex +
template <int dim, typename BaseClass>
void DoFHexAccessor<dim,BaseClass>::
-distribute_local_to_global (const dFMatrix &local_source,
- dSMatrix &global_destination) const {
+distribute_local_to_global (const FullMatrix<double> &local_source,
+ SparseMatrix<double> &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
Assert (local_source.m() == (8*dof_handler->get_fe().dofs_per_vertex +
template <>
void
-DoFCellAccessor<1>::get_dof_values (const dVector &values,
- dVector &local_values) const {
+DoFCellAccessor<1>::get_dof_values (const Vector<double> &values,
+ Vector<double> &local_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (local_values.size() == dof_handler->get_fe().total_dofs,
template <>
void
-DoFCellAccessor<1>::set_dof_values (const dVector &local_values,
- dVector &values) const {
+DoFCellAccessor<1>::set_dof_values (const Vector<double> &local_values,
+ Vector<double> &values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (local_values.size() == dof_handler->get_fe().total_dofs,
template <>
void
-DoFCellAccessor<2>::get_dof_values (const dVector &values,
- dVector &local_values) const {
+DoFCellAccessor<2>::get_dof_values (const Vector<double> &values,
+ Vector<double> &local_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (local_values.size() == dof_handler->get_fe().total_dofs,
template <>
void
-DoFCellAccessor<2>::set_dof_values (const dVector &local_values,
- dVector &values) const {
+DoFCellAccessor<2>::set_dof_values (const Vector<double> &local_values,
+ Vector<double> &values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (local_values.size() == dof_handler->get_fe().total_dofs,
template <>
void
-DoFCellAccessor<3>::get_dof_values (const dVector &values,
- dVector &local_values) const {
+DoFCellAccessor<3>::get_dof_values (const Vector<double> &values,
+ Vector<double> &local_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (local_values.size() == dof_handler->get_fe().total_dofs,
template <>
void
-DoFCellAccessor<3>::set_dof_values (const dVector &local_values,
- dVector &values) const {
+DoFCellAccessor<3>::set_dof_values (const Vector<double> &local_values,
+ Vector<double> &values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (local_values.size() == dof_handler->get_fe().total_dofs,
template <int dim>
void
-DoFCellAccessor<dim>::get_interpolated_dof_values (const dVector &values,
- dVector &interpolated_values) const {
+DoFCellAccessor<dim>::get_interpolated_dof_values (const Vector<double> &values,
+ Vector<double> &interpolated_values) const {
const unsigned int total_dofs = dof_handler->get_fe().total_dofs;
Assert (dof_handler != 0, ExcInvalidObject());
else
// otherwise clobber them from the children
{
- dVector tmp1(total_dofs);
- dVector tmp2(total_dofs);
+ Vector<double> tmp1(total_dofs);
+ Vector<double> tmp2(total_dofs);
interpolated_values.clear ();
template <int dim>
void
-DoFCellAccessor<dim>::set_dof_values_by_interpolation (const dVector &local_values,
- dVector &values) const {
+DoFCellAccessor<dim>::set_dof_values_by_interpolation (const Vector<double> &local_values,
+ Vector<double> &values) const {
const unsigned int total_dofs = dof_handler->get_fe().total_dofs;
Assert (dof_handler != 0, ExcInvalidObject());
else
// otherwise distribute them to the children
{
- dVector tmp(total_dofs);
+ Vector<double> tmp(total_dofs);
for (unsigned int child=0; child<GeometryInfo<dim>::children_per_cell;
++child)
#include <grid/dof_constraints.h>
-#include <lac/dsmatrix.h>
-#include <lac/dvector.h>
+#include <lac/sparsematrix.h>
+#include <lac/vector.h>
#include <iostream>
#include <algorithm>
#include <numeric>
-void ConstraintMatrix::condense (const dSMatrixStruct &uncondensed,
- dSMatrixStruct &condensed) const {
+void ConstraintMatrix::condense (const SparseMatrixStruct &uncondensed,
+ SparseMatrixStruct &condensed) const {
Assert (sorted == true, ExcMatrixNotClosed());
Assert (uncondensed.is_compressed() == true, ExcMatrixNotClosed());
Assert (uncondensed.n_rows() == uncondensed.n_cols(),
-void ConstraintMatrix::condense (dSMatrixStruct &sparsity) const {
+void ConstraintMatrix::condense (SparseMatrixStruct &sparsity) const {
Assert (sorted == true, ExcMatrixNotClosed());
Assert (sparsity.is_compressed() == false, ExcMatrixIsClosed());
Assert (sparsity.n_rows() == sparsity.n_cols(),
-void ConstraintMatrix::condense (const dSMatrix &uncondensed,
- dSMatrix &condensed) const {
- const dSMatrixStruct &uncondensed_struct = uncondensed.get_sparsity_pattern ();
+void ConstraintMatrix::condense (const SparseMatrix<double> &uncondensed,
+ SparseMatrix<double> &condensed) const {
+ const SparseMatrixStruct &uncondensed_struct = uncondensed.get_sparsity_pattern ();
Assert (sorted == true, ExcMatrixNotClosed());
Assert (uncondensed_struct.is_compressed() == true, ExcMatrixNotClosed());
-void ConstraintMatrix::condense (dSMatrix &uncondensed) const {
- const dSMatrixStruct &sparsity = uncondensed.get_sparsity_pattern ();
+void ConstraintMatrix::condense (SparseMatrix<double> &uncondensed) const {
+ const SparseMatrixStruct &sparsity = uncondensed.get_sparsity_pattern ();
Assert (sorted == true, ExcMatrixNotClosed());
Assert (sparsity.is_compressed() == true, ExcMatrixNotClosed());
-void ConstraintMatrix::condense (const dVector &uncondensed,
- dVector &condensed) const {
+void ConstraintMatrix::condense (const Vector<double> &uncondensed,
+ Vector<double> &condensed) const {
Assert (sorted == true, ExcMatrixNotClosed());
Assert (condensed.size()+n_constraints() == uncondensed.size(),
ExcWrongDimension());
-void ConstraintMatrix::condense (dVector &vec) const {
+void ConstraintMatrix::condense (Vector<double> &vec) const {
Assert (sorted == true, ExcMatrixNotClosed());
if (lines.size() == 0)
-void ConstraintMatrix::distribute (const dVector &condensed,
- dVector &uncondensed) const {
+void ConstraintMatrix::distribute (const Vector<double> &condensed,
+ Vector<double> &uncondensed) const {
Assert (sorted == true, ExcMatrixNotClosed());
Assert (condensed.size()+n_constraints() == uncondensed.size(),
ExcWrongDimension());
-void ConstraintMatrix::distribute (dVector &vec) const {
+void ConstraintMatrix::distribute (Vector<double> &vec) const {
Assert (sorted == true, ExcMatrixNotClosed());
vector<ConstraintLine>::const_iterator next_constraint = lines.begin();
#include <grid/tria.h>
#include <grid/geometry_info.h>
#include <fe/fe.h>
-#include <lac/dsmatrix.h>
-#include <lac/dvector.h>
+#include <lac/sparsematrix.h>
+#include <lac/vector.h>
+#include <lac/vector.h>
#include <map>
#include <set>
const bool use_constraints,
const vector<int> &starting_points) {
// make the connection graph
- dSMatrixStruct sparsity (n_dofs(), max_couplings_between_dofs());
+ SparseMatrixStruct sparsity (n_dofs(), max_couplings_between_dofs());
make_sparsity_pattern (sparsity);
if (use_constraints)
template <int dim>
-void DoFHandler<dim>::make_sparsity_pattern (dSMatrixStruct &sparsity) const {
+void DoFHandler<dim>::make_sparsity_pattern (SparseMatrixStruct &sparsity) const {
Assert (selected_fe != 0, ExcNoFESelected());
Assert (sparsity.n_rows() == n_dofs(),
ExcDifferentDimensions (sparsity.n_rows(), n_dofs()));
#if deal_II_dimension == 1
template <>
-void DoFHandler<1>::make_boundary_sparsity_pattern (const vector<int> &,
- dSMatrixStruct &) const {
+void DoFHandler<1>::make_boundary_sparsity_pattern (const vector<int> &,
+ SparseMatrixStruct &) const {
Assert (selected_fe != 0, ExcNoFESelected());
Assert (false, ExcInternalError());
};
template <>
-void DoFHandler<1>::make_boundary_sparsity_pattern (const FunctionMap &,
- const vector<int> &,
- dSMatrixStruct &) const {
+void DoFHandler<1>::make_boundary_sparsity_pattern (const FunctionMap &,
+ const vector<int> &,
+ SparseMatrixStruct &) const {
Assert (selected_fe != 0, ExcNoFESelected());
Assert (false, ExcInternalError());
};
template <int dim>
-void DoFHandler<dim>::make_boundary_sparsity_pattern (const vector<int> &dof_to_boundary_mapping,
- dSMatrixStruct &sparsity) const {
+void DoFHandler<dim>::make_boundary_sparsity_pattern (const vector<int> &dof_to_boundary_mapping,
+ SparseMatrixStruct &sparsity) const {
Assert (selected_fe != 0, ExcNoFESelected());
Assert (dof_to_boundary_mapping.size() == n_dofs(), ExcInternalError());
Assert (sparsity.n_rows() == n_boundary_dofs(),
template <int dim>
-void DoFHandler<dim>::make_boundary_sparsity_pattern (const FunctionMap &boundary_indicators,
- const vector<int> &dof_to_boundary_mapping,
- dSMatrixStruct &sparsity) const {
+void DoFHandler<dim>::make_boundary_sparsity_pattern (const FunctionMap &boundary_indicators,
+ const vector<int> &dof_to_boundary_mapping,
+ SparseMatrixStruct &sparsity) const {
Assert (selected_fe != 0, ExcNoFESelected());
Assert (dof_to_boundary_mapping.size() == n_dofs(), ExcInternalError());
Assert (boundary_indicators.find(255) == boundary_indicators.end(),
template <int dim>
void DoFHandler<dim>::make_transfer_matrix (const DoFHandler<dim> &transfer_from,
- dSMatrixStruct &transfer_pattern) const {
+ SparseMatrixStruct &transfer_pattern) const {
Assert (tria->n_cells(0) == transfer_from.tria->n_cells(0),
ExcGridsDoNotMatch());
Assert (*selected_fe == *transfer_from.selected_fe,
template <int dim>
void DoFHandler<dim>::make_transfer_matrix (const DoFHandler<dim> &transfer_from,
- dSMatrix &transfer_matrix) const {
+ SparseMatrix<double> &transfer_matrix) const {
cell_iterator old_cell = transfer_from.begin(0),
new_cell = begin(0);
unsigned int n_coarse_cells = tria->n_cells(0);
template <int dim>
void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterator &old_cell,
const typename DoFHandler<dim>::cell_iterator &new_cell,
- dSMatrixStruct &transfer_pattern) const {
+ SparseMatrixStruct &transfer_pattern) const {
if (!new_cell->active() && !old_cell->active())
// both cells have children; go deeper
for (unsigned int child=0; child<GeometryInfo<dim>::children_per_cell; ++child)
template <int dim>
void DoFHandler<dim>::transfer_cell (const typename DoFHandler<dim>::cell_iterator &old_cell,
const typename DoFHandler<dim>::cell_iterator &new_cell,
- dSMatrix &transfer_matrix) const {
+ SparseMatrix<double> &transfer_matrix) const {
if (!new_cell->active() && !old_cell->active())
// both cells have children; go deeper
for (unsigned int child=0; child<GeometryInfo<dim>::children_per_cell; ++child)
// copy dofs one-by-one
for (unsigned int j=0; j<old_dofs.size(); ++j)
- // use the dSMatrix:: as a workaround
+ // use the SparseMatrix<double>:: as a workaround
// for a bug in egcs
- transfer_matrix.dSMatrix::set (new_dofs[j], old_dofs[j], 1.0);
+ transfer_matrix.SparseMatrix<double>::set (new_dofs[j], old_dofs[j], 1.0);
}
else
if (!new_cell->active() && old_cell->active())
for (unsigned int k=0; k<selected_fe->total_dofs; ++k)
for (unsigned int j=0; j<selected_fe->total_dofs; ++j)
if (selected_fe->prolongate(c)(k,j) != 0.)
- // use the dSMatrix::
+ // use the SparseMatrix<double>::
// as a workaround
// for a bug in egcs
- transfer_matrix.dSMatrix::set (child_dof_indices[k],
- old_dof_indices[j],
- selected_fe->prolongate(c)(k,j));
+ transfer_matrix.SparseMatrix<double>::set (child_dof_indices[k],
+ old_dof_indices[j],
+ selected_fe->prolongate(c)(k,j));
};
} else {
// old cell has children, new one has not
for (unsigned int k=0; k<selected_fe->total_dofs; ++k)
for (unsigned int j=0; j<selected_fe->total_dofs; ++j)
if (selected_fe->restrict(c)(k,j) != 0.)
- // use the dSMatrix:: as
+ // use the SparseMatrix<double>:: as
// a workaround
// for a bug in egcs
- transfer_matrix.dSMatrix::set (new_dof_indices[k],
- child_dof_indices[j],
- selected_fe->restrict(c)(k,j));
+ transfer_matrix.SparseMatrix<double>::set (new_dof_indices[k],
+ child_dof_indices[j],
+ selected_fe->restrict(c)(k,j));
};
};
};
+
template <int dim>
-void DoFHandler<dim>::distribute_cell_to_dof_vector (const dVector &cell_data,
- dVector &dof_data) const {
+template <typename Number>
+void DoFHandler<dim>::distribute_cell_to_dof_vector (const Vector<Number> &cell_data,
+ Vector<double> &dof_data) const {
Assert (cell_data.size()==tria->n_active_cells(),
ExcWrongSize (cell_data.size(), tria->n_active_cells()));
/*-------------- Explicit Instantiations -------------------------------*/
template class DoFHandler<deal_II_dimension>;
+template
+void
+DoFHandler<deal_II_dimension>::distribute_cell_to_dof_vector (const Vector<float> &cell_data,
+ Vector<double> &dof_data) const;
+
+template
+void
+DoFHandler<deal_II_dimension>::distribute_cell_to_dof_vector (const Vector<double> &cell_data,
+ Vector<double> &dof_data) const;
+
#include <grid/tria_iterator.templates.h>
#include <fe/fe.h>
-#include <lac/dvector.h>
-#include <lac/dfmatrix.h>
-#include <lac/dsmatrix.h>
+#include <lac/vector.h>
+#include <lac/fullmatrix.h>
+#include <lac/sparsematrix.h>
template <>
void
-MGDoFCellAccessor<1>::get_mg_dof_values (const dVector &values,
- dVector &dof_values) const {
+MGDoFCellAccessor<1>::get_mg_dof_values (const Vector<double> &values,
+ Vector<double> &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
template <>
void
-MGDoFCellAccessor<2>::get_mg_dof_values (const dVector &values,
- dVector &dof_values) const {
+MGDoFCellAccessor<2>::get_mg_dof_values (const Vector<double> &values,
+ Vector<double> &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
#include <grid/tria.h>
#include <grid/geometry_info.h>
#include <fe/fe.h>
-#include <lac/dsmatrix.h>
+#include <lac/sparsematrix.h>
#include <algorithm>
template <int dim>
void MGDoFHandler<dim>::make_sparsity_pattern (const unsigned int level,
- dSMatrixStruct &sparsity) const {
+ SparseMatrixStruct &sparsity) const {
Assert (selected_fe != 0, ExcNoFESelected());
Assert (sparsity.n_rows() == n_dofs(level),
ExcDifferentDimensions (sparsity.n_rows(), n_dofs(level)));
const RenumberingMethod method,
const vector<int> &starting_points) {
// make the connection graph
- dSMatrixStruct sparsity (n_dofs(level), max_couplings_between_dofs());
+ SparseMatrixStruct sparsity (n_dofs(level), max_couplings_between_dofs());
make_sparsity_pattern (level, sparsity);
int n_dofs = sparsity.n_rows();
template <int dim>
-const dFMatrix &
+const FullMatrix<double> &
FiniteElementBase<dim>::restrict (const unsigned int child) const {
Assert (child<GeometryInfo<dim>::children_per_cell, ExcInvalidIndex(child));
return restriction[child];
template <int dim>
-const dFMatrix &
+const FullMatrix<double> &
FiniteElementBase<dim>::prolongate (const unsigned int child) const {
Assert (child<GeometryInfo<dim>::children_per_cell, ExcInvalidIndex(child));
return prolongation[child];
template <int dim>
-const dFMatrix &
+const FullMatrix<double> &
FiniteElementBase<dim>::constraints () const {
if (dim==1)
Assert ((interface_constraints.m()==1) && (interface_constraints.n()==1),
const bool compute_support_points,
vector<Point<1> > &q_points,
const bool compute_q_points,
- const dFMatrix &,
+ const FullMatrix<double> &,
const vector<vector<Tensor<1,1> > > &,
const Boundary<1> &boundary) const {
Assert ((!compute_jacobians) || (jacobians.size() == unit_points.size()),
const bool ,
vector<Point<1> > &,
const bool,
- const dFMatrix &,
+ const FullMatrix<double> &,
const vector<vector<Tensor<1,1> > > &,
const Boundary<1> &) const {
Assert (false, ExcNotImplemented());
const bool ,
vector<Point<1> > &,
const bool,
- const dFMatrix &,
+ const FullMatrix<double> &,
const vector<vector<Tensor<1,1> > > &,
const Boundary<1> &) const {
Assert (false, ExcNotImplemented());
const bool,
vector<Point<dim> > &,
const bool,
- const dFMatrix &,
+ const FullMatrix<double> &,
const vector<vector<Tensor<1,dim> > > &,
const Boundary<dim> &) const {
Assert (false, ExcPureFunctionCalled());
const bool compute_face_jacobians,
vector<Point<dim> > &normal_vectors,
const bool compute_normal_vectors,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_gradients_transform,
const Boundary<dim> &boundary) const {
Assert (jacobians.size() == unit_points.size(),
const bool compute_face_jacobians,
vector<Point<dim> > &normal_vectors,
const bool compute_normal_vectors,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_gradients_transform,
const Boundary<dim> &boundary) const {
Assert (jacobians.size() == unit_points.size(),
template <>
void FECrissCross<1>::get_local_mass_matrix (const DoFHandler<1>::cell_iterator &,
const Boundary<1> &,
- dFMatrix &) const {
+ FullMatrix<double> &) const {
Assert (false, ExcNotUseful());
};
const bool ,
vector<Point<1> > &,
const bool ,
- const dFMatrix &,
+ const FullMatrix<double> &,
const vector<vector<Tensor<1,1> > > &,
const Boundary<1> &) const {
Assert (false, ExcNotUseful());
template <>
void FECrissCross<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator &cell,
const Boundary<2> &,
- dFMatrix &mass_matrix) const {
+ FullMatrix<double> &mass_matrix) const {
Assert (mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(mass_matrix.n(),total_dofs));
Assert (mass_matrix.m() == total_dofs,
const bool compute_support_points,
vector<Point<dim> > &q_points,
const bool compute_q_points,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &/*shape_grad_transform*/,
const Boundary<dim> &boundary) const {
Assert (jacobians.size() == unit_points.size(),
template <>
void FECubicSub<1>::get_local_mass_matrix (const DoFHandler<1>::cell_iterator &cell,
const Boundary<1> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FECubicSub<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator &cell,
const Boundary<2> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FECubicSub<3>::get_local_mass_matrix (const DoFHandler<3>::cell_iterator &,
const Boundary<3> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <int dim>
-const dFMatrix &
+const FullMatrix<double> &
FEDGLinear<dim>::restrict (const unsigned int child) const {
Assert (false, ExcNotImplemented());
return restriction[child];
template <int dim>
-const dFMatrix &
+const FullMatrix<double> &
FEDGQuadraticSub<dim>::restrict (const unsigned int child) const {
Assert (false, ExcNotImplemented());
return restriction[child];
template <int dim>
-const dFMatrix &
+const FullMatrix<double> &
FEDGCubicSub<dim>::restrict (const unsigned int child) const {
Assert (false, ExcNotImplemented());
return restriction[child];
template <int dim>
-const dFMatrix &
+const FullMatrix<double> &
FEDGQuarticSub<dim>::restrict (const unsigned int child) const {
Assert (false, ExcNotImplemented());
return restriction[child];
template <int dim>
void FEDGConstant<dim>::get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &cell,
const Boundary<dim> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <int dim>
-const dFMatrix &
+const FullMatrix<double> &
FEDGConstant<dim>::restrict (const unsigned int child) const {
Assert (false, ExcNotImplemented());
return restriction[child];
template <>
void FELinear<1>::get_local_mass_matrix (const DoFHandler<1>::cell_iterator &cell,
const Boundary<1> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FELinear<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator &cell,
const Boundary<2> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FELinear<3>::get_local_mass_matrix (const DoFHandler<3>::cell_iterator &,
const Boundary<3> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FEQuadraticSub<1>::get_local_mass_matrix (const DoFHandler<1>::cell_iterator &cell,
const Boundary<1> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FEQuadraticSub<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator &cell,
const Boundary<2> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FEQuadraticSub<3>::get_local_mass_matrix (const DoFHandler<3>::cell_iterator &,
const Boundary<3> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FEQuarticSub<1>::get_local_mass_matrix (const DoFHandler<1>::cell_iterator &cell,
const Boundary<1> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
template <>
void FEQuarticSub<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator &cell,
const Boundary<2> &,
- dFMatrix &local_mass_matrix) const {
+ FullMatrix<double> &local_mass_matrix) const {
Assert (local_mass_matrix.n() == total_dofs,
ExcWrongFieldDimension(local_mass_matrix.n(),total_dofs));
Assert (local_mass_matrix.m() == total_dofs,
const bool compute_support_points,
vector<Point<1> > &q_points,
const bool compute_q_points,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,1> > > &shape_gradients_transform,
const Boundary<1> &boundary) const {
// simply pass down
const bool compute_support_points,
vector<Point<dim> > &q_points,
const bool compute_q_points,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &/*shape_grad_transform*/,
const Boundary<dim> &boundary) const
{
if (compute_jacobians)
{
- dFMatrix M(dim,dim);
+ FullMatrix<double> M(dim,dim);
for (unsigned int l=0; l<n_points; ++l)
{
M.clear ();
#include <base/quadrature_lib.h>
#include <grid/tria_iterator.h>
#include <grid/dof_accessor.h>
- #include <lac/dvector.h>
+ #include <lac/vector.h>
int main () {
Triangulation<2> tria;
fevalues.reinit (dof.begin_active(),b);
- dVector val(4);
+ Vector<double> val(4);
val(2) = 1;
vector<Point<2> > grads(4);
/* $Id$ */
-/* Copyright W. Bangerth, G. Kanschat University of Heidelberg, 1990 */
+/* Copyright W. Bangerth, G. Kanschat University of Heidelberg, 1999 */
#include <fe/fe_system.h>
template <int dim>
void FESystem<dim>::get_local_mass_matrix (const DoFHandler<dim>::cell_iterator &/*cell*/,
const Boundary<dim> &/*boundary*/,
- dFMatrix &/*local_mass_matrix*/) const
+ FullMatrix<double> &/*local_mass_matrix*/) const
{
Assert(false, ExcNotImplemented());
/*
// first get the local mass matrix for
// the base object
- dFMatrix base_mass_matrix (base_element->total_dofs,
+ FullMatrix<double> base_mass_matrix (base_element->total_dofs,
base_element->total_dofs);
base_element->get_local_mass_matrix (cell, boundary, base_mass_matrix);
const bool compute_support_points,
vector<Point<dim> > &q_points,
const bool compute_q_points,
- const dFMatrix &shape_values_transform,
+ const FullMatrix<double> &shape_values_transform,
const vector<vector<Tensor<1,dim> > > &shape_grad_transform,
const Boundary<dim> &boundary) const
{
#include <grid/tria_accessor.h>
#include <grid/tria_boundary.h>
#include <grid/dof_accessor.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
/*------------------------------- FEValuesBase ---------------------------*/
n_quadrature_points (n_q_points),
total_dofs (total_dofs),
n_transform_functions (n_transform_functions),
- shape_values (n_values_arrays, dFMatrix(total_dofs, n_q_points)),
+ shape_values (n_values_arrays, FullMatrix<double>(total_dofs, n_q_points)),
shape_gradients (total_dofs, vector<Tensor<1,dim> >(n_q_points)),
shape_2nd_derivatives (total_dofs, vector<Tensor<2,dim> >(n_q_points)),
weights (n_q_points, 0),
jacobi_matrices (n_q_points, Tensor<2,dim>()),
jacobi_matrices_grad (n_q_points, Tensor<3,dim>()),
shape_values_transform (n_values_arrays,
- dFMatrix(n_transform_functions,
- n_quadrature_points)),
+ FullMatrix<double>(n_transform_functions,
+ n_quadrature_points)),
selected_dataset (0),
update_flags (update_flags),
fe(&fe) {};
template <int dim>
-void FEValuesBase<dim>::get_function_values (const dVector &fe_function,
- vector<double> &values) const
+void FEValuesBase<dim>::get_function_values (const Vector<double> &fe_function,
+ vector<double> &values) const
{
Assert (fe->n_components == 1,
ExcWrongNoOfComponents());
// get function values of dofs
// on this cell
- dVector dof_values (total_dofs);
+ Vector<double> dof_values (total_dofs);
present_cell->get_dof_values (fe_function, dof_values);
// initialize with zero
template <int dim>
-void FEValuesBase<dim>::get_function_values (const dVector &fe_function,
+void FEValuesBase<dim>::get_function_values (const Vector<double> &fe_function,
vector< vector<double> > &values) const
{
Assert (fe->n_components == values.size(),
// get function values of dofs
// on this cell
- dVector dof_values (total_dofs);
+ Vector<double> dof_values (total_dofs);
present_cell->get_dof_values (fe_function, dof_values);
// initialize with zero
template <int dim>
-void FEValuesBase<dim>::get_function_grads (const dVector &fe_function,
+void FEValuesBase<dim>::get_function_grads (const Vector<double> &fe_function,
vector<Tensor<1,dim> > &gradients) const {
Assert (gradients.size() == n_quadrature_points,
ExcWrongVectorSize(gradients.size(), n_quadrature_points));
// get function values of dofs
// on this cell
- dVector dof_values (total_dofs);
+ Vector<double> dof_values (total_dofs);
present_cell->get_dof_values (fe_function, dof_values);
// initialize with zero
template <int dim>
-void FEValuesBase<dim>::get_function_2nd_derivatives (const dVector &fe_function,
+void FEValuesBase<dim>::get_function_2nd_derivatives (const Vector<double> &fe_function,
vector<Tensor<2,dim> > &second_derivatives) const {
Assert (second_derivatives.size() == n_quadrature_points,
ExcWrongVectorSize(second_derivatives.size(), n_quadrature_points));
// get function values of dofs
// on this cell
- dVector dof_values (total_dofs);
+ Vector<double> dof_values (total_dofs);
present_cell->get_dof_values (fe_function, dof_values);
// initialize with zero
#include <grid/tria_iterator.h>
#include <grid/geometry_info.h>
#include <basic/magic_numbers.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <iostream>
#include <algorithm>
#include <numeric>
template <int dim>
-void Triangulation<dim>::refine (const dVector &criteria,
- const double threshold) {
+void Triangulation<dim>::refine (const Vector<float> &criteria,
+ const double threshold) {
Assert (criteria.size() == n_active_cells(),
ExcInvalidVectorSize(criteria.size(), n_active_cells()));
template <int dim>
-void Triangulation<dim>::coarsen (const dVector &criteria,
- const double threshold) {
+void Triangulation<dim>::coarsen (const Vector<float> &criteria,
+ const double threshold) {
Assert (criteria.size() == n_active_cells(),
ExcInvalidVectorSize(criteria.size(), n_active_cells()));
template <int dim>
-void Triangulation<dim>::refine_and_coarsen_fixed_number (const dVector &criteria,
- const double top_fraction,
- const double bottom_fraction) {
+void Triangulation<dim>::refine_and_coarsen_fixed_number (const Vector<float> &criteria,
+ const double top_fraction,
+ const double bottom_fraction) {
// correct number of cells is
// checked in #refine#
Assert ((top_fraction>=0) && (top_fraction<=1), ExcInvalidParameterValue());
const int coarsen_cells = max(static_cast<int>(bottom_fraction*criteria.size()),
1);
- dVector tmp(criteria);
+ Vector<float> tmp(criteria);
nth_element (tmp.begin(), tmp.begin()+refine_cells,
tmp.end(),
greater<double>());
template <int dim>
void
-Triangulation<dim>::refine_and_coarsen_fixed_fraction (const dVector &criteria,
- const double top_fraction,
- const double bottom_fraction) {
+Triangulation<dim>::refine_and_coarsen_fixed_fraction (const Vector<float> &criteria,
+ const double top_fraction,
+ const double bottom_fraction) {
// correct number of cells is
// checked in #refine#
Assert ((top_fraction>=0) && (top_fraction<=1), ExcInvalidParameterValue());
// error, which is what we have to sum
// up and compare with
// #fraction_of_error*total_error#.
- dVector tmp(criteria);
+ Vector<float> tmp(criteria);
const double total_error = tmp.l1_norm();
-
- dVector partial_sums(criteria.size());
+
+ Vector<float> partial_sums(criteria.size());
sort (tmp.begin(), tmp.end(), greater<double>());
partial_sum (tmp.begin(), tmp.end(), partial_sums.begin());
// compute thresholds
- dVector::const_iterator p;
+ Vector<float>::const_iterator p;
double top_threshold, bottom_threshold;
p = lower_bound (partial_sums.begin(), partial_sums.end(),
top_fraction*total_error);
#include <grid/tria_iterator.templates.h>
#include <fe/fe.h>
-#include <lac/dvector.h>
-#include <lac/dfmatrix.h>
-#include <lac/dsmatrix.h>
+#include <lac/vector.h>
+#include <lac/fullmatrix.h>
+#include <lac/sparsematrix.h>
template <>
void
-MGDoFCellAccessor<1>::get_mg_dof_values (const dVector &values,
- dVector &dof_values) const {
+MGDoFCellAccessor<1>::get_mg_dof_values (const Vector<double> &values,
+ Vector<double> &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
template <>
void
-MGDoFCellAccessor<2>::get_mg_dof_values (const dVector &values,
- dVector &dof_values) const {
+MGDoFCellAccessor<2>::get_mg_dof_values (const Vector<double> &values,
+ Vector<double> &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
#include <grid/tria.h>
#include <grid/geometry_info.h>
#include <fe/fe.h>
-#include <lac/dsmatrix.h>
+#include <lac/sparsematrix.h>
#include <algorithm>
template <int dim>
void MGDoFHandler<dim>::make_sparsity_pattern (const unsigned int level,
- dSMatrixStruct &sparsity) const {
+ SparseMatrixStruct &sparsity) const {
Assert (selected_fe != 0, ExcNoFESelected());
Assert (sparsity.n_rows() == n_dofs(level),
ExcDifferentDimensions (sparsity.n_rows(), n_dofs(level)));
const RenumberingMethod method,
const vector<int> &starting_points) {
// make the connection graph
- dSMatrixStruct sparsity (n_dofs(level), max_couplings_between_dofs());
+ SparseMatrixStruct sparsity (n_dofs(level), max_couplings_between_dofs());
make_sparsity_pattern (level, sparsity);
int n_dofs = sparsity.n_rows();
// Copyright Guido Kanschat, Universitaet Heidelberg, 1999
#include <numerics/multigrid.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
void
-MultiGrid::vmult(dVector& dst, const dVector& src) const
+MultiGrid::vmult(Vector<double>& dst, const Vector<double>& src) const
{
dst = 0.;
void
-MultiGrid::precondition(dVector& dst, const dVector& src) const
+MultiGrid::precondition(Vector<double>& dst, const Vector<double>& src) const
{
copy_to_mg(s,src);
copy_to_mg(d,dst);
void MultiGrid::post_smooth(unsigned level,
- dVector& dst, const dVector& src,
+ Vector<double>& dst, const Vector<double>& src,
unsigned steps) const
{
smooth(level, dst, src, steps);
#include <grid/tria_iterator.h>
#include <grid/tria_iterator.templates.h>
#include <fe/fe.h>
-#include <lac/dfmatrix.h>
-#include <lac/dvector.h>
-#include <lac/dsmatrix.h>
+#include <lac/fullmatrix.h>
+#include <lac/vector.h>
+#include <lac/sparsematrix.h>
template <int dim>
template <int dim>
-void Equation<dim>::assemble (dFMatrix &,
- dVector &,
+void Equation<dim>::assemble (FullMatrix<double> &,
+ Vector<double> &,
const FEValues<dim> &,
const typename DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void Equation<dim>::assemble (dFMatrix &,
+void Equation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const typename DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void Equation<dim>::assemble (dVector &,
+void Equation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const typename DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
AssemblerData<dim>::AssemblerData (const DoFHandler<dim> &dof,
const bool assemble_matrix,
const bool assemble_rhs,
- dSMatrix &matrix,
- dVector &rhs_vector,
+ SparseMatrix<double> &matrix,
+ Vector<double> &rhs_vector,
const Quadrature<dim> &quadrature,
const FiniteElement<dim> &fe,
const UpdateFlags &update_flags,
const AssemblerData<dim> *local_data) :
DoFCellAccessor<dim> (tria,level,index, &local_data->dof),
cell_matrix (dof_handler->get_fe().total_dofs),
- cell_vector (dVector(dof_handler->get_fe().total_dofs)),
+ cell_vector (Vector<double>(dof_handler->get_fe().total_dofs)),
assemble_matrix (local_data->assemble_matrix),
assemble_rhs (local_data->assemble_rhs),
matrix(local_data->matrix),
#include <base/function.h>
#include <fe/fe.h>
#include <base/quadrature.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <lac/solver_cg.h>
#include <lac/vector_memory.h>
Assert ((tria!=0) && (dof_handler!=0), ExcNoTriaSelected());
SolverControl control(4000, 1e-16);
- PrimitiveVectorMemory<dVector> memory;
- SolverCG<dSMatrix,dVector> cg(control,memory);
+ PrimitiveVectorMemory<Vector<double> > memory;
+ SolverCG<SparseMatrix<double>,Vector<double> > cg(control,memory);
// solve
cg.solve (system_matrix, solution, right_hand_side);
#include <base/quadrature_lib.h>
#include <fe/fe_values.h>
#include <fe/fe.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <map>
#include <iostream>
template <int dim>
-DataOut<dim>::DataEntry::DataEntry (const dVector *data,
+DataOut<dim>::DataEntry::DataEntry (const Vector<double> *data,
const string name,
const string units) :
data(data), name(name), units(units) {};
template <int dim>
-void DataOut<dim>::add_data_vector (const dVector &vec,
+void DataOut<dim>::add_data_vector (const Vector<double> &vec,
const string &name,
const string &units) {
Assert (dofs != 0, ExcNoDoFHandlerSelected ());
#include <grid/tria_iterator.h>
#include <grid/dof_accessor.h>
#include <grid/geometry_info.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
+#include <lac/vector.h>
#include <numeric>
#include <algorithm>
const FiniteElement<1> &,
const Boundary<1> &,
const FunctionMap &,
- const dVector &,
- dVector &,
+ const Vector<double> &,
+ Vector<float> &,
const Function<1> *) {
Assert(false, ExcNotImplemented());
};
const FiniteElement<dim> &fe,
const Boundary<dim> &boundary,
const FunctionMap &neumann_bc,
- const dVector &solution,
- dVector &error,
+ const Vector<double> &solution,
+ Vector<float> &error,
const Function<dim> *coefficient) {
Assert (neumann_bc.find(255) == neumann_bc.end(),
ExcInvalidBoundaryIndicator());
FEFaceValues<1> &,
FEFaceValues<1> &,
FaceIntegrals &,
- const dVector &,
+ const Vector<double> &,
const Function<1> *) {
Assert (false, ExcInternalError());
};
FEFaceValues<1> &,
FESubfaceValues<1> &,
FaceIntegrals &,
- const dVector &,
+ const Vector<double> &,
const Function<1> *) {
Assert (false, ExcInternalError());
};
FEFaceValues<dim> &fe_face_values_cell,
FEFaceValues<dim> &fe_face_values_neighbor,
FaceIntegrals &face_integrals,
- const dVector &solution,
+ const Vector<double> &solution,
const Function<dim> *coefficient) {
const DoFHandler<dim>::face_iterator face = cell->face(face_no);
FEFaceValues<dim> &fe_face_values,
FESubfaceValues<dim> &fe_subface_values,
FaceIntegrals &face_integrals,
- const dVector &solution,
+ const Vector<double> &solution,
const Function<dim> *coefficient) {
const DoFHandler<dim>::cell_iterator neighbor = cell->neighbor(face_no);
Assert (neighbor.state() == valid, ExcInternalError());
#include <fe/fe_values.h>
#include <numerics/matrices.h>
#include <numerics/assembler.h>
-#include <lac/dvector.h>
-#include <lac/dsmatrix.h>
+#include <lac/vector.h>
+#include <lac/sparsematrix.h>
#include <algorithm>
#include <set>
void MatrixCreator<dim>::create_mass_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const Function<dim> * const a) {
const FiniteElement<dim> &fe = dof.get_fe();
- dVector dummy; // no entries, should give an error if accessed
+ Vector<double> dummy; // no entries, should give an error if accessed
UpdateFlags update_flags = update_JxW_values;
if (a != 0)
update_flags = UpdateFlags (update_flags | update_q_points);
void MatrixCreator<dim>::create_mass_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const Function<dim> &rhs,
- dVector &rhs_vector,
+ Vector<double> &rhs_vector,
const Function<dim> * const a) {
const FiniteElement<dim> &fe = dof.get_fe();
template <int dim>
void MatrixCreator<dim>::create_mass_matrix (const DoFHandler<dim> &dof,
const Boundary<dim> &boundary,
- dSMatrix &matrix) {
+ SparseMatrix<double> &matrix) {
const FiniteElement<dim> &fe = dof.get_fe();
const unsigned int total_dofs = fe.total_dofs;
- dFMatrix local_mass_matrix (total_dofs, total_dofs);
- vector<int> dofs_on_this_cell (total_dofs);
+ FullMatrix<double> local_mass_matrix (total_dofs, total_dofs);
+ vector<int> dofs_on_this_cell (total_dofs);
DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
endc = dof.end();
void MatrixCreator<1>::create_boundary_mass_matrix (const DoFHandler<1> &,
const Quadrature<0> &,
const Boundary<1> &,
- dSMatrix &,
+ SparseMatrix<double> &,
const FunctionMap &,
- dVector &,
+ Vector<double> &,
vector<int> &,
const Function<1> *) {
Assert (false, ExcNotImplemented());
void MatrixCreator<dim>::create_boundary_mass_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim-1> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const FunctionMap &rhs,
- dVector &rhs_vector,
+ Vector<double> &rhs_vector,
vector<int> &dof_to_boundary_mapping,
const Function<dim> *a) {
const FiniteElement<dim> &fe = dof.get_fe();
const unsigned int dofs_per_cell = fe.total_dofs,
dofs_per_face = fe.dofs_per_face;
- dFMatrix cell_matrix(dofs_per_cell, dofs_per_cell);
- dVector cell_vector(dofs_per_cell);
+ FullMatrix<double> cell_matrix(dofs_per_cell, dofs_per_cell);
+ Vector<double> cell_vector(dofs_per_cell);
UpdateFlags update_flags = UpdateFlags (update_JxW_values | update_q_points);
fe_values.reinit (cell, face, boundary);
- const dFMatrix &values = fe_values.get_shape_values ();
- const vector<double> &weights = fe_values.get_JxW_values ();
- vector<double> rhs_values (fe_values.n_quadrature_points);
+ const FullMatrix<double> &values = fe_values.get_shape_values ();
+ const vector<double> &weights = fe_values.get_JxW_values ();
+ vector<double> rhs_values (fe_values.n_quadrature_points);
rhs.find(cell->face(face)->boundary_indicator())
->second->value_list (fe_values.get_quadrature_points(), rhs_values);
void MatrixCreator<dim>::create_laplace_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const Function<dim> * const a) {
const FiniteElement<dim> &fe = dof.get_fe();
- dVector dummy; // no entries, should give an error if accessed
+ Vector<double> dummy; // no entries, should give an error if accessed
UpdateFlags update_flags = UpdateFlags(update_gradients |
update_JxW_values);
if (a != 0)
void MatrixCreator<dim>::create_laplace_matrix (const DoFHandler<dim> &dof,
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
- dSMatrix &matrix,
+ SparseMatrix<double> &matrix,
const Function<dim> &rhs,
- dVector &rhs_vector,
+ Vector<double> &rhs_vector,
const Function<dim> * const a) {
const FiniteElement<dim> &fe = dof.get_fe();
template <int dim>
void MatrixTools<dim>::apply_boundary_values (const map<int,double> &boundary_values,
- dSMatrix &matrix,
- dVector &solution,
- dVector &right_hand_side) {
+ SparseMatrix<double> &matrix,
+ Vector<double> &solution,
+ Vector<double> &right_hand_side) {
Assert (matrix.n() == matrix.m(),
ExcDimensionsDontMatch(matrix.n(), matrix.m()));
Assert (matrix.n() == right_hand_side.size(),
map<int,double>::const_iterator dof = boundary_values.begin(),
endd = boundary_values.end();
const unsigned int n_dofs = matrix.m();
- const dSMatrixStruct &sparsity = matrix.get_sparsity_pattern();
+ const SparseMatrixStruct &sparsity = matrix.get_sparsity_pattern();
const unsigned int *sparsity_rowstart = sparsity.get_rowstart_indices();
const int *sparsity_colnums = sparsity.get_column_numbers();
break;
};
- // use the dSMatrix:: as a workaround
+ // use the SparseMatrix<double>:: as a workaround
// for a bug in egcs
- matrix.dSMatrix::set(dof_number, dof_number,
- first_diagonal_entry);
+ matrix.SparseMatrix<double>::set(dof_number, dof_number,
+ first_diagonal_entry);
new_rhs = right_hand_side(dof_number)
= dof->second * first_diagonal_entry;
};
template <int dim>
-void MassMatrix<dim>::assemble (dFMatrix &cell_matrix,
- const FEValues<dim> &fe_values,
+void MassMatrix<dim>::assemble (FullMatrix<double> &cell_matrix,
+ const FEValues<dim> &fe_values,
const typename DoFHandler<dim>::cell_iterator &) const {
const unsigned int total_dofs = fe_values.total_dofs,
n_q_points = fe_values.n_quadrature_points;
ExcWrongSize(cell_matrix.m(), total_dofs));
Assert (cell_matrix.all_zero(), ExcObjectNotEmpty());
- const dFMatrix &values = fe_values.get_shape_values ();
- const vector<double> &weights = fe_values.get_JxW_values ();
+ const FullMatrix<double> &values = fe_values.get_shape_values ();
+ const vector<double> &weights = fe_values.get_JxW_values ();
if (coefficient != 0)
template <int dim>
-void MassMatrix<dim>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void MassMatrix<dim>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &) const {
Assert (right_hand_side != 0, ExcNoRHSSelected());
Assert (cell_matrix.all_zero(), ExcObjectNotEmpty());
Assert (rhs.all_zero(), ExcObjectNotEmpty());
- const dFMatrix &values = fe_values.get_shape_values ();
- const vector<double> &weights = fe_values.get_JxW_values ();
- vector<double> rhs_values (fe_values.n_quadrature_points);
+ const FullMatrix<double> &values = fe_values.get_shape_values ();
+ const vector<double> &weights = fe_values.get_JxW_values ();
+ vector<double> rhs_values (fe_values.n_quadrature_points);
right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values);
if (coefficient != 0)
template <int dim>
-void MassMatrix<dim>::assemble (dVector &rhs,
+void MassMatrix<dim>::assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &) const {
Assert (right_hand_side != 0, ExcNoRHSSelected());
Assert (rhs.size() == total_dofs, ExcWrongSize(rhs.size(), total_dofs));
Assert (rhs.all_zero(), ExcObjectNotEmpty());
- const dFMatrix &values = fe_values.get_shape_values ();
- const vector<double> &weights = fe_values.get_JxW_values ();
- vector<double> rhs_values(fe_values.n_quadrature_points);
+ const FullMatrix<double> &values = fe_values.get_shape_values ();
+ const vector<double> &weights = fe_values.get_JxW_values ();
+ vector<double> rhs_values(fe_values.n_quadrature_points);
right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values);
for (unsigned int point=0; point<n_q_points; ++point)
template <int dim>
-void LaplaceMatrix<dim>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
- const FEValues<dim> &fe_values,
+void LaplaceMatrix<dim>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
+ const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &) const {
Assert (right_hand_side != 0, ExcNoRHSSelected());
Assert (rhs.all_zero(), ExcObjectNotEmpty());
const vector<vector<Tensor<1,dim> > >&gradients = fe_values.get_shape_grads ();
- const dFMatrix &values = fe_values.get_shape_values ();
+ const FullMatrix<double> &values = fe_values.get_shape_values ();
vector<double> rhs_values(fe_values.n_quadrature_points);
const vector<double> &weights = fe_values.get_JxW_values ();
right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values);
template <int dim>
-void LaplaceMatrix<dim>::assemble (dFMatrix &cell_matrix,
+void LaplaceMatrix<dim>::assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &) const {
const unsigned int total_dofs = fe_values.total_dofs,
template <int dim>
-void LaplaceMatrix<dim>::assemble (dVector &rhs,
+void LaplaceMatrix<dim>::assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &) const {
Assert (right_hand_side != 0, ExcNoRHSSelected());
Assert (rhs.size() == total_dofs, ExcWrongSize(rhs.size(), total_dofs));
Assert (rhs.all_zero(), ExcObjectNotEmpty());
- const dFMatrix &values = fe_values.get_shape_values ();
- const vector<double> &weights = fe_values.get_JxW_values ();
+ const FullMatrix<double> &values = fe_values.get_shape_values ();
+ const vector<double> &weights = fe_values.get_JxW_values ();
vector<double> rhs_values(fe_values.n_quadrature_points);
right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values);
template<int dim> void
MatrixCreator<dim>::create_interpolation_matrix(const FiniteElement<dim> &high,
const FiniteElement<dim> &low,
- dFMatrix& result)
+ FullMatrix<double>& result)
{
result.reinit (low.total_dofs, high.total_dofs);
// Copyright Guido Kanschat, Universitaet Heidelberg, 1999
#include <numerics/multigrid.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
void
-MultiGrid::vmult(dVector& dst, const dVector& src) const
+MultiGrid::vmult(Vector<double>& dst, const Vector<double>& src) const
{
dst = 0.;
void
-MultiGrid::precondition(dVector& dst, const dVector& src) const
+MultiGrid::precondition(Vector<double>& dst, const Vector<double>& src) const
{
copy_to_mg(s,src);
copy_to_mg(d,dst);
void MultiGrid::post_smooth(unsigned level,
- dVector& dst, const dVector& src,
+ Vector<double>& dst, const Vector<double>& src,
unsigned steps) const
{
smooth(level, dst, src, steps);
#include <numerics/assembler.h>
#include <numerics/vectors.h>
#include <numerics/matrices.h>
-#include <lac/dvector.h>
-#include <lac/dsmatrix.h>
+#include <lac/vector.h>
+#include <lac/sparsematrix.h>
#include <lac/solver_cg.h>
#include <lac/vector_memory.h>
void VectorTools<dim>::interpolate (const DoFHandler<dim> &dof,
const Boundary<dim> &boundary,
const Function<dim> &function,
- dVector &vec)
+ Vector<double> &vec)
{
const FiniteElement<dim> &fe = dof.get_fe();
template <int dim> void
VectorTools<dim>::interpolate(const DoFHandler<dim> &high_dof,
const DoFHandler<dim> &low_dof,
- const dFMatrix &transfer,
- const dVector &high,
- dVector &low)
+ const FullMatrix<double> &transfer,
+ const Vector<double> &high,
+ Vector<double> &low)
{
- dVector cell_high(high_dof.get_fe().total_dofs);
- dVector cell_low(low_dof.get_fe().total_dofs);
+ Vector<double> cell_high(high_dof.get_fe().total_dofs);
+ Vector<double> cell_low(low_dof.get_fe().total_dofs);
DoFHandler<dim>::active_cell_iterator h = high_dof.begin_active();
DoFHandler<dim>::active_cell_iterator l = low_dof.begin_active();
const Boundary<1> &,
const Quadrature<1> &,
const Function<1> &,
- dVector &,
+ Vector<double> &,
const bool ,
const Quadrature<0> &,
const bool ) {
const Boundary<dim> &boundary,
const Quadrature<dim> &q,
const Function<dim> &function,
- dVector &vec,
+ Vector<double> &vec,
const bool enforce_zero_boundary,
const Quadrature<dim-1> &q_boundary,
const bool project_to_boundary_first) {
// set up mass matrix and right hand side
vec.reinit (dof.n_dofs());
- dSMatrixStruct sparsity(dof.n_dofs(),
- dof.n_dofs(),
- dof.max_couplings_between_dofs());
+ SparseMatrixStruct sparsity(dof.n_dofs(),
+ dof.n_dofs(),
+ dof.max_couplings_between_dofs());
dof.make_sparsity_pattern (sparsity);
constraints.condense (sparsity);
- dSMatrix mass_matrix (sparsity);
- dVector tmp (mass_matrix.n());
+ SparseMatrix<double> mass_matrix (sparsity);
+ Vector<double> tmp (mass_matrix.n());
MatrixCreator<dim>::create_mass_matrix (dof, boundary, mass_matrix);
VectorTools<dim>::create_right_hand_side (dof, q, boundary,
function, tmp);
mass_matrix, vec, tmp);
SolverControl control(1000,1e-16);
- PrimitiveVectorMemory<dVector> memory;
- SolverCG<dSMatrix,dVector> cg(control,memory);
+ PrimitiveVectorMemory<Vector<double> > memory;
+ SolverCG<SparseMatrix<double>,Vector<double> > cg(control,memory);
// solve
cg.solve (mass_matrix, vec, tmp);
const Quadrature<dim> &q,
const Boundary<dim> &boundary,
const Function<dim> &rhs,
- dVector &rhs_vector) {
+ Vector<double> &rhs_vector) {
const FiniteElement<dim> &fe = dof.get_fe();
UpdateFlags update_flags = UpdateFlags(update_q_points |
update_JxW_values);
- dSMatrix dummy;
+ SparseMatrix<double> dummy;
const AssemblerData<dim> data (dof,
false, true,
dummy, rhs_vector,
dof.map_dof_to_boundary_indices (boundary_functions, dof_to_boundary_mapping);
// set up sparsity structure
- dSMatrixStruct sparsity(dof.n_boundary_dofs(boundary_functions),
- dof.max_couplings_between_boundary_dofs());
+ SparseMatrixStruct sparsity(dof.n_boundary_dofs(boundary_functions),
+ dof.max_couplings_between_boundary_dofs());
dof.make_boundary_sparsity_pattern (boundary_functions, dof_to_boundary_mapping,
sparsity);
// make mass matrix and right hand side
- dSMatrix mass_matrix(sparsity);
- dVector rhs(sparsity.n_rows());
+ SparseMatrix<double> mass_matrix(sparsity);
+ Vector<double> rhs(sparsity.n_rows());
MatrixTools<dim>::create_boundary_mass_matrix (dof, q, boundary,
Assert (dim<3, ExcNotImplemented());
- dVector boundary_projection (rhs.size());
+ Vector<double> boundary_projection (rhs.size());
SolverControl control(1000, 1e-16);
- PrimitiveVectorMemory<dVector> memory;
- SolverCG<dSMatrix,dVector> cg(control,memory);
+ PrimitiveVectorMemory<Vector<double> > memory;
+ SolverCG<SparseMatrix<double>,Vector<double> > cg(control,memory);
// solve
cg.solve (mass_matrix, boundary_projection, rhs);
template <int dim>
void VectorTools<dim>::integrate_difference (const DoFHandler<dim> &dof,
- const dVector &fe_function,
+ const Vector<double> &fe_function,
const Function<dim> &exact_solution,
- dVector &difference,
+ Vector<float> &difference,
const Quadrature<dim> &q,
const NormType &norm,
const Boundary<dim> &boundary) {
};
};
-template <int dim>
-void VectorTools<dim>::integrate_difference (const DoFHandler<dim> &dof,
- const dVector &fe_function,
- const TensorFunction<1,dim> &exact_solution,
- dVector &difference,
- const Quadrature<dim> &q,
- const FiniteElement<dim> &fe,
- const NormType &norm,
- const Boundary<dim> &boundary) {
- Assert (fe == dof.get_fe(), ExcInvalidFE());
+
+
+// template <int dim>
+// void VectorTools<dim>::integrate_difference (const DoFHandler<dim> &dof,
+// const Vector<double> &fe_function,
+// const TensorFunction<1,dim>&exact_solution,
+// Vector<float> &difference,
+// const Quadrature<dim> &q,
+// const FiniteElement<dim> &fe,
+// const NormType &norm,
+// const Boundary<dim> &boundary) {
+// Assert (fe == dof.get_fe(), ExcInvalidFE());
- difference.reinit (dof.get_tria().n_active_cells());
+// difference.reinit (dof.get_tria().n_active_cells());
- UpdateFlags update_flags = UpdateFlags (update_q_points |
- update_JxW_values);
- if ((norm==H1_seminorm) || (norm==H1_norm))
- update_flags = UpdateFlags (update_flags | update_gradients);
- FEValues<dim> fe_values(fe, q, update_flags);
+// UpdateFlags update_flags = UpdateFlags (update_q_points |
+// update_JxW_values);
+// if ((norm==H1_seminorm) || (norm==H1_norm))
+// update_flags = UpdateFlags (update_flags | update_gradients);
+// FEValues<dim> fe_values(fe, q, update_flags);
- // loop over all cells
- DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
- endc = dof.end();
- for (unsigned int index=0; cell != endc; ++cell, ++index)
- {
- double diff=0;
- // initialize for this cell
- fe_values.reinit (cell, boundary);
-
- switch (norm)
- {
- case mean:
- case L1_norm:
- case L2_norm:
- case Linfty_norm:
- case H1_norm:
- {
- // we need the finite element
- // function \psi at the different
- // integration points. Compute
- // it like this:
- // \psi(x_j)=\sum_i v_i \phi_i(x_j)
- // with v_i the nodal values of the
- // fe_function and \phi_i(x_j) the
- // matrix of the trial function
- // values at the integration point
- // x_j. Then the vector
- // of the \psi(x_j) is v*Phi with
- // v being the vector of nodal
- // values on this cell and Phi
- // the matrix.
- //
- // we then need the difference:
- // reference_function(x_j)-\psi_j
- // and assign that to the vector
- // \psi.
- const unsigned int n_q_points = q.n_quadrature_points;
- vector<double> psi (n_q_points);
-
- // in praxi: first compute
- // exact fe_function vector
- exact_solution.value_list (fe_values.get_quadrature_points(),
- psi);
- // then subtract finite element
- // fe_function
- if (true)
- {
- vector< vector<double> > function_values (fe.n_components, vector<double>(n_q_points, 0));
- fe_values.get_function_values (fe_function, function_values);
-
- transform (psi.begin(), psi.end(),
- function_values.begin(),
- psi.begin(),
- minus<double>());
- };
-
- // for L1_norm and Linfty_norm:
- // take absolute
- // value, for the L2_norm take
- // square of psi
- switch (norm)
- {
- case mean:
- break;
- case L1_norm:
- case Linfty_norm:
- transform (psi.begin(), psi.end(),
- psi.begin(), ptr_fun(fabs));
- break;
- case L2_norm:
- case H1_norm:
- transform (psi.begin(), psi.end(),
- psi.begin(), ptr_fun(sqr));
- break;
- default:
- Assert (false, ExcNotImplemented());
- };
-
- // ok, now we have the integrand,
- // let's compute the integral,
- // which is
- // sum_j psi_j JxW_j
- // (or |psi_j| or |psi_j|^2
- switch (norm)
- {
- case mean:
- case L1_norm:
- diff = inner_product (psi.begin(), psi.end(),
- fe_values.get_JxW_values().begin(),
- 0.0);
- break;
- case L2_norm:
- case H1_norm:
- diff = sqrt(inner_product (psi.begin(), psi.end(),
- fe_values.get_JxW_values().begin(),
- 0.0));
- break;
- case Linfty_norm:
- diff = *max_element (psi.begin(), psi.end());
- break;
- default:
- Assert (false, ExcNotImplemented());
- };
-
- // note: the H1_norm uses the result
- // of the L2_norm and control goes
- // over to the next case statement!
- if (norm != H1_norm)
- break;
- };
-
- case H1_seminorm:
- {
- // note: the computation of the
- // H1_norm starts at the previous
- // case statement, but continues
- // here!
-
- // for H1_norm: re-square L2_norm.
- diff = sqr(diff);
-
- // same procedure as above, but now
- // psi is a vector of gradients
- const unsigned int n_q_points = q.n_quadrature_points;
- vector<Tensor<1,dim> > psi (n_q_points);
-
- // in praxi: first compute
- // exact fe_function vector
- exact_solution.gradient_list (fe_values.get_quadrature_points(),
- psi);
+// // loop over all cells
+// DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
+// endc = dof.end();
+// for (unsigned int index=0; cell != endc; ++cell, ++index)
+// {
+// double diff=0;
+// // initialize for this cell
+// fe_values.reinit (cell, boundary);
+
+// switch (norm)
+// {
+// case mean:
+// case L1_norm:
+// case L2_norm:
+// case Linfty_norm:
+// case H1_norm:
+// {
+// // we need the finite element
+// // function \psi at the different
+// // integration points. Compute
+// // it like this:
+// // \psi(x_j)=\sum_i v_i \phi_i(x_j)
+// // with v_i the nodal values of the
+// // fe_function and \phi_i(x_j) the
+// // matrix of the trial function
+// // values at the integration point
+// // x_j. Then the vector
+// // of the \psi(x_j) is v*Phi with
+// // v being the vector of nodal
+// // values on this cell and Phi
+// // the matrix.
+// //
+// // we then need the difference:
+// // reference_function(x_j)-\psi_j
+// // and assign that to the vector
+// // \psi.
+// const unsigned int n_q_points = q.n_quadrature_points;
+// vector<double> psi (n_q_points);
+
+// // in praxi: first compute
+// // exact fe_function vector
+// exact_solution.value_list (fe_values.get_quadrature_points(),
+// psi);
+// // then subtract finite element
+// // fe_function
+// if (true)
+// {
+// vector< vector<double> > function_values (fe.n_components, vector<double>(n_q_points, 0));
+// fe_values.get_function_values (fe_function, function_values);
+
+// transform (psi.begin(), psi.end(),
+// function_values.begin(),
+// psi.begin(),
+// minus<double>());
+// };
+
+// // for L1_norm and Linfty_norm:
+// // take absolute
+// // value, for the L2_norm take
+// // square of psi
+// switch (norm)
+// {
+// case mean:
+// break;
+// case L1_norm:
+// case Linfty_norm:
+// transform (psi.begin(), psi.end(),
+// psi.begin(), ptr_fun(fabs));
+// break;
+// case L2_norm:
+// case H1_norm:
+// transform (psi.begin(), psi.end(),
+// psi.begin(), ptr_fun(sqr));
+// break;
+// default:
+// Assert (false, ExcNotImplemented());
+// };
+
+// // ok, now we have the integrand,
+// // let's compute the integral,
+// // which is
+// // sum_j psi_j JxW_j
+// // (or |psi_j| or |psi_j|^2
+// switch (norm)
+// {
+// case mean:
+// case L1_norm:
+// diff = inner_product (psi.begin(), psi.end(),
+// fe_values.get_JxW_values().begin(),
+// 0.0);
+// break;
+// case L2_norm:
+// case H1_norm:
+// diff = sqrt(inner_product (psi.begin(), psi.end(),
+// fe_values.get_JxW_values().begin(),
+// 0.0));
+// break;
+// case Linfty_norm:
+// diff = *max_element (psi.begin(), psi.end());
+// break;
+// default:
+// Assert (false, ExcNotImplemented());
+// };
+
+// // note: the H1_norm uses the result
+// // of the L2_norm and control goes
+// // over to the next case statement!
+// if (norm != H1_norm)
+// break;
+// };
+
+// case H1_seminorm:
+// {
+// // note: the computation of the
+// // H1_norm starts at the previous
+// // case statement, but continues
+// // here!
+
+// // for H1_norm: re-square L2_norm.
+// diff = sqr(diff);
+
+// // same procedure as above, but now
+// // psi is a vector of gradients
+// const unsigned int n_q_points = q.n_quadrature_points;
+// vector<Tensor<1,dim> > psi (n_q_points);
+
+// // in praxi: first compute
+// // exact fe_function vector
+// exact_solution.gradient_list (fe_values.get_quadrature_points(),
+// psi);
- // then subtract finite element
- // fe_function
- if (true)
- {
- vector<Tensor<1,dim> > function_grads (n_q_points, Tensor<1,dim>());
- fe_values.get_function_grads (fe_function, function_grads);
-
- transform (psi.begin(), psi.end(),
- function_grads.begin(),
- psi.begin(),
- minus<Tensor<1,dim> >());
- };
- // take square of integrand
- vector<double> psi_square (psi.size(), 0.0);
- for (unsigned int i=0; i<n_q_points; ++i)
- psi_square[i] = sqr_point(psi[i]);
-
- // add seminorm to L_2 norm or
- // to zero
- diff += inner_product (psi_square.begin(), psi_square.end(),
- fe_values.get_JxW_values().begin(),
- 0.0);
- diff = sqrt(diff);
-
- break;
- };
+// // then subtract finite element
+// // fe_function
+// if (true)
+// {
+// vector<Tensor<1,dim> > function_grads (n_q_points, Tensor<1,dim>());
+// fe_values.get_function_grads (fe_function, function_grads);
+
+// transform (psi.begin(), psi.end(),
+// function_grads.begin(),
+// psi.begin(),
+// minus<Tensor<1,dim> >());
+// };
+// // take square of integrand
+// vector<double> psi_square (psi.size(), 0.0);
+// for (unsigned int i=0; i<n_q_points; ++i)
+// psi_square[i] = sqr_point(psi[i]);
+
+// // add seminorm to L_2 norm or
+// // to zero
+// diff += inner_product (psi_square.begin(), psi_square.end(),
+// fe_values.get_JxW_values().begin(),
+// 0.0);
+// diff = sqrt(diff);
+
+// break;
+// };
- default:
- Assert (false, ExcNotImplemented());
- };
+// default:
+// Assert (false, ExcNotImplemented());
+// };
- // append result of this cell
- // to the end of the vector
- difference(index) = diff;
- };
-};
+// // append result of this cell
+// // to the end of the vector
+// difference(index) = diff;
+// };
+// };
template VectorTools<deal_II_dimension>;
#include <numerics/base.h>
#include <numerics/assembler.h>
#include <numerics/vectors.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <map>
#include <fstream>
Equation<dim>(1),
right_hand_side (rhs) {};
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
protected:
template <>
-void PoissonEquation<2>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<2>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<2> &fe_values,
const DoFHandler<2>::cell_iterator &) const {
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void PoissonEquation<dim>::assemble (dVector &,
+void PoissonEquation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
solve ();
Solution<dim> sol;
- dVector l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell;
- dVector h1_seminorm_error_per_cell, h1_error_per_cell;
+ Vector<float> l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell;
+ Vector<float> h1_seminorm_error_per_cell, h1_error_per_cell;
cout << " Calculating L1 error... ";
VectorTools<dim>::integrate_difference (*dof_handler,
if (dof->n_dofs()<=5000)
{
- dVector l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof;
- dVector h1_seminorm_error_per_dof, h1_error_per_dof;
+ Vector<double> l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof;
+ Vector<double> h1_seminorm_error_per_dof, h1_error_per_dof;
dof->distribute_cell_to_dof_vector (l1_error_per_cell, l1_error_per_dof);
dof->distribute_cell_to_dof_vector (l2_error_per_cell, l2_error_per_dof);
dof->distribute_cell_to_dof_vector (linfty_error_per_cell,
h1_seminorm_error_per_dof);
dof->distribute_cell_to_dof_vector (h1_error_per_cell, h1_error_per_dof);
-// dVector projected_solution;
+// Vector<double> projected_solution;
// ConstraintMatrix constraints;
// constraints.close ();
// VectorTools<dim>::project (*dof, constraints, *fe,
#include <grid/tria_boundary.h>
#include <grid/tria_iterator.h>
#include <grid/tria_accessor.h>
-#include <lac/dsmatrix.h>
+#include <lac/sparsematrix.h>
#include <base/parameter_handler.h>
#include <grid/dof_constraints.h>
cout << " Renumbering degrees of freedom..." << endl;
dof->renumber_dofs (Cuthill_McKee, false);
- dSMatrixStruct sparsity (dof->n_dofs(),
- dof->max_couplings_between_dofs());
+ SparseMatrixStruct sparsity (dof->n_dofs(),
+ dof->max_couplings_between_dofs());
dof->make_sparsity_pattern (sparsity);
#include <numerics/assembler.h>
#include <numerics/vectors.h>
#include <numerics/error_estimator.h>
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <map.h>
#include <fstream.h>
right_hand_side (rhs),
coefficient (coefficient) {};
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
protected:
template <>
-void PoissonEquation<2>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<2>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<2> &fe_values,
const DoFHandler<2>::cell_iterator &) const {
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void PoissonEquation<dim>::assemble (dVector &,
+void PoissonEquation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
cout << " Solving..." << endl;
solve ();
- dVector l2_error_per_cell, linfty_error_per_cell, h1_error_per_cell;
- dVector estimated_error_per_cell;
+ Vector<float> l2_error_per_cell, linfty_error_per_cell, h1_error_per_cell;
+ Vector<float> estimated_error_per_cell;
QGauss3<dim> q;
cout << " Calculating L2 error... ";
cout << estimated_error_per_cell.l2_norm() << endl;
estimated_error.push_back (estimated_error_per_cell.l2_norm());
- dVector l2_error_per_dof, linfty_error_per_dof;
- dVector h1_error_per_dof, estimated_error_per_dof;
- dVector error_ratio;
+ Vector<double> l2_error_per_dof, linfty_error_per_dof;
+ Vector<double> h1_error_per_dof, estimated_error_per_dof;
+ Vector<double> error_ratio;
dof->distribute_cell_to_dof_vector (l2_error_per_cell, l2_error_per_dof);
dof->distribute_cell_to_dof_vector (linfty_error_per_cell,
linfty_error_per_dof);
#include <numerics/vectors.h>
#include <numerics/matrices.h>
-#include <lac/dvector.h>
-#include <lac/dsmatrix.h>
+#include <lac/vector.h>
+#include <lac/sparsematrix.h>
#include <lac/solver_cg.h>
#include <lac/vector_memory.h>
Equation<dim>(1),
right_hand_side (rhs) {};
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
protected:
/**
* Sparsity pattern of the system matrix.
*/
- dSMatrixStruct system_sparsity;
+ SparseMatrixStruct system_sparsity;
/**
* System matrix.
*/
- dSMatrix system_matrix;
+ SparseMatrix<double> system_matrix;
/**
* Vector storing the right hand side.
*/
- dVector right_hand_side;
+ Vector<double> right_hand_side;
/**
* Solution vector.
*/
- dVector solution;
+ Vector<double> solution;
/**
* List of constraints introduced by
template <>
-void PoissonEquation<2>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<2>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<2> &fe_values,
const DoFHandler<2>::cell_iterator &) const {
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void PoissonEquation<dim>::assemble (dVector &,
+void PoissonEquation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
Assert ((tria!=0) && (dof!=0), ExcNoTriaSelected());
SolverControl control(4000, 1e-16);
- PrimitiveVectorMemory<dVector> memory;
- SolverCG<dSMatrix,dVector> cg(control,memory);
+ PrimitiveVectorMemory<Vector<double> > memory;
+ SolverCG<SparseMatrix<double>,Vector<double> > cg(control,memory);
// solve
cg.solve (system_matrix, solution, right_hand_side);
solve ();
Solution<dim> sol;
- dVector l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell;
- dVector h1_seminorm_error_per_cell, h1_error_per_cell;
+ Vector<float> l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell;
+ Vector<float> h1_seminorm_error_per_cell, h1_error_per_cell;
cout << " Calculating L1 error... ";
VectorTools<dim>::integrate_difference (*dof,
if (dof->DoFHandler<dim>::n_dofs()<=5000)
{
- dVector l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof;
- dVector h1_seminorm_error_per_dof, h1_error_per_dof;
+ Vector<double> l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof;
+ Vector<double> h1_seminorm_error_per_dof, h1_error_per_dof;
dof->distribute_cell_to_dof_vector (l1_error_per_cell, l1_error_per_dof);
dof->distribute_cell_to_dof_vector (l2_error_per_cell, l2_error_per_dof);
dof->distribute_cell_to_dof_vector (linfty_error_per_cell,
h1_seminorm_error_per_dof);
dof->distribute_cell_to_dof_vector (h1_error_per_cell, h1_error_per_dof);
-// dVector projected_solution;
+// Vector<double> projected_solution;
// ConstraintMatrix constraints;
// constraints.close ();
// VectorTools<dim>::project (*dof, constraints, *fe,
#include "poisson.h"
-#include <lac/dvector.h>
+#include <lac/vector.h>
#if deal_II_dimension == 1
template <>
-void PoissonEquation<1>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<1>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<1> &fe_values,
const DoFHandler<1>::cell_iterator &) const {
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
//#if deal_II_dimension >= 2
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &) const {
const vector<vector<Tensor<1,dim> > >&gradients = fe_values.get_shape_grads ();
- const dFMatrix &values = fe_values.get_shape_values ();
+ const FullMatrix<double> &values = fe_values.get_shape_values ();
vector<double> rhs_values (fe_values.n_quadrature_points);
const vector<double> &weights = fe_values.get_JxW_values ();
template <int dim>
-void PoissonEquation<dim>::assemble (dFMatrix &,
+void PoissonEquation<dim>::assemble (FullMatrix<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
template <int dim>
-void PoissonEquation<dim>::assemble (dVector &,
+void PoissonEquation<dim>::assemble (Vector<double> &,
const FEValues<dim> &,
const DoFHandler<dim>::cell_iterator &) const {
Assert (false, ExcPureVirtualFunctionCalled());
#include "poisson.h"
-#include <lac/dvector.h>
+#include <lac/vector.h>
#include <base/quadrature_lib.h>
#include <numerics/base.h>
#include <numerics/assembler.h>
-#include <lac/dsmatrix.h>
+#include <lac/sparsematrix.h>
#include <map.h>
#include <fstream.h>
#include <cmath>
-extern "C" {
-# include <stdlib.h>
-}
+#include <cstdlib>
+
Equation<dim>(1),
right_hand_side (rhs) {};
- virtual void assemble (dFMatrix &cell_matrix,
- dVector &rhs,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
+ Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dFMatrix &cell_matrix,
+ virtual void assemble (FullMatrix<double> &cell_matrix,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
- virtual void assemble (dVector &rhs,
+ virtual void assemble (Vector<double> &rhs,
const FEValues<dim> &fe_values,
const DoFHandler<dim>::cell_iterator &cell) const;
protected:
#include "poisson.h"
-#include <lac/dvector.h>
+#include <lac/vector.h>