From: Wolfgang Bangerth Date: Thu, 1 Apr 2004 17:00:13 +0000 (+0000) Subject: Run the thing in parallel -- yeay! X-Git-Tag: v8.0.0~15402 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5d5b66935bc4ff5c4eb47dde429df52eb9c49494;p=dealii.git Run the thing in parallel -- yeay! git-svn-id: https://svn.dealii.org/trunk@8946 0785d39b-7218-0410-832d-ea1e28bc413d --- 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::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 void ElasticProblem::assemble_system () { - // move to front + // xxx move to front std::map boundary_values; VectorTools::interpolate_boundary_values (dof_handler, 0, @@ -270,84 +270,86 @@ void ElasticProblem::assemble_system () typename DoFHandler::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; iget_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_pointget_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::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::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::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 elastic_problem_2d; + ElasticProblem<3> elastic_problem_2d; elastic_problem_2d.run (); }