[cig-commits] [commit] master: fix denormalize pressure in parallel (f83fde5)

[cig_noreply at geodynamics.org]

Repository : https://github.com/geodynamics/aspect

On branch  : master
Link       : https://github.com/geodynamics/aspect/compare/21aac238409facd328313998aa9ab9dddb2f16a7...c3b8728148f634b279a3a9e6c4013e696081cf12

>---------------------------------------------------------------

commit f83fde5dfc838a50f608ea8afc215511dfeb5aba
Author: Timo Heister <timo.heister at gmail.com>
Date:   Fri Jun 13 08:54:28 2014 -0400

    fix denormalize pressure in parallel
    
    We were not handing a ghosted vector to denormalize pressure when using
    a direct solver, which caused a runtime error.


>---------------------------------------------------------------

f83fde5dfc838a50f608ea8afc215511dfeb5aba
 source/simulator/helper_functions.cc | 5 +++--
 source/simulator/solver.cc           | 5 +++--
 2 files changed, 6 insertions(+), 4 deletions(-)

diff --git a/source/simulator/helper_functions.cc b/source/simulator/helper_functions.cc
index 4c9bc79..1630bbb 100644
--- a/source/simulator/helper_functions.cc
+++ b/source/simulator/helper_functions.cc
@@ -725,9 +725,10 @@ namespace aspect
           {
             // velocity and pressure are in the same block, so we have to modify the values manually
 
+            const unsigned int block_idx = introspection.block_indices.pressure;
             LinearAlgebra::BlockVector distributed_vector (introspection.index_sets.stokes_partitioning,
                                                            mpi_communicator);
-            distributed_vector = vector;
+            distributed_vector.block(block_idx) = vector.block(block_idx);
 
             std::vector<types::global_dof_index> local_dof_indices (finite_element.dofs_per_cell);
             typename DoFHandler<dim>::active_cell_iterator
@@ -751,7 +752,7 @@ namespace aspect
                     }
                 }
             distributed_vector.compress(VectorOperation::insert);
-            vector = distributed_vector;
+            vector.block(block_idx) = distributed_vector.block(block_idx);
           }
       }
     else
diff --git a/source/simulator/solver.cc b/source/simulator/solver.cc
index d8a6579..11a8197 100644
--- a/source/simulator/solver.cc
+++ b/source/simulator/solver.cc
@@ -471,8 +471,9 @@ namespace aspect
         // nonlinear residual (see initial_residual below).
         // TODO: if there was an easy way to know if the caller needs the
         // initial residual we could skip all of this stuff.
-        distributed_stokes_solution.block(0) = current_linearization_point.block(0);
-        denormalize_pressure (distributed_stokes_solution);
+        solution.block(0) = current_linearization_point.block(0);
+        denormalize_pressure (solution);
+        distributed_stokes_solution.block(0) = solution.block(0);
         current_constraints.set_zero (distributed_stokes_solution);
 
         // Undo the pressure scaling: