From 45514f7f0e04fecd50bdc29549a3355fb9047883 Mon Sep 17 00:00:00 2001 From: wolf Date: Wed, 24 Mar 2004 23:24:32 +0000 Subject: [PATCH] Inch forward. git-svn-id: https://svn.dealii.org/trunk@8870 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-17/step-17.cc | 169 +++++++++++++++++++--------- 1 file changed, 115 insertions(+), 54 deletions(-) diff --git a/deal.II/examples/step-17/step-17.cc b/deal.II/examples/step-17/step-17.cc index 7e827bed64..0ffdc325c4 100644 --- a/deal.II/examples/step-17/step-17.cc +++ b/deal.II/examples/step-17/step-17.cc @@ -19,7 +19,8 @@ // xxx #include -#include +#include +#include #include #include @@ -32,6 +33,10 @@ #include #include #include + + // xxx +#include + #include #include #include @@ -71,14 +76,21 @@ class ElasticProblem ConstraintMatrix hanging_node_constraints; // xxx no sparsity - PETScWrappers::SparseMatrix system_matrix; + PETScWrappers::MPI::SparseMatrix system_matrix; - PETScWrappers::Vector solution; - PETScWrappers::Vector system_rhs; + PETScWrappers::MPI::Vector solution; + PETScWrappers::MPI::Vector system_rhs; + MPI_Comm mpi_communicator; + // xxx const unsigned int n_partitions; const unsigned int this_partition; + + unsigned int local_dofs; + + static unsigned int get_n_partitions (const MPI_Comm &mpi_communicator); + static unsigned int get_this_partition (const MPI_Comm &mpi_communicator); }; @@ -145,6 +157,31 @@ void RightHandSide::vector_value_list (const std::vector > &poin + // xxx +template +unsigned int +ElasticProblem::get_n_partitions (const MPI_Comm &mpi_communicator) +{ + // xxx int vs uint + int n_jobs; + MPI_Comm_size (mpi_communicator, &n_jobs); + + return n_jobs; +} + + + +template +unsigned int +ElasticProblem::get_this_partition (const MPI_Comm &mpi_communicator) +{ + int rank; + MPI_Comm_rank (mpi_communicator, &rank); + + return rank; +} + + template ElasticProblem::ElasticProblem () @@ -152,8 +189,9 @@ ElasticProblem::ElasticProblem () dof_handler (triangulation), fe (FE_Q(1), dim), // xxx - n_partitions (8), - this_partition (0) + mpi_communicator (MPI_COMM_WORLD), + n_partitions (get_n_partitions(mpi_communicator)), + this_partition (get_this_partition(mpi_communicator)) {} @@ -168,19 +206,29 @@ ElasticProblem::~ElasticProblem () template void ElasticProblem::setup_system () { + // xxx dof_handler.distribute_dofs (fe); + DoFRenumbering::subdomain_wise (dof_handler); + + local_dofs + = DoFTools::count_dofs_with_subdomain_association (dof_handler, + this_partition); + + hanging_node_constraints.clear (); DoFTools::make_hanging_node_constraints (dof_handler, hanging_node_constraints); hanging_node_constraints.close (); // no sparsity pattern - system_matrix.reinit (dof_handler.n_dofs(), + system_matrix.reinit (mpi_communicator, + dof_handler.n_dofs(), dof_handler.n_dofs(), + local_dofs, dof_handler.max_couplings_between_dofs()); - solution.reinit (dof_handler.n_dofs()); - system_rhs.reinit (dof_handler.n_dofs()); + solution.reinit (dof_handler.n_dofs(), local_dofs, mpi_communicator); + system_rhs.reinit (dof_handler.n_dofs(), local_dofs, mpi_communicator); } @@ -334,12 +382,14 @@ void ElasticProblem::refine_grid () { Vector estimated_error_per_cell (triangulation.n_active_cells()); + // XXX typename FunctionMap::type neumann_boundary; - KellyErrorEstimator::estimate (dof_handler, - QGauss2(), - neumann_boundary, - solution, - estimated_error_per_cell); + Assert (false, ExcNotImplemented()); +// KellyErrorEstimator::estimate (dof_handler, +// QGauss2(), +// neumann_boundary, +// solution, +// estimated_error_per_cell); GridRefinement::refine_and_coarsen_fixed_number (triangulation, estimated_error_per_cell, @@ -355,40 +405,45 @@ void ElasticProblem::refine_grid () template void ElasticProblem::output_results (const unsigned int cycle) const { - std::string filename = "solution-"; - filename += ('0' + cycle); - Assert (cycle < 10, ExcInternalError()); - - filename += ".gmv"; - std::ofstream output (filename.c_str()); - - DataOut data_out; - data_out.attach_dof_handler (dof_handler); - - - - std::vector solution_names; - switch (dim) + // xxx + PETScWrappers::Vector global_solution; + global_solution = solution; + + if (this_partition == 0) { - case 1: - solution_names.push_back ("displacement"); - break; - case 2: - solution_names.push_back ("x_displacement"); - solution_names.push_back ("y_displacement"); - break; - case 3: - solution_names.push_back ("x_displacement"); - solution_names.push_back ("y_displacement"); - solution_names.push_back ("z_displacement"); - break; - default: - Assert (false, ExcInternalError()); - }; - - data_out.add_data_vector (solution, solution_names); - data_out.build_patches (); - data_out.write_gmv (output); + std::string filename = "solution-"; + filename += ('0' + cycle); + Assert (cycle < 10, ExcInternalError()); + + filename += ".gmv"; + std::ofstream output (filename.c_str()); + + DataOut data_out; + data_out.attach_dof_handler (dof_handler); + + std::vector solution_names; + switch (dim) + { + case 1: + solution_names.push_back ("displacement"); + break; + case 2: + solution_names.push_back ("x_displacement"); + solution_names.push_back ("y_displacement"); + break; + case 3: + solution_names.push_back ("x_displacement"); + solution_names.push_back ("y_displacement"); + solution_names.push_back ("z_displacement"); + break; + default: + Assert (false, ExcInternalError()); + }; + + data_out.add_data_vector (global_solution, solution_names); + data_out.build_patches (); + data_out.write_gmv (output); + } } @@ -398,7 +453,9 @@ void ElasticProblem::run () { for (unsigned int cycle=0; cycle<12; ++cycle) { - std::cout << "Cycle " << cycle << ':' << std::endl; + // xxx + if (this_partition == 0) + std::cout << "Cycle " << cycle << ':' << std::endl; if (cycle == 0) { @@ -411,15 +468,19 @@ void ElasticProblem::run () else refine_grid (); - std::cout << " Number of active cells: " - << triangulation.n_active_cells() - << std::endl; + // xxx + if (this_partition == 0) + std::cout << " Number of active cells: " + << triangulation.n_active_cells() + << std::endl; setup_system (); - std::cout << " Number of degrees of freedom: " - << dof_handler.n_dofs() - << std::endl; + // xxx + if (this_partition == 0) + std::cout << " Number of degrees of freedom: " + << dof_handler.n_dofs() + << std::endl; assemble_system (); solve (); -- 2.39.5