From a374488a4ef348e7c6dea0eb2eb0b867fbfcc954 Mon Sep 17 00:00:00 2001
From: wolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Thu, 1 Apr 2004 17:00:13 +0000
Subject: [PATCH] Run the thing in parallel -- yeay!

git-svn-id: https://svn.dealii.org/trunk@8946 0785d39b-7218-0410-832d-ea1e28bc413d
---
 deal.II/examples/step-17/Makefile   |   4 +-
 deal.II/examples/step-17/step-17.cc | 183 +++++++++++++++-------------
 2 files changed, 101 insertions(+), 86 deletions(-)

diff --git a/deal.II/examples/step-17/Makefile b/deal.II/examples/step-17/Makefile
index 2d7deb8ece..57de7de067 100644
--- a/deal.II/examples/step-17/Makefile
+++ b/deal.II/examples/step-17/Makefile
@@ -76,10 +76,10 @@ else
 #
 # You may need to augment the lists of libraries when compiling your
 # program for other dimensions, or when using third party libraries
-libs.g   = $(lib-deal2-2d.g) \
+libs.g   = $(lib-deal2-3d.g) \
 	   $(lib-lac.g)      \
            $(lib-base.g)
-libs.o   = $(lib-deal2-2d.o) \
+libs.o   = $(lib-deal2-3d.o) \
 	   $(lib-lac.o)      \
            $(lib-base.o)
 
diff --git a/deal.II/examples/step-17/step-17.cc b/deal.II/examples/step-17/step-17.cc
index 0ffdc325c4..d87c5a32d4 100644
--- a/deal.II/examples/step-17/step-17.cc
+++ b/deal.II/examples/step-17/step-17.cc
@@ -227,15 +227,15 @@ void ElasticProblem<dim>::setup_system ()
                         local_dofs,
                         dof_handler.max_couplings_between_dofs());
 
-  solution.reinit (dof_handler.n_dofs(), local_dofs, mpi_communicator);
-  system_rhs.reinit (dof_handler.n_dofs(), local_dofs, mpi_communicator);
+  solution.reinit (mpi_communicator, dof_handler.n_dofs(), local_dofs);
+  system_rhs.reinit (mpi_communicator, dof_handler.n_dofs(), local_dofs);
 }
 
 
 template <int dim>
 void ElasticProblem<dim>::assemble_system () 
 {
-                                   // move to front
+                                   // xxx move to front
   std::map<unsigned int,double> boundary_values;
   VectorTools::interpolate_boundary_values (dof_handler,
 					    0,
@@ -270,84 +270,86 @@ void ElasticProblem<dim>::assemble_system ()
   typename DoFHandler<dim>::active_cell_iterator cell = dof_handler.begin_active(),
 						 endc = dof_handler.end();
   for (; cell!=endc; ++cell)
-    {
-      cell_matrix.clear ();
-      cell_rhs.clear ();
+                                     // xxx
+    if (cell->subdomain_id() == this_partition)
+      {
+        cell_matrix.clear ();
+        cell_rhs.clear ();
 
-      fe_values.reinit (cell);
+        fe_values.reinit (cell);
       
-      lambda.value_list (fe_values.get_quadrature_points(), lambda_values);
-      mu.value_list     (fe_values.get_quadrature_points(), mu_values);
+        lambda.value_list (fe_values.get_quadrature_points(), lambda_values);
+        mu.value_list     (fe_values.get_quadrature_points(), mu_values);
 
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-	{
-	  const unsigned int 
-	    component_i = fe.system_to_component_index(i).first;
+        for (unsigned int i=0; i<dofs_per_cell; ++i)
+          {
+            const unsigned int 
+              component_i = fe.system_to_component_index(i).first;
 	  
-	  for (unsigned int j=0; j<dofs_per_cell; ++j) 
-	    {
-	      const unsigned int 
-		component_j = fe.system_to_component_index(j).first;
+            for (unsigned int j=0; j<dofs_per_cell; ++j) 
+              {
+                const unsigned int 
+                  component_j = fe.system_to_component_index(j).first;
 	      
-	      for (unsigned int q_point=0; q_point<n_q_points;
-		   ++q_point)
-		{
-		  cell_matrix(i,j) 
-		    += 
-		    (
-		      (fe_values.shape_grad(i,q_point)[component_i] *
-		       fe_values.shape_grad(j,q_point)[component_j] *
-		       lambda_values[q_point])
-		      +
-		      (fe_values.shape_grad(i,q_point)[component_j] *
-		       fe_values.shape_grad(j,q_point)[component_i] *
-		       mu_values[q_point])
-		      +
-		      ((component_i == component_j) ?
-		       (fe_values.shape_grad(i,q_point) *
-			fe_values.shape_grad(j,q_point) *
-			mu_values[q_point])  :
-		       0)
-		    )
-		    *
-		    fe_values.JxW(q_point);
-		};
-	    };
-	};
-
-      right_hand_side.vector_value_list (fe_values.get_quadrature_points(),
-					 rhs_values);
-      for (unsigned int i=0; i<dofs_per_cell; ++i)
-	{
-	  const unsigned int 
-	    component_i = fe.system_to_component_index(i).first;
+                for (unsigned int q_point=0; q_point<n_q_points;
+                     ++q_point)
+                  {
+                    cell_matrix(i,j) 
+                      += 
+                      (
+                        (fe_values.shape_grad(i,q_point)[component_i] *
+                         fe_values.shape_grad(j,q_point)[component_j] *
+                         lambda_values[q_point])
+                        +
+                        (fe_values.shape_grad(i,q_point)[component_j] *
+                         fe_values.shape_grad(j,q_point)[component_i] *
+                         mu_values[q_point])
+                        +
+                        ((component_i == component_j) ?
+                         (fe_values.shape_grad(i,q_point) *
+                          fe_values.shape_grad(j,q_point) *
+                          mu_values[q_point])  :
+                         0)
+                        )
+                      *
+                      fe_values.JxW(q_point);
+                  };
+              };
+          };
+
+        right_hand_side.vector_value_list (fe_values.get_quadrature_points(),
+                                           rhs_values);
+        for (unsigned int i=0; i<dofs_per_cell; ++i)
+          {
+            const unsigned int 
+              component_i = fe.system_to_component_index(i).first;
 	  
-	  for (unsigned int q_point=0; q_point<n_q_points; ++q_point)
-	    cell_rhs(i) += fe_values.shape_value(i,q_point) *
-			   rhs_values[q_point](component_i) *
-			   fe_values.JxW(q_point);
-	};
-
-      cell->get_dof_indices (local_dof_indices);
-
-                                       //xxx
-      MatrixTools::local_apply_boundary_values (boundary_values,
-                                                local_dof_indices,
-                                                cell_matrix,
-                                                cell_rhs,
-                                                false);
-
-                                       // xxx
-      hanging_node_constraints
-        .distribute_local_to_global (cell_matrix,
-                                     local_dof_indices,
-                                     system_matrix);
-
-      hanging_node_constraints
-        .distribute_local_to_global (cell_rhs,
-                                     local_dof_indices,
-                                     system_rhs);
-    }
+            for (unsigned int q_point=0; q_point<n_q_points; ++q_point)
+              cell_rhs(i) += fe_values.shape_value(i,q_point) *
+                             rhs_values[q_point](component_i) *
+                             fe_values.JxW(q_point);
+          };
+
+        cell->get_dof_indices (local_dof_indices);
+
+                                         //xxx
+        MatrixTools::local_apply_boundary_values (boundary_values,
+                                                  local_dof_indices,
+                                                  cell_matrix,
+                                                  cell_rhs,
+                                                  false);
+
+                                         // xxx
+        hanging_node_constraints
+          .distribute_local_to_global (cell_matrix,
+                                       local_dof_indices,
+                                       system_matrix);
+
+        hanging_node_constraints
+          .distribute_local_to_global (cell_rhs,
+                                       local_dof_indices,
+                                       system_rhs);
+      }
 
                                    //xxx no condense necessary, no apply_b_v
                                    //either
@@ -365,14 +367,19 @@ void ElasticProblem<dim>::solve ()
 {
                                    // xxx
   SolverControl           solver_control (1000, 1e-10);
-  PETScWrappers::SolverCG cg (solver_control);
-
-  PETScWrappers::PreconditionSSOR preconditioner(system_matrix, 1.2);
+  PETScWrappers::SolverCG cg (solver_control,
+                              mpi_communicator);
 
+  PETScWrappers::PreconditionBlockJacobi preconditioner(system_matrix);
+  
   cg.solve (system_matrix, solution, system_rhs,
 	    preconditioner);
-
   hanging_node_constraints.distribute (solution);
+
+  if (this_partition == 0)
+    std::cout << "   Solver converged in "
+              << solver_control.last_step()
+              << " iterations." << std::endl;
 }
 
 
@@ -460,7 +467,7 @@ void ElasticProblem<dim>::run ()
       if (cycle == 0)
 	{
 	  GridGenerator::hyper_cube (triangulation, -1, 1);
-	  triangulation.refine_global (2);
+	  triangulation.refine_global (4);
 
                                            // xxx
           GridTools::partition_triangulation (n_partitions, triangulation);
@@ -478,9 +485,17 @@ void ElasticProblem<dim>::run ()
 
                                        // xxx
       if (this_partition == 0)
-        std::cout << "   Number of degrees of freedom: "
-                  << dof_handler.n_dofs()
-                  << std::endl;
+        {
+          std::cout << "   Number of degrees of freedom: "
+                    << dof_handler.n_dofs()
+                    << " (by partition:";
+          for (unsigned int partition=0; partition<n_partitions; ++partition)
+            std::cout << (partition==0 ? ' ' : '+')
+                      << (DoFTools::
+                          count_dofs_with_subdomain_association (dof_handler,
+                                                                 partition));
+          std::cout << ")" << std::endl;
+        }
       
       assemble_system ();
       solve ();
@@ -499,7 +514,7 @@ int main (int argc, char **argv)
 
                                        // xxx localize scope
       {
-        ElasticProblem<2> elastic_problem_2d;
+        ElasticProblem<3> elastic_problem_2d;
         elastic_problem_2d.run ();
       }
 
-- 
2.39.5