]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Always use compress() with an argument, remove unneccesary compress() calls in algori...
authorTimo Heister <timo.heister@gmail.com>
Wed, 13 Feb 2013 15:26:53 +0000 (15:26 +0000)
committerTimo Heister <timo.heister@gmail.com>
Wed, 13 Feb 2013 15:26:53 +0000 (15:26 +0000)
git-svn-id: https://svn.dealii.org/trunk@28378 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/source/fe/fe_tools.cc
deal.II/source/lac/constraint_matrix.cc
deal.II/source/numerics/matrix_tools.cc

index 8707becfd3e783b8e9da69925817616133207517..26ed3a42fa02dd22feb85a5f149ee4766192e3c9 100644 (file)
@@ -1417,8 +1417,8 @@ namespace FETools
     // be cell2
     Assert (cell2 == endc2, ExcInternalError());
 
-    u2.compress();
-    touch_count.compress();
+    u2.compress(VectorOperation::add);
+    touch_count.compress(VectorOperation::add);
 
     // if we work on parallel distributed
     // vectors, we have to ensure, that we only
@@ -1439,12 +1439,9 @@ namespace FETools
         }
 
     // finish the work on parallel vectors
-    u2.compress();
+    u2.compress(VectorOperation::insert);
     // Apply hanging node constraints.
     constraints.distribute(u2);
-
-    // and finally update ghost values
-    u2.update_ghost_values();
   }
 
 
@@ -1522,10 +1519,8 @@ namespace FETools
         }
 
     // if we work on a parallel PETSc vector
-    // we have to finish the work and
-    // update ghost values
-    u1_interpolated.compress();
-    u1_interpolated.update_ghost_values();
+    // we have to finish the work
+    u1_interpolated.compress(VectorOperation::insert);
   }
 
 
@@ -1625,10 +1620,8 @@ namespace FETools
         };
 
     // if we work on a parallel PETSc vector
-    // we have to finish the work and
-    // update ghost values
-    u1_interpolated.compress();
-    u1_interpolated.update_ghost_values();
+    // we have to finish the work
+    u1_interpolated.compress(VectorOperation::insert);
   }
 
 
@@ -1771,8 +1764,7 @@ namespace FETools
     // if we work on a parallel PETSc vector
     // we have to finish the work and
     // update ghost values
-    u1_difference.compress();
-    u1_difference.update_ghost_values();
+    u1_difference.compress(VectorOperation::insert);
   }
 
 
@@ -1796,12 +1788,6 @@ namespace FETools
       {
         back_interpolate(dof1, constraints1, u1, dof2, constraints2, u1_difference);
         u1_difference.sadd(-1, u1);
-
-        // if we work on a parallel PETSc vector
-        // we have to finish the work and
-        // update ghost values
-        u1_difference.compress();
-        u1_difference.update_ghost_values();
       }
   }
 
index a03a4791499e09f623bf77fc4039a21798185173..bba9512788f2079e5fae3afe235c3ebafb04130d 100644 (file)
@@ -2106,7 +2106,7 @@ ConstraintMatrix::distribute (PETScWrappers::MPI::Vector &vec) const
       vec(it->line) = new_value;
     }
 
-  vec.compress ();
+  vec.compress (VectorOperation::insert);
 }
 
 
index d48a26831cd5e84d5621fe88178d3c5a77edb0c1..2282d6cd41ad38bd63aa32705cfd41e2b6a194a6 100644 (file)
@@ -2384,13 +2384,6 @@ namespace MatrixTools
         Assert (local_range == solution.local_range(),
                 ExcInternalError());
 
-
-        // we have to read and write from this
-        // matrix (in this order). this will only
-        // work if we compress the matrix first,
-        // done here
-        matrix.compress ();
-
         // determine the first nonzero diagonal
         // entry from within the part of the
         // matrix that we can see. if we can't
@@ -2427,23 +2420,6 @@ namespace MatrixTools
         // treated in the other functions.
         matrix.clear_rows (constrained_rows, average_nonzero_diagonal_entry);
 
-        // the next thing is to set right hand
-        // side to the wanted value. there's one
-        // drawback: if we write to individual
-        // vector elements, then we have to do
-        // that on all processors. however, some
-        // processors may not need to set
-        // anything because their chunk of
-        // matrix/rhs do not contain any boundary
-        // nodes. therefore, rather than using
-        // individual calls, we use one call for
-        // all elements, thereby making sure that
-        // all processors call this function,
-        // even if some only have an empty set of
-        // elements to set
-        right_hand_side.compress ();
-        solution.compress ();
-
         std::vector<unsigned int> indices;
         std::vector<PetscScalar>  solution_values;
         for (std::map<unsigned int,double>::const_iterator
@@ -2467,8 +2443,8 @@ namespace MatrixTools
 
         // clean up
         matrix.compress ();
-        solution.compress ();
-        right_hand_side.compress ();
+        solution.compress (VectorOperation::insert);
+        right_hand_side.compress (VectorOperation::insert);
       }
     }
   }
@@ -2656,27 +2632,6 @@ namespace MatrixTools
         // matrix classes in deal.II.
         matrix.clear_rows (constrained_rows, average_nonzero_diagonal_entry);
 
-        // the next thing is to set right
-        // hand side to the wanted
-        // value. there's one drawback:
-        // if we write to individual
-        // vector elements, then we have
-        // to do that on all
-        // processors. however, some
-        // processors may not need to set
-        // anything because their chunk
-        // of matrix/rhs do not contain
-        // any boundary nodes. therefore,
-        // rather than using individual
-        // calls, we use one call for all
-        // elements, thereby making sure
-        // that all processors call this
-        // function, even if some only
-        // have an empty set of elements
-        // to set
-        right_hand_side.compress ();
-        solution.compress ();
-
         std::vector<unsigned int> indices;
         std::vector<TrilinosScalar>  solution_values;
         for (std::map<unsigned int,double>::const_iterator
@@ -2700,8 +2655,8 @@ namespace MatrixTools
 
         // clean up
         matrix.compress ();
-        solution.compress ();
-        right_hand_side.compress ();
+        solution.compress (VectorOperation::insert);
+        right_hand_side.compress (VectorOperation::insert);
       }
 
 

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.