]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Made the collective add and set functions a bit faster.
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 28 Nov 2008 10:24:51 +0000 (10:24 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 28 Nov 2008 10:24:51 +0000 (10:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@17767 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_constraints.templates.h
deal.II/lac/include/lac/block_matrix_base.h
deal.II/lac/include/lac/petsc_matrix_base.h
deal.II/lac/include/lac/sparse_matrix.h
deal.II/lac/include/lac/trilinos_sparse_matrix.h

index 22f1ba11b0c4625136cfa484f8bb157a0ff5ccdb..befc30a16bd22f9abcfd370a2d6b0e937f07e6df 100644 (file)
@@ -709,8 +709,11 @@ distribute_local_to_global (const FullMatrix<double>        &local_matrix,
                                       // order to obtain a sufficient size
                                       // for the scratch array.
       n_max_entries_per_row += n_local_dofs;
-      column_indices.resize(n_max_entries_per_row);
-      column_values.resize(n_max_entries_per_row);
+      if (column_indices.size() < n_max_entries_per_row)
+        {
+         column_indices.resize(n_max_entries_per_row);
+         column_values.resize(n_max_entries_per_row);
+       }
 
                                        // now distribute entries row by row
       for (unsigned int i=0; i<n_local_dofs; ++i)
index e31ee10ad020ec96c8858e82757a80b44262d6f3..748bd12c3e61828177a337bf6f65b7d429a129e6 100644 (file)
@@ -1861,30 +1861,46 @@ BlockMatrixBase<MatrixType>::set (const unsigned int  row,
                                  const bool          elide_zero_values)
 {
                                   // Resize scratch arrays
-  column_indices.resize (this->n_block_cols());
-  column_values.resize (this->n_block_cols());
-  counter_within_block.resize (this->n_block_cols());
-
-                                  // Resize sub-arrays to n_cols. This is a
-                                  // bit wasteful, but we resize only a few
-                                  // times (then the maximum row length
-                                  // won't increase that much any more).
-  for (unsigned int i=0; i<this->n_block_cols(); ++i)
+  if (column_indices.size() < this->n_block_cols())
+    {
+      column_indices.resize (this->n_block_cols());
+      column_values.resize (this->n_block_cols());
+      counter_within_block.resize (this->n_block_cols());
+    }
+
+                                  // Resize sub-arrays to n_cols. This
+                                  // is a bit wasteful, but we resize
+                                  // only a few times (then the maximum
+                                  // row length won't increase that
+                                  // much any more). At least we know
+                                  // that all arrays are going to be of
+                                  // the same size, so we can check
+                                  // whether the size of one is large
+                                  // enough before actually going
+                                  // through all of them.
+  if (column_indices[0].size() < n_cols)
     {
-      column_indices[i].resize(n_cols);
-      column_values[i].resize(n_cols);
-      counter_within_block[i] = 0;
+      for (unsigned int i=0; i<this->n_block_cols(); ++i)
+        {
+         column_indices[i].resize(n_cols);
+         column_values[i].resize(n_cols);
+       }
     }
 
-                                  // Go through the column indices to find
-                                  // out which portions of the values
-                                  // should be written into which block
-                                  // matrix. We need to restore all the
-                                  // data, since we can't be sure that the
-                                  // data of one block is stored
-                                  // contiguously (in fact, it will be
-                                  // intermixed when the data comes from an
-                                  // element matrix).
+                                  // Reset the number of added elements
+                                  // in each block to zero.
+  for (unsigned int i=0; i<this->n_block_cols(); ++i)
+    counter_within_block[i] = 0;
+
+                                  // Go through the column indices to
+                                  // find out which portions of the
+                                  // values should be set in which
+                                  // block of the matrix. We need to
+                                  // touch all the data, since we can't
+                                  // be sure that the data of one block
+                                  // is stored contiguously (in fact,
+                                  // indices will be intermixed when it
+                                  // comes from an element matrix).
   for (unsigned int j=0; j<n_cols; ++j)
     {
       double value = values[j];
@@ -2037,31 +2053,48 @@ BlockMatrixBase<MatrixType>::add (const unsigned int   row,
                                   // TODO: Look over this to find out
                                   // whether we can do that more
                                   // efficiently.
+
                                   // Resize scratch arrays
-  column_indices.resize (this->n_block_cols());
-  column_values.resize (this->n_block_cols());
-  counter_within_block.resize (this->n_block_cols());
-
-                                  // Resize sub-arrays to n_cols. This is a
-                                  // bit wasteful, but we resize only a few
-                                  // times (then the maximum row length
-                                  // won't increase that much any more).
-  for (unsigned int i=0; i<this->n_block_cols(); ++i)
+  if (column_indices.size() < this->n_block_cols())
     {
-      column_indices[i].resize(n_cols);
-      column_values[i].resize(n_cols);
-      counter_within_block[i] = 0;
+      column_indices.resize (this->n_block_cols());
+      column_values.resize (this->n_block_cols());
+      counter_within_block.resize (this->n_block_cols());
     }
 
-                                  // Go through the column indices to find
-                                  // out which portions of the values
-                                  // should be written into which block
-                                  // matrix. We need to restore all the
-                                  // data, since we can't be sure that the
-                                  // data of one block is stored
-                                  // contiguously (in fact, it will be
-                                  // intermixed when the data comes from an
-                                  // element matrix).
+                                  // Resize sub-arrays to n_cols. This
+                                  // is a bit wasteful, but we resize
+                                  // only a few times (then the maximum
+                                  // row length won't increase that
+                                  // much any more). At least we know
+                                  // that all arrays are going to be of
+                                  // the same size, so we can check
+                                  // whether the size of one is large
+                                  // enough before actually going
+                                  // through all of them.
+  if (column_indices[0].size() < n_cols)
+    {
+      for (unsigned int i=0; i<this->n_block_cols(); ++i)
+        {
+         column_indices[i].resize(n_cols);
+         column_values[i].resize(n_cols);
+       }
+    }
+
+                                  // Reset the number of added elements
+                                  // in each block to zero.
+  for (unsigned int i=0; i<this->n_block_cols(); ++i)
+    counter_within_block[i] = 0;
+
+                                  // Go through the column indices to
+                                  // find out which portions of the
+                                  // values should be written into
+                                  // which block of the matrix. We need
+                                  // to touch all the data, since we
+                                  // can't be sure that the data of one
+                                  // block is stored contiguously (in
+                                  // fact, data will be intermixed when
+                                  // it comes from an element matrix).
   for (unsigned int j=0; j<n_cols; ++j)
     {
       double value = values[j];
index 14c92050805dffc01b81f25603c81f7e532ddd15..a6ed6d523a81af0465a0978d2fd444e4853ca739 100644 (file)
@@ -1445,8 +1445,11 @@ namespace PETScWrappers
       {
                                   // Otherwise, extract nonzero values in
                                   // each row and get the respective index.
-       column_indices.resize(n_cols);
-       column_values.resize(n_cols);
+       if (column_indices.size() < n_cols)
+         {
+           column_indices.resize(n_cols);
+           column_values.resize(n_cols);
+         }
 
        n_columns = 0;
        for (unsigned int j=0; j<n_cols; ++j)
@@ -1609,8 +1612,11 @@ namespace PETScWrappers
       {
                                   // Otherwise, extract nonzero values in
                                   // each row and get the respective index.
-       column_indices.resize(n_cols);
-       column_values.resize(n_cols);
+       if (column_indices.size() < n_cols)
+         {
+           column_indices.resize(n_cols);
+           column_values.resize(n_cols);
+         }
 
        n_columns = 0;
        for (unsigned int j=0; j<n_cols; ++j)
index 218e4bc4c9b55120fd41d63454bd66c3293715b0..c3c81d11f824c776fe4b00f856ef6ed07541bc5b 100644 (file)
@@ -2252,8 +2252,11 @@ SparseMatrix<number>::add (const unsigned int  row,
 
   unsigned int n_columns = 0;
 
-  global_indices.resize(n_cols);
-  column_values.resize(n_cols);
+  if (global_indices.size() < n_cols)
+    {
+      global_indices.resize(n_cols);
+      column_values.resize(n_cols);
+    }
 
                                   // First, search all the indices to find
                                   // out which values we actually need to
index 5ac56828d500044496652add7b7712a75de37d31..2d25253f119a9b30995a6f2733c86e6a60132edc 100755 (executable)
@@ -2170,8 +2170,11 @@ namespace TrilinosWrappers
                                   // Otherwise, extract nonzero values in
                                   // each row and get the respective
                                   // indices.
-       column_indices.resize(n_cols);
-       column_values.resize(n_cols);
+       if (column_indices.size() < n_cols)
+         {
+           column_indices.resize(n_cols);
+           column_values.resize(n_cols);
+         }
 
        n_columns = 0;
        for (unsigned int j=0; j<n_cols; ++j)
@@ -2395,8 +2398,11 @@ namespace TrilinosWrappers
       {
                                   // Otherwise, extract nonzero values in
                                   // each row and the corresponding index.
-       column_indices.resize(n_cols);
-       column_values.resize(n_cols);
+       if (column_indices.size() < n_cols)
+         {
+           column_indices.resize(n_cols);
+           column_values.resize(n_cols);
+         }
 
        n_columns = 0;
        for (unsigned int j=0; j<n_cols; ++j)

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.