From d4aa4a761d17bb55ca533ef3ba0ceef97c407a52 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 18 Mar 1998 15:41:11 +0000 Subject: [PATCH] Cleanup and small changes. git-svn-id: https://svn.dealii.org/trunk@78 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/numerics/assembler.h | 23 +++------ deal.II/deal.II/include/numerics/base.h | 17 +++--- deal.II/deal.II/source/numerics/assembler.cc | 54 ++++++++------------ deal.II/deal.II/source/numerics/base.cc | 31 ++++------- 4 files changed, 47 insertions(+), 78 deletions(-) diff --git a/deal.II/deal.II/include/numerics/assembler.h b/deal.II/deal.II/include/numerics/assembler.h index a7d2e036c1..bd2192629f 100644 --- a/deal.II/deal.II/include/numerics/assembler.h +++ b/deal.II/deal.II/include/numerics/assembler.h @@ -24,7 +24,7 @@ template class ProblemBase; finite element discretisation of one or more equations. Equation objects need only provide functions which set up the cell - matrices and the right hand side(s). These are then automatically inserted + matrices and the cell right hand side. These are then automatically inserted into the global matrices and vectors. */ template @@ -39,12 +39,11 @@ class Equation { /** * Virtual function which assembles the - * cell matrix and a specific (user - * selectable) number of right hand sides + * cell matrix and the right hand side * on a given cell. */ virtual void assemble (dFMatrix &cell_matrix, - vector &rhs, + dVector &rhs, const FEValues &fe_values, const Triangulation::cell_iterator &cell) const; @@ -58,9 +57,9 @@ class Equation { /** * Virtual function which only assembles - * the right hand side(s) on a given cell. + * the right hand side on a given cell. */ - virtual void assemble (vector &rhs, + virtual void assemble (dVector &rhs, const FEValues &fe_values, const Triangulation::cell_iterator &cell) const; @@ -104,7 +103,6 @@ struct AssemblerData { AssemblerData (const DoFHandler &dof, const bool assemble_matrix, const bool assemble_rhs, - const unsigned int n_rhs, ProblemBase &problem, const Quadrature &quadrature, const FiniteElement &fe); @@ -119,10 +117,6 @@ struct AssemblerData { * and right hand sides. */ const bool assemble_matrix, assemble_rhs; - /** - * How many right hand sides are there. - */ - const unsigned int n_rhs; /** * Pointer to the #ProblemBase# object * of which the global matrices and @@ -198,10 +192,9 @@ class Assembler : public DoFCellAccessor { dFMatrix cell_matrix; /** - * Have room for a certain number of - * right hand sides. + * Right hand side local to cell. */ - vector cell_vectors; + dVector cell_vector; /** * Store whether to assemble the @@ -211,7 +204,7 @@ class Assembler : public DoFCellAccessor { /** * Store whether to assemble the - * right hand sides. + * right hand side. */ bool assemble_rhs; diff --git a/deal.II/deal.II/include/numerics/base.h b/deal.II/deal.II/include/numerics/base.h index 0a98cc5b5c..2a5f1e94f4 100644 --- a/deal.II/deal.II/include/numerics/base.h +++ b/deal.II/deal.II/include/numerics/base.h @@ -43,8 +43,8 @@ template class DataOut; \item Loop over all cells and assemble matrix and vectors using the given quadrature formula and the equation object which contains the weak formulation of the equation. - \item Condense the system matrix with the constraints induced by hanging - nodes. + \item Condense the system matrix and right hand side with the constraints + induced by hanging nodes. \end{itemize} */ template @@ -55,12 +55,10 @@ class ProblemBase { * degree of freedom object during the * lifetime of this object. The dof * object must refer to the given - * triangulation. The number of right hand - * sides being assembled defaults to one. + * triangulation. */ ProblemBase (Triangulation *tria, - DoFHandler *dof_handler, - const unsigned int n_rhs=1); + DoFHandler *dof_handler); /** * Initiate the process of assemblage of @@ -111,12 +109,9 @@ class ProblemBase { dSMatrix system_matrix; /** - * Pointers to the right hand sides to the - * problem. Usually, one only needs one - * right hand side, but simetimes one - * wants more. + * Vector storing the right hand side. */ - vector right_hand_sides; + dVector right_hand_side; /** * Solution vector. diff --git a/deal.II/deal.II/source/numerics/assembler.cc b/deal.II/deal.II/source/numerics/assembler.cc index 3b99a37c69..0153753146 100644 --- a/deal.II/deal.II/source/numerics/assembler.cc +++ b/deal.II/deal.II/source/numerics/assembler.cc @@ -19,7 +19,7 @@ Equation::Equation (const unsigned int n_equations) : void Equation<1>::assemble (dFMatrix &, - vector &, + dVector &, const FEValues<1> &, const Triangulation<1>::cell_iterator &) const { Assert (false, ExcPureVirtualFunctionCalled()); @@ -28,7 +28,7 @@ void Equation<1>::assemble (dFMatrix &, void Equation<2>::assemble (dFMatrix &, - vector &, + dVector &, const FEValues<2> &, const Triangulation<2>::cell_iterator &) const { Assert (false, ExcPureVirtualFunctionCalled()); @@ -52,7 +52,7 @@ void Equation<2>::assemble (dFMatrix &, -void Equation<1>::assemble (vector &, +void Equation<1>::assemble (dVector &, const FEValues<1> &, const Triangulation<1>::cell_iterator &) const { Assert (false, ExcPureVirtualFunctionCalled()); @@ -60,7 +60,7 @@ void Equation<1>::assemble (vector &, -void Equation<2>::assemble (vector &, +void Equation<2>::assemble (dVector &, const FEValues<2> &, const Triangulation<2>::cell_iterator &) const { Assert (false, ExcPureVirtualFunctionCalled()); @@ -73,13 +73,14 @@ template AssemblerData::AssemblerData (const DoFHandler &dof, const bool assemble_matrix, const bool assemble_rhs, - const unsigned int n_rhs, ProblemBase &problem, const Quadrature &quadrature, const FiniteElement &fe) : - dof(dof), assemble_matrix(assemble_matrix), - assemble_rhs(assemble_rhs), n_rhs(n_rhs), - problem(problem), quadrature(quadrature), + dof(dof), + assemble_matrix(assemble_matrix), + assemble_rhs(assemble_rhs), + problem(problem), + quadrature(quadrature), fe(fe) {}; @@ -93,8 +94,7 @@ Assembler::Assembler (Triangulation *tria, DoFCellAccessor (tria,level,index, &((AssemblerData*)local_data)->dof), cell_matrix (dof_handler->get_selected_fe().total_dofs), - cell_vectors (((AssemblerData*)local_data)->n_rhs, - dVector(dof_handler->get_selected_fe().total_dofs)), + 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), @@ -102,21 +102,14 @@ Assembler::Assembler (Triangulation *tria, fe_values (((AssemblerData*)local_data)->fe, ((AssemblerData*)local_data)->quadrature) { - Assert (((AssemblerData*)local_data)->n_rhs != 0, - ExcInvalidData()); 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 (problem.right_hand_sides.size() == cell_vectors.size(), - ExcInvalidData()); Assert (((AssemblerData*)local_data)->fe == dof_handler->get_selected_fe(), ExcInvalidData()); - - - for (unsigned int i=0; in() == dof_handler->n_dofs(), - ExcInvalidData()); + Assert ((unsigned int)problem.right_hand_side.n() == dof_handler->n_dofs(), + ExcInvalidData()); }; @@ -128,24 +121,22 @@ void Assembler::assemble (const Equation &equation) { present_level, present_index), fe); - + + // clear cell matrix if (assemble_matrix) - // clear cell matrix for (unsigned int i=0; iget_selected_fe().total_dofs; ++i) for (unsigned int j=0; jget_selected_fe().total_dofs; ++j) cell_matrix(i,j) = 0; + // clear cell vector if (assemble_rhs) - // clear cell vector - for (unsigned int vec=0; vecget_selected_fe().total_dofs; ++j) - cell_vectors[vec](j) = 0; + cell_vector.clear (); // fill cell matrix and vector if required if (assemble_matrix && assemble_rhs) - equation.assemble (cell_matrix, cell_vectors, fe_values, + equation.assemble (cell_matrix, cell_vector, fe_values, Triangulation::cell_iterator(tria, present_level, present_index)); @@ -157,7 +148,7 @@ void Assembler::assemble (const Equation &equation) { present_index)); else if (assemble_rhs) - equation.assemble (cell_vectors, fe_values, + equation.assemble (cell_vector, fe_values, Triangulation::cell_iterator(tria, present_level, present_index)); @@ -169,17 +160,16 @@ void Assembler::assemble (const Equation &equation) { vector dofs; dof_indices (dofs); + // distribute cell matrix if (assemble_matrix) - // distribute cell matrix for (unsigned int i=0; iget_selected_fe().total_dofs; ++i) for (unsigned int j=0; jget_selected_fe().total_dofs; ++j) problem.system_matrix.add(dofs[i], dofs[j], cell_matrix(i,j)); + // distribute cell vector if (assemble_rhs) - // distribute cell vector - for (unsigned int vec=0; vecget_selected_fe().total_dofs; ++j) - (*problem.right_hand_sides[vec])(dofs[j]) += cell_vectors[vec](j); + for (unsigned int j=0; jget_selected_fe().total_dofs; ++j) + problem.right_hand_side(dofs[j]) += cell_vector(j); }; diff --git a/deal.II/deal.II/source/numerics/base.cc b/deal.II/deal.II/source/numerics/base.cc index 722129f720..106fdddfaf 100644 --- a/deal.II/deal.II/source/numerics/base.cc +++ b/deal.II/deal.II/source/numerics/base.cc @@ -6,11 +6,12 @@ #include #include - -#include "../../../mia/vectormemory.h" #include "../../../mia/control.h" +#include "../../../mia/vectormemory.h" #include "../../../mia/cg.h" + + extern TriaActiveIterator<1,CellAccessor<1> > __dummy1233; // for gcc2.7 extern TriaActiveIterator<2,CellAccessor<2> > __dummy1234; @@ -18,21 +19,13 @@ extern TriaActiveIterator<2,CellAccessor<2> > __dummy1234; template ProblemBase::ProblemBase (Triangulation *tria, - DoFHandler *dof, - const unsigned int n_rhs) : + DoFHandler *dof) : tria(tria), dof_handler(dof), system_sparsity(1,1,1), // dummy initialisation, is later reinit'd - system_matrix(), // dummy initialisation, is later reinit'd - right_hand_sides (n_rhs, (dVector*)0) + system_matrix() // dummy initialisation, is later reinit'd { Assert (tria == &dof->get_tria(), ExcDofAndTriaDontMatch()); - - for (unsigned int i=0; i::assemble (const Equation &equation, system_sparsity.reinit (dof_handler->n_dofs(), dof_handler->n_dofs(), dof_handler->max_couplings_between_dofs()); - for (unsigned int i=0; ireinit (dof_handler->n_dofs()); + right_hand_side.reinit (dof_handler->n_dofs()); // make up sparsity pattern and // compress with constraints @@ -62,7 +54,6 @@ void ProblemBase::assemble (const Equation &equation, // create assembler AssemblerData data (*dof_handler, true, true, //assemble matrix and rhs - right_hand_sides.size(), *this, quadrature, fe); @@ -80,9 +71,9 @@ void ProblemBase::assemble (const Equation &equation, // condense system matrix in-place constraints.condense (system_matrix); - // condense right hand sides in-place - for (unsigned int i=0; i::solve () { double tolerance = 1.e-16; Control control1(max_iter,tolerance); - PrimitiveVectorMemory memory(right_hand_sides[0]->n()); + PrimitiveVectorMemory memory(right_hand_side.n()); CG cg(control1,memory); // solve - cg (system_matrix, solution, *right_hand_sides[0]); + cg (system_matrix, solution, right_hand_side); // distribute solution constraints.distribute (solution); }; -- 2.39.5