]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add a function make_sparsity_pattern which uses a mask; this is useful only for syste...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 20 Mar 1999 14:59:44 +0000 (14:59 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 20 Mar 1999 14:59:44 +0000 (14:59 +0000)
git-svn-id: https://svn.dealii.org/trunk@1031 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_handler.h
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/source/dofs/dof_handler.cc

index 626b52476a073342a582fd8a778d5f38ee352634..55f5a0d5f549788dac4ea447e5e38d40eb720af2 100644 (file)
@@ -633,6 +633,50 @@ class DoFHandler : public DoFDimensionInfo<dim> {
                                      */
     void make_sparsity_pattern (SparseMatrixStruct &) const; 
 
+                                    /**
+                                     * This function does mistly the same as
+                                     * the above one, but it is specialized for
+                                     * mixed finite elements and allows to
+                                     * specify which variables couple in which
+                                     * equation. For example, if wanted to solve
+                                     * the Stokes equations,
+                                     * \begin{verbatim}
+                                     *    -\Delta \vec u + \nabla p = 0,
+                                     *    \div u                    = 0
+                                     * \end{verbatim}
+                                     * in, say, two space dimensions, using
+                                     * nonstabilized Q2/Q2/Q1 mixed elements
+                                     * (i.e. using the #FESystem# class), then
+                                     * you don't want all degrees of freedom
+                                     * to couple in each equation. You rather
+                                     * may want to give the following pattern
+                                     * of couplings:
+                                     * \begin{verbatim}
+                                     *    1  0  1
+                                     *    0  1  1
+                                     *    1  1  0
+                                     * \end{verbatim}
+                                     * where "1" indicates that two variables
+                                     * (i.e. components of the #FESystem#)
+                                     * couple in the respective equation, and
+                                     * a "0" means no coupling, in which case
+                                     * it is not necessary to allocate space
+                                     * in the matrix structure. Obviously, the
+                                     * mask refers to components of the
+                                     * composed #FESystem#, rather than to the
+                                     * degrees of freedom contained in there.
+                                     *
+                                     * This function is designed to accept
+                                     * a mask, like the one shown above,
+                                     * through the #mask# parameter, which
+                                     * contains boolean values. It builds
+                                     * the matrix structure just like the
+                                     * previous function, but does not create
+                                     * elements if not specified by the mask.
+                                     */
+    void make_sparsity_pattern (const vector<vector<bool> > &mask,
+                               SparseMatrixStruct          &) const;
+    
                                     /**
                                      * Write the sparsity structure of the
                                      * matrix composed of the basis functions
@@ -1499,7 +1543,14 @@ class DoFHandler : public DoFDimensionInfo<dim> {
                    int,
                    << "The given list of new dof indices is not consecutive: "
                    << "the index " << arg1 << " does not exist.");
-
+                                    /**
+                                     * Exception
+                                     */
+    DeclException2 (ExcInvalidMaskDimension,
+                   int, int,
+                   << "The dimension of the mask " << arg1 << " does not match"
+                   << " the number of components in the finite element object.");
+    
   protected:
     
                                     /**
index 459711ffe970852f27b4d72d933689aaa2ae52e2..2e463e4aa62ded0b476690bcc38dad7105040073 100644 (file)
@@ -250,7 +250,7 @@ class FiniteElementBase :
                                      * Compute system index from components.
                                      */
     unsigned int component_to_system_index (unsigned int component,
-                                       unsigned int component_index) const;
+                                           unsigned int component_index) const;
   
                                     /**
                                      * Compute component and index from system index.
index 074bd35434588e0d9181ebb72a8f8a6296e2ef7b..b515b8494f613b60caf39e83308cf273642fdec9 100644 (file)
@@ -2202,6 +2202,55 @@ void DoFHandler<dim>::make_sparsity_pattern (SparseMatrixStruct &sparsity) const
 };
 
 
+
+template <int dim>
+void
+DoFHandler<dim>::make_sparsity_pattern (const vector<vector<bool> > &mask,
+                                       SparseMatrixStruct          &sparsity) const
+{
+  const unsigned int total_dofs = selected_fe->total_dofs;
+
+  Assert (selected_fe != 0, ExcNoFESelected());
+  Assert (sparsity.n_rows() == n_dofs(),
+         ExcDifferentDimensions (sparsity.n_rows(), n_dofs()));
+  Assert (sparsity.n_cols() == n_dofs(),
+         ExcDifferentDimensions (sparsity.n_cols(), n_dofs()));
+  Assert (mask.size() == selected_fe->n_components,
+         ExcInvalidMaskDimension(mask.size(), selected_fe->n_components));
+  for (unsigned int i=0; i<mask.size(); ++i)
+    Assert (mask[i].size() == selected_fe->n_components,
+           ExcInvalidMaskDimension(mask[i].size(), selected_fe->n_components));
+
+                                  // first build a mask for each dof,
+                                  // not like the one given which represents
+                                  // components
+  vector<vector<bool> > dof_mask(total_dofs,
+                                vector<bool>(total_dofs, false));
+  for (unsigned int i=0; i<total_dofs; ++i)
+    for (unsigned int j=0; j<total_dofs; ++j)
+      dof_mask[i][j] = mask
+                      [selected_fe->system_to_component_index(i).first]
+                      [selected_fe->system_to_component_index(j).first];
+  
+  
+  vector<int> dofs_on_this_cell(total_dofs);
+  active_cell_iterator cell = begin_active(),
+                      endc = end();
+  for (; cell!=endc; ++cell) 
+    {
+      cell->get_dof_indices (dofs_on_this_cell);
+                                      // make sparsity pattern for this cell
+      for (unsigned int i=0; i<total_dofs; ++i)
+       for (unsigned int j=0; j<total_dofs; ++j)
+         if (dof_mask[i][j] == true)
+           sparsity.add (dofs_on_this_cell[i],
+                         dofs_on_this_cell[j]);
+    };
+};
+
+
+
+
 #if deal_II_dimension == 1
 
 template <>

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.