]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make the assembler class more flexible (why didn't I think of that beforehand? The...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 14 Apr 1998 14:47:32 +0000 (14:47 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 14 Apr 1998 14:47:32 +0000 (14:47 +0000)
git-svn-id: https://svn.dealii.org/trunk@171 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/assembler.h
deal.II/deal.II/source/numerics/assembler.cc
deal.II/deal.II/source/numerics/base.cc

index 44af7879eb023bec15cccfac147573852ba2551a..2bfc408a1b4ae35fe1490b1806c896220bac0829 100644 (file)
@@ -14,7 +14,6 @@
 // forward declarations
 class dFMatrix;
 class dVector;
-template <int dim> class ProblemBase;
 
 
 
@@ -103,7 +102,8 @@ struct AssemblerData {
     AssemblerData (const DoFHandler<dim>    &dof,
                   const bool                assemble_matrix,
                   const bool                assemble_rhs,
-                  ProblemBase<dim>         &problem,
+                  dSMatrix                 &matrix,
+                  dVector                  &rhs_vector,
                   const Quadrature<dim>    &quadrature,
                   const FiniteElement<dim> &fe,
                   const FEValues<dim>::UpdateStruct &update_flags);
@@ -113,22 +113,31 @@ struct AssemblerData {
                                      * to be used to iterate on.
                                      */
     const DoFHandler<dim>  &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<dim>       &problem;
+    dVector                &rhs_vector;
+    
                                     /**
                                      * Pointer to a quadrature object to be
                                      * used for this assemblage process.
                                      */
     const Quadrature<dim>  &quadrature;
+    
                                     /**
                                      * Use this FE type for the assemblage
                                      * process. The FE object must match that
@@ -136,6 +145,7 @@ struct AssemblerData {
                                      * object.
                                      */
     const FiniteElement<dim> &fe;
+    
                                     /**
                                      * Store which of the fields of the
                                      * FEValues object need to be reinitialized
@@ -216,11 +226,16 @@ class Assembler : public DoFCellAccessor<dim> {
     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<dim> &problem;
+    dVector                &rhs_vector;
 
                                     /**
                                      * Pointer to the finite element used for
index 591e2de10d1f3c8ab519d4ef53d439e1180bc0d8..66411ddf05526ef79c181b0d3e917a184d57bc77 100644 (file)
@@ -1,12 +1,11 @@
 /* $Id$ */
 
 #include <numerics/assembler.h>
-#include <numerics/base.h>
 #include <grid/tria_iterator.h>
 #include <grid/tria_iterator.templates.h>
 #include <lac/dfmatrix.h>
 #include <lac/dvector.h>
-
+#include <lac/dsmatrix.h>
 
 
 template <int dim>
@@ -45,17 +44,19 @@ void Equation<dim>::assemble (dVector           &,
 
 
 template <int dim>
-AssemblerData<dim>::AssemblerData (const DoFHandler<dim>   &dof,
-                                  const bool               assemble_matrix,
-                                  const bool               assemble_rhs,
-                                  ProblemBase<dim>         &problem,
+AssemblerData<dim>::AssemblerData (const DoFHandler<dim>    &dof,
+                                  const bool                assemble_matrix,
+                                  const bool                assemble_rhs,
+                                  dSMatrix                 &matrix,
+                                  dVector                  &rhs_vector,
                                   const Quadrature<dim>    &quadrature,
                                   const FiniteElement<dim> &fe,
                                   const FEValues<dim>::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<dim>::Assembler (Triangulation<dim> *tria,
                cell_vector (dVector(dof_handler->get_selected_fe().total_dofs)),
                assemble_matrix (((AssemblerData<dim>*)local_data)->assemble_matrix),
                assemble_rhs (((AssemblerData<dim>*)local_data)->assemble_rhs),
-               problem (((AssemblerData<dim>*)local_data)->problem),
+               matrix(((AssemblerData<dim>*)local_data)->matrix),
+               rhs_vector(((AssemblerData<dim>*)local_data)->rhs_vector),
                fe(((AssemblerData<dim>*)local_data)->fe),
                fe_values (((AssemblerData<dim>*)local_data)->fe,
                           ((AssemblerData<dim>*)local_data)->quadrature,
                           ((AssemblerData<dim>*)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<dim>*)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<dim>::assemble (const Equation<dim> &equation) {
   if (assemble_matrix)
     for (unsigned int i=0; i<n_dofs; ++i)
       for (unsigned int j=0; j<n_dofs; ++j)
-       problem.system_matrix.add(dofs[i], dofs[j], cell_matrix(i,j));
+       matrix.add(dofs[i], dofs[j], cell_matrix(i,j));
 
                                   // distribute cell vector
   if (assemble_rhs)
     for (unsigned int j=0; j<n_dofs; ++j)
-      problem.right_hand_side(dofs[j]) += cell_vector(j);
+      rhs_vector(dofs[j]) += cell_vector(j);
 };
 
                
index b96bb8ba933abc6d7563319a917c62bd8c04d605..d672c0c4d0058b1696176befe52f216e6b2cbd51 100644 (file)
@@ -100,7 +100,8 @@ void ProblemBase<dim>::assemble (const Equation<dim>               &equation,
                                   // create assembler
   AssemblerData<dim> data (*dof_handler,
                           true, true, //assemble matrix and rhs
-                          *this,
+                          system_matrix,
+                          right_hand_side,
                           quadrature,
                           fe,
                           update_flags);

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

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.