]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add function add_entries_local_to_global for rectangular matrices.
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 25 Aug 2010 14:08:55 +0000 (14:08 +0000)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 25 Aug 2010 14:08:55 +0000 (14:08 +0000)
git-svn-id: https://svn.dealii.org/trunk@21704 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/constraint_matrix.h
deal.II/lac/include/lac/constraint_matrix.templates.h
deal.II/lac/source/constraint_matrix.cc

index b3296ce25c02a018dcfa3a87831b8fd5bddf1776..a175c86a65627c407b81eeaecd8b30a066f6377d 100644 (file)
@@ -1359,9 +1359,9 @@ class ConstraintMatrix : public Subscriptor
                                 MatrixType                      &global_matrix) const;
 
                                      /**
-                                      * Does the same as the function above
-                                      * but can treat
-                                      * non quadratic matrices.
+                                      * Does the same as the function
+                                      * above but can treat non
+                                      * quadratic matrices.
                                       */
     template <typename MatrixType>
     void
@@ -1474,6 +1474,20 @@ class ConstraintMatrix : public Subscriptor
     template <typename SparsityType>
     void
     add_entries_local_to_global (const std::vector<unsigned int> &local_dof_indices,
+                                SparsityType                    &sparsity_pattern,
+                                const bool                       keep_constrained_entries = true,
+                                const Table<2,bool>             &dof_mask = default_empty_table) const;
+
+                                    /**
+                                     * Similar to the other function,
+                                     * but for non-quadratic sparsity
+                                     * patterns.
+                                     */
+
+    template <typename SparsityType>
+    void
+    add_entries_local_to_global (const std::vector<unsigned int> &row_indices,
+                                const std::vector<unsigned int> &col_indices,
                                 SparsityType                    &sparsity_pattern,
                                 const bool                       keep_constrained_entries = true,
                                 const Table<2,bool>             &dof_mask = default_empty_table) const;
@@ -1894,12 +1908,27 @@ class ConstraintMatrix : public Subscriptor
                                      * function.
                                      *
                                      * Creates a list of affected
-                                     * rows for distribution.
+                                     * global rows for distribution,
+                                     * including the local rows where
+                                     * the entries come from.
                                      */
     void
     make_sorted_row_list (const std::vector<unsigned int> &local_dof_indices,
                          internals::GlobalRowsFromLocal  &global_rows) const;
 
+                                    /**
+                                     * Internal helper function for
+                                     * add_entries_local_to_global
+                                     * function.
+                                     *
+                                     * Creates a list of affected
+                                     * rows for distribution without
+                                     * any additional information.
+                                     */
+    void
+    make_sorted_row_list (const std::vector<unsigned int> &local_dof_indices,
+                         std::vector<unsigned int>       &active_dofs) const;
+
                                     /**
                                      * Internal helper function for
                                      * distribute_local_to_global function.
index 0b6aae5a0ab0aaba6481536c6e00f388c518395a..980933b62890115091fe8af9835bcc6d4bc3ee76 100644 (file)
@@ -2119,6 +2119,54 @@ make_sorted_row_list (const std::vector<unsigned int> &local_dof_indices,
 
 
 
+inline
+void
+ConstraintMatrix::
+make_sorted_row_list (const std::vector<unsigned int> &local_dof_indices,
+                     std::vector<unsigned int>       &active_dofs) const
+{
+  const unsigned int n_local_dofs = local_dof_indices.size();
+  unsigned int added_rows = 0;
+  for (unsigned int i = 0; i<n_local_dofs; ++i)
+    {
+      if (is_constrained(local_dof_indices[i]) == false)
+       {
+         active_dofs[added_rows++] = local_dof_indices[i];
+         continue;
+       }
+
+      active_dofs[n_local_dofs-i+added_rows-1] = i;
+    }
+  std::sort (active_dofs.begin(), active_dofs.begin()+added_rows);
+
+  const unsigned int n_constrained_dofs = n_local_dofs-added_rows;
+  for (unsigned int i=n_constrained_dofs; i>0; --i)
+    {
+      const unsigned int local_row = active_dofs.back();
+      active_dofs.pop_back();
+      const unsigned int global_row = local_dof_indices[local_row];
+      const ConstraintLine & position =
+       lines[lines_cache[calculate_line_index(global_row)]];
+      for (unsigned int q=0; q<position.entries.size(); ++q)
+       {
+         const unsigned int new_index = position.entries[q].first;
+         if (active_dofs[active_dofs.size()-i] < new_index)
+           active_dofs.insert(active_dofs.end()-i+1,new_index);
+         else
+           {
+             std::vector<unsigned int>::iterator it =
+               std::lower_bound(active_dofs.begin(),
+                                active_dofs.end()-i+1,
+                                new_index);
+             if (*it != new_index)
+               active_dofs.insert(it, new_index);
+           }
+       }
+    }
+}
+
+
+
 template <typename SparsityType>
 inline
 void
@@ -2598,6 +2646,56 @@ add_entries_local_to_global (const std::vector<unsigned int> &local_dof_indices,
 
 
 
+template <typename SparsityType>
+void
+ConstraintMatrix::
+add_entries_local_to_global (const std::vector<unsigned int> &row_indices,
+                            const std::vector<unsigned int> &col_indices,
+                            SparsityType                    &sparsity_pattern,
+                            const bool                       keep_constrained_entries,
+                            const Table<2,bool>             &dof_mask) const
+{
+  const unsigned int n_local_rows = row_indices.size();
+  const unsigned int n_local_cols = col_indices.size();
+  bool dof_mask_is_active = false;
+  Assert (keep_constrained_entries == false, ExcNotImplemented());
+  if (dof_mask.n_rows() == n_local_rows && dof_mask.n_cols() == n_local_cols)
+    dof_mask_is_active = true;
+
+                                  // if the dof mask is not active, all we
+                                  // have to do is to add some indices in a
+                                  // matrix format. To do this, we first
+                                  // create an array of all the indices
+                                  // that are to be added. these indices
+                                  // are the local dof indices plus some
+                                  // indices that come from constraints.
+  if (dof_mask_is_active == false)
+    {
+      std::vector<unsigned int> actual_row_indices (n_local_rows);
+      std::vector<unsigned int> actual_col_indices (n_local_cols);
+      make_sorted_row_list (row_indices, actual_row_indices);
+      make_sorted_row_list (col_indices, actual_col_indices);
+      const unsigned int n_actual_rows = actual_row_indices.size();
+
+                                      // now add the indices we collected above
+                                      // to the sparsity pattern. Very easy
+                                      // here - just add the same array to all
+                                      // the rows...
+      for (unsigned int i=0; i<n_actual_rows; ++i)
+       sparsity_pattern.add_entries(actual_row_indices[i],
+                                    actual_col_indices.begin(),
+                                    actual_col_indices.end(),
+                                    true);
+      return;
+    }
+
+                               // TODO: implement this
+  Assert (false, ExcNotImplemented());
+}
+
+
+
+
 template <typename SparsityType>
 void
 ConstraintMatrix::
index 040a854e4dde31c65f41be1db34581eba9e06c63..689e63be6c5a6438503590b546496b09e901f702 100644 (file)
@@ -2511,7 +2511,13 @@ BLOCK_MATRIX_VECTOR_FUNCTIONS(TrilinosWrappers::BlockSparseMatrix, TrilinosWrapp
     SparsityType &,                    \
     const bool,                        \
     const Table<2,bool> &, \
-    internal::bool2type<false>) const
+    internal::bool2type<false>) const; \
+  template void ConstraintMatrix::add_entries_local_to_global<SparsityType> (\
+    const std::vector<unsigned int> &, \
+    const std::vector<unsigned int> &, \
+    SparsityType &,                    \
+    const bool,                        \
+    const Table<2,bool> &) const
 #define BLOCK_SPARSITY_FUNCTIONS(SparsityType) \
   template void ConstraintMatrix::add_entries_local_to_global<SparsityType> (\
     const std::vector<unsigned int> &, \

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.