From f4f23563dd9db95c1f4df161d6e01c29bf94f0d7 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 14 Apr 1998 14:47:32 +0000 Subject: [PATCH] Make the assembler class more flexible (why didn't I think of that beforehand? The old way it was pretty unusable...) git-svn-id: https://svn.dealii.org/trunk@171 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/numerics/assembler.h | 35 ++++++++++++++------ deal.II/deal.II/source/numerics/assembler.cc | 30 ++++++++--------- deal.II/deal.II/source/numerics/base.cc | 3 +- 3 files changed, 42 insertions(+), 26 deletions(-) diff --git a/deal.II/deal.II/include/numerics/assembler.h b/deal.II/deal.II/include/numerics/assembler.h index 44af7879eb..2bfc408a1b 100644 --- a/deal.II/deal.II/include/numerics/assembler.h +++ b/deal.II/deal.II/include/numerics/assembler.h @@ -14,7 +14,6 @@ // forward declarations class dFMatrix; class dVector; -template class ProblemBase; @@ -103,7 +102,8 @@ struct AssemblerData { AssemblerData (const DoFHandler &dof, const bool assemble_matrix, const bool assemble_rhs, - ProblemBase &problem, + dSMatrix &matrix, + dVector &rhs_vector, const Quadrature &quadrature, const FiniteElement &fe, const FEValues::UpdateStruct &update_flags); @@ -113,22 +113,31 @@ struct AssemblerData { * to be used to iterate on. */ const DoFHandler &dof; + /** * Flags whether to assemble the matrix * and right hand sides. */ const bool assemble_matrix, assemble_rhs; + + /** + * Pointer to the matrix to be assembled + * by this object. + */ + dSMatrix &matrix; + /** - * Pointer to the #ProblemBase# object - * of which the global matrices and - * vectors are to be assembled. + * Pointer to the vector to be assembled + * by this object. */ - ProblemBase &problem; + dVector &rhs_vector; + /** * Pointer to a quadrature object to be * used for this assemblage process. */ const Quadrature &quadrature; + /** * Use this FE type for the assemblage * process. The FE object must match that @@ -136,6 +145,7 @@ struct AssemblerData { * object. */ const FiniteElement &fe; + /** * Store which of the fields of the * FEValues object need to be reinitialized @@ -216,11 +226,16 @@ class Assembler : public DoFCellAccessor { bool assemble_rhs; /** - * Pointer to the problem class object - * of which we are to use the matrices - * and vectors. + * Pointer to the matrix to be assembled + * by this object. + */ + dSMatrix &matrix; + + /** + * Pointer to the vector to be assembled + * by this object. */ - ProblemBase &problem; + dVector &rhs_vector; /** * Pointer to the finite element used for diff --git a/deal.II/deal.II/source/numerics/assembler.cc b/deal.II/deal.II/source/numerics/assembler.cc index 591e2de10d..66411ddf05 100644 --- a/deal.II/deal.II/source/numerics/assembler.cc +++ b/deal.II/deal.II/source/numerics/assembler.cc @@ -1,12 +1,11 @@ /* $Id$ */ #include -#include #include #include #include #include - +#include template @@ -45,17 +44,19 @@ void Equation::assemble (dVector &, template -AssemblerData::AssemblerData (const DoFHandler &dof, - const bool assemble_matrix, - const bool assemble_rhs, - ProblemBase &problem, +AssemblerData::AssemblerData (const DoFHandler &dof, + const bool assemble_matrix, + const bool assemble_rhs, + dSMatrix &matrix, + dVector &rhs_vector, const Quadrature &quadrature, const FiniteElement &fe, const FEValues::UpdateStruct &update_flags) : dof(dof), assemble_matrix(assemble_matrix), assemble_rhs(assemble_rhs), - problem(problem), + matrix(matrix), + rhs_vector(rhs_vector), quadrature(quadrature), fe(fe), update_flags(update_flags) {}; @@ -74,19 +75,18 @@ Assembler::Assembler (Triangulation *tria, cell_vector (dVector(dof_handler->get_selected_fe().total_dofs)), assemble_matrix (((AssemblerData*)local_data)->assemble_matrix), assemble_rhs (((AssemblerData*)local_data)->assemble_rhs), - problem (((AssemblerData*)local_data)->problem), + matrix(((AssemblerData*)local_data)->matrix), + rhs_vector(((AssemblerData*)local_data)->rhs_vector), fe(((AssemblerData*)local_data)->fe), fe_values (((AssemblerData*)local_data)->fe, ((AssemblerData*)local_data)->quadrature, ((AssemblerData*)local_data)->update_flags) { - Assert ((unsigned int)problem.system_matrix.m() == dof_handler->n_dofs(), - ExcInvalidData()); - Assert ((unsigned int)problem.system_matrix.n() == dof_handler->n_dofs(), - ExcInvalidData()); + Assert (matrix.m() == dof_handler->n_dofs(), ExcInvalidData()); + Assert (matrix.n() == dof_handler->n_dofs(), ExcInvalidData()); Assert (((AssemblerData*)local_data)->fe == dof_handler->get_selected_fe(), ExcInvalidData()); - Assert ((unsigned int)problem.right_hand_side.n() == dof_handler->n_dofs(), + Assert (rhs_vector.n() == dof_handler->n_dofs(), ExcInvalidData()); }; @@ -143,12 +143,12 @@ void Assembler::assemble (const Equation &equation) { if (assemble_matrix) for (unsigned int i=0; i::assemble (const Equation &equation, // create assembler AssemblerData data (*dof_handler, true, true, //assemble matrix and rhs - *this, + system_matrix, + right_hand_side, quadrature, fe, update_flags); -- 2.39.5