From: kronbichler Date: Thu, 4 Sep 2008 09:28:15 +0000 (+0000) Subject: Made relevant changes for running this program in parallel. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1e6f2d8f94a6a214d89d9460dcdc8c37bd8d66ce;p=dealii-svn.git Made relevant changes for running this program in parallel. git-svn-id: https://svn.dealii.org/trunk@16742 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/examples/step-32/step-32.cc b/deal.II/examples/step-32/step-32.cc index 7e8da2e2d7..661ce22747 100644 --- a/deal.II/examples/step-32/step-32.cc +++ b/deal.II/examples/step-32/step-32.cc @@ -20,6 +20,7 @@ #include #include #include +#include #include #include @@ -55,10 +56,16 @@ #include #include -#include +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI +# include +#else +# include +#endif + #include #include +#include #include @@ -325,7 +332,13 @@ class BoussinesqFlowProblem const double old_time_step); +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + Epetra_MpiComm trilinos_communicator; +#else Epetra_SerialComm trilinos_communicator; +#endif + + ConditionalOStream pcout; Triangulation triangulation; @@ -377,6 +390,11 @@ class BoussinesqFlowProblem template BoussinesqFlowProblem::BoussinesqFlowProblem () : +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + trilinos_communicator (MPI_COMM_WORLD), +#endif + pcout (std::cout, trilinos_communicator.MyPID()==0), + triangulation (Triangulation::maximum_smoothing), stokes_degree (1), @@ -551,11 +569,15 @@ void BoussinesqFlowProblem::setup_dofs () { std::vector stokes_block_component (dim+1,0); stokes_block_component[dim] = 1; - + + GridTools::partition_triangulation (trilinos_communicator.NumProc(), + triangulation); + { stokes_dof_handler.distribute_dofs (stokes_fe); DoFRenumbering::Cuthill_McKee (stokes_dof_handler); DoFRenumbering::component_wise (stokes_dof_handler, stokes_block_component); + DoFRenumbering::subdomain_wise (stokes_dof_handler); stokes_constraints.clear (); DoFTools::make_hanging_node_constraints (stokes_dof_handler, @@ -570,6 +592,7 @@ void BoussinesqFlowProblem::setup_dofs () { temperature_dof_handler.distribute_dofs (temperature_fe); DoFRenumbering::Cuthill_McKee (temperature_dof_handler); + DoFRenumbering::subdomain_wise (temperature_dof_handler); temperature_constraints.clear (); DoFTools::make_hanging_node_constraints (temperature_dof_handler, @@ -585,17 +608,17 @@ void BoussinesqFlowProblem::setup_dofs () n_p = stokes_dofs_per_block[1], n_T = temperature_dof_handler.n_dofs(); - std::cout << "Number of active cells: " - << triangulation.n_active_cells() - << " (on " - << triangulation.n_levels() - << " levels)" - << std::endl - << "Number of degrees of freedom: " - << n_u + n_p + n_T - << " (" << n_u << '+' << n_p << '+'<< n_T <<')' - << std::endl - << std::endl; + pcout << "Number of active cells: " + << triangulation.n_active_cells() + << " (on " + << triangulation.n_levels() + << " levels)" + << std::endl + << "Number of degrees of freedom: " + << n_u + n_p + n_T + << " (" << n_u << '+' << n_p << '+'<< n_T <<')' + << std::endl + << std::endl; stokes_partitioner.clear(); { @@ -607,31 +630,35 @@ void BoussinesqFlowProblem::setup_dofs () { stokes_matrix.clear (); - BlockCompressedSetSparsityPattern csp (2,2); - - csp.block(0,0).reinit (n_u, n_u); - csp.block(0,1).reinit (n_u, n_p); - csp.block(1,0).reinit (n_p, n_u); - csp.block(1,1).reinit (n_p, n_p); + if (trilinos_communicator.MyPID() == 0) + { + BlockCompressedSetSparsityPattern csp (2,2); + + csp.block(0,0).reinit (n_u, n_u); + csp.block(0,1).reinit (n_u, n_p); + csp.block(1,0).reinit (n_p, n_u); + csp.block(1,1).reinit (n_p, n_p); - csp.collect_sizes (); + csp.collect_sizes (); - Table<2,DoFTools::Coupling> coupling (dim+1, dim+1); + Table<2,DoFTools::Coupling> coupling (dim+1, dim+1); - for (unsigned int c=0; c::setup_dofs () Mp_preconditioner.reset (); stokes_preconditioner_matrix.clear (); - BlockCompressedSetSparsityPattern csp (2,2); - - csp.block(0,0).reinit (n_u, n_u); - csp.block(0,1).reinit (n_u, n_p); - csp.block(1,0).reinit (n_p, n_u); - csp.block(1,1).reinit (n_p, n_p); - - csp.collect_sizes (); - - Table<2,DoFTools::Coupling> coupling (dim+1, dim+1); - for (unsigned int c=0; c coupling (dim+1, dim+1); + for (unsigned int c=0; c::setup_dofs () temperature_stiffness_matrix.clear (); temperature_matrix.clear (); - CompressedSetSparsityPattern csp (n_T, n_T); - DoFTools::make_sparsity_pattern (temperature_dof_handler, csp); - temperature_constraints.condense (csp); - - SparsityPattern temperature_sparsity_pattern; - temperature_sparsity_pattern.copy_from (csp); - - temperature_matrix.reinit (temperature_partitioner, - temperature_sparsity_pattern); - temperature_mass_matrix.reinit (temperature_partitioner, - temperature_sparsity_pattern); - temperature_stiffness_matrix.reinit (temperature_partitioner, - temperature_sparsity_pattern); + if (trilinos_communicator.MyPID() == 0) + { + CompressedSetSparsityPattern csp (n_T, n_T); + DoFTools::make_sparsity_pattern (temperature_dof_handler, csp); + temperature_constraints.condense (csp); + + SparsityPattern temperature_sparsity_pattern; + temperature_sparsity_pattern.copy_from (csp); + + temperature_matrix.reinit (temperature_partitioner, + temperature_sparsity_pattern); + temperature_mass_matrix.reinit (temperature_partitioner, + temperature_sparsity_pattern); + temperature_stiffness_matrix.reinit (temperature_partitioner, + temperature_sparsity_pattern); + } + temperature_matrix.compress(); + temperature_mass_matrix.compress(); + temperature_stiffness_matrix.compress(); } stokes_solution.reinit (stokes_partitioner); @@ -729,31 +766,33 @@ BoussinesqFlowProblem::assemble_stokes_preconditioner () cell = stokes_dof_handler.begin_active(), endc = stokes_dof_handler.end(); for (; cell!=endc; ++cell) - { - stokes_fe_values.reinit (cell); - local_matrix = 0; + if (cell->subdomain_id() == (unsigned int)trilinos_communicator.MyPID()) + { + stokes_fe_values.reinit (cell); + local_matrix = 0; - for (unsigned int q=0; qget_dof_indices (local_dof_indices); + stokes_constraints.distribute_local_to_global (local_matrix, + local_dof_indices, + stokes_preconditioner_matrix); + } - cell->get_dof_indices (local_dof_indices); - stokes_constraints.distribute_local_to_global (local_matrix, - local_dof_indices, - stokes_preconditioner_matrix); - } stokes_preconditioner_matrix.compress(); } @@ -765,11 +804,11 @@ BoussinesqFlowProblem::build_stokes_preconditioner () { if (rebuild_stokes_preconditioner == false) return; - - std::cout << " Rebuilding Stokes preconditioner..." << std::flush; - + + pcout << " Rebuilding Stokes preconditioner..." << std::flush; + assemble_stokes_preconditioner (); - + Amg_preconditioner = boost::shared_ptr (new TrilinosWrappers::PreconditionAMG()); @@ -779,13 +818,13 @@ BoussinesqFlowProblem::build_stokes_preconditioner () DoFTools::extract_constant_modes (stokes_dof_handler, velocity_components, null_space); Amg_preconditioner->initialize(stokes_preconditioner_matrix.block(0,0), - true, true, null_space, false); + true, true, 5e-2, null_space, false); Mp_preconditioner = boost::shared_ptr (new TrilinosWrappers::PreconditionSSOR( stokes_preconditioner_matrix.block(1,1),1.2)); - std::cout << std::endl; + pcout << std::endl; rebuild_stokes_preconditioner = false; } @@ -796,7 +835,7 @@ BoussinesqFlowProblem::build_stokes_preconditioner () template void BoussinesqFlowProblem::assemble_stokes_system () { - std::cout << " Assembling..." << std::flush; + pcout << " Assembling..." << std::flush; if (rebuild_stokes_matrix == true) stokes_matrix=0; @@ -857,87 +896,89 @@ void BoussinesqFlowProblem::assemble_stokes_system () temperature_cell = temperature_dof_handler.begin_active(); for (; cell!=endc; ++cell, ++temperature_cell) - { - stokes_fe_values.reinit (cell); - temperature_fe_values.reinit (temperature_cell); - - local_matrix = 0; - local_rhs = 0; - - temperature_fe_values.get_function_values (old_temperature_solution, old_temperature_values); - - for (unsigned int q=0; q gravity = ( (dim == 2) ? (Point (0,1)) : - (Point (0,0,1)) ); - for (unsigned int i=0; i::faces_per_cell; - ++face_no) - if (cell->at_boundary(face_no)) + if (cell->subdomain_id() == (unsigned int)trilinos_communicator.MyPID()) + { + stokes_fe_values.reinit (cell); + temperature_fe_values.reinit (temperature_cell); + + local_matrix = 0; + local_rhs = 0; + + temperature_fe_values.get_function_values (old_temperature_solution, + old_temperature_values); + + for (unsigned int q=0; q - phi_i_u = stokes_fe_face_values[velocities].value (i, q); - - local_rhs(i) += -(phi_i_u * - stokes_fe_face_values.normal_vector(q) * - boundary_values[q] * - stokes_fe_face_values.JxW(q)); - } - } - - cell->get_dof_indices (local_dof_indices); - - if (rebuild_stokes_matrix == true) - stokes_constraints.distribute_local_to_global (local_matrix, + for (unsigned int j=0; j gravity = ( (dim == 2) ? (Point (0,1)) : + (Point (0,0,1)) ); + for (unsigned int i=0; i::faces_per_cell; + ++face_no) + if (cell->at_boundary(face_no)) + { + stokes_fe_face_values.reinit (cell, face_no); + + pressure_boundary_values + .value_list (stokes_fe_face_values.get_quadrature_points(), + boundary_values); + + for (unsigned int q=0; q + phi_i_u = stokes_fe_face_values[velocities].value (i, q); + + local_rhs(i) += -(phi_i_u * + stokes_fe_face_values.normal_vector(q) * + boundary_values[q] * + stokes_fe_face_values.JxW(q)); + } + } + + cell->get_dof_indices (local_dof_indices); + + if (rebuild_stokes_matrix == true) + stokes_constraints.distribute_local_to_global (local_matrix, + local_dof_indices, + stokes_matrix); + + stokes_constraints.distribute_local_to_global (local_rhs, local_dof_indices, - stokes_matrix); - - stokes_constraints.distribute_local_to_global (local_rhs, - local_dof_indices, - stokes_rhs); - } + stokes_rhs); + } stokes_matrix.compress(); stokes_rhs.compress(); rebuild_stokes_matrix = false; - std::cout << std::endl; + pcout << std::endl; } @@ -977,44 +1018,48 @@ void BoussinesqFlowProblem::assemble_temperature_matrix () cell = temperature_dof_handler.begin_active(), endc = temperature_dof_handler.end(); for (; cell!=endc; ++cell) - { - local_mass_matrix = 0; - local_stiffness_matrix = 0; - - temperature_fe_values.reinit (cell); - - for (unsigned int q=0; qsubdomain_id() == (unsigned int)trilinos_communicator.MyPID()) + { + local_mass_matrix = 0; + local_stiffness_matrix = 0; + + temperature_fe_values.reinit (cell); + + for (unsigned int q=0; qget_dof_indices (local_dof_indices); - - temperature_constraints.distribute_local_to_global (local_mass_matrix, - local_dof_indices, - temperature_mass_matrix); - temperature_constraints.distribute_local_to_global (local_stiffness_matrix, - local_dof_indices, - temperature_stiffness_matrix); - } + + for (unsigned int i=0; iget_dof_indices (local_dof_indices); + temperature_constraints.distribute_local_to_global (local_mass_matrix, + local_dof_indices, + temperature_mass_matrix); + temperature_constraints.distribute_local_to_global (local_stiffness_matrix, + local_dof_indices, + temperature_stiffness_matrix); + } + + temperature_mass_matrix.compress(); + temperature_stiffness_matrix.compress(); + rebuild_temperature_matrices = false; } @@ -1038,6 +1083,7 @@ void BoussinesqFlowProblem::assemble_temperature_system () temperature_matrix.copy_from (temperature_mass_matrix); temperature_matrix.add (time_step, temperature_stiffness_matrix); } + temperature_matrix.compress(); temperature_rhs = 0; @@ -1087,119 +1133,122 @@ void BoussinesqFlowProblem::assemble_temperature_system () stokes_cell = stokes_dof_handler.begin_active(); for (; cell!=endc; ++cell, ++stokes_cell) - { - local_rhs = 0; - - temperature_fe_values.reinit (cell); - stokes_fe_values.reinit (stokes_cell); - - temperature_fe_values.get_function_values (old_temperature_solution, - old_temperature_values); - temperature_fe_values.get_function_values (old_old_temperature_solution, - old_old_temperature_values); - - temperature_fe_values.get_function_gradients (old_temperature_solution, - old_temperature_grads); - temperature_fe_values.get_function_gradients (old_old_temperature_solution, - old_old_temperature_grads); - - temperature_fe_values.get_function_hessians (old_temperature_solution, - old_temperature_hessians); - temperature_fe_values.get_function_hessians (old_old_temperature_solution, - old_old_temperature_hessians); - - temperature_right_hand_side.value_list (temperature_fe_values.get_quadrature_points(), - gamma_values); - - stokes_fe_values.get_function_values (stokes_solution, - present_stokes_values); - - const double nu - = compute_viscosity (old_temperature_values, - old_old_temperature_values, - old_temperature_grads, - old_old_temperature_grads, - old_temperature_hessians, - old_old_temperature_hessians, - present_stokes_values, - gamma_values, - global_u_infty, - global_T_range.second - global_T_range.first, - global_Omega_diameter, cell->diameter(), - old_time_step); - - for (unsigned int q=0; q old_grad_T = old_temperature_grads[q]; - const Tensor<1,dim> old_old_grad_T = old_old_temperature_grads[q]; - - - Tensor<1,dim> present_u; - for (unsigned int d=0; dsubdomain_id() == (unsigned int)trilinos_communicator.MyPID() ) + { + local_rhs = 0; + + temperature_fe_values.reinit (cell); + stokes_fe_values.reinit (stokes_cell); + + temperature_fe_values.get_function_values (old_temperature_solution, + old_temperature_values); + temperature_fe_values.get_function_values (old_old_temperature_solution, + old_old_temperature_values); + + temperature_fe_values.get_function_gradients (old_temperature_solution, + old_temperature_grads); + temperature_fe_values.get_function_gradients (old_old_temperature_solution, + old_old_temperature_grads); + + temperature_fe_values.get_function_hessians (old_temperature_solution, + old_temperature_hessians); + temperature_fe_values.get_function_hessians (old_old_temperature_solution, + old_old_temperature_hessians); + + temperature_right_hand_side.value_list (temperature_fe_values.get_quadrature_points(), + gamma_values); + + stokes_fe_values.get_function_values (stokes_solution, + present_stokes_values); + + const double nu + = compute_viscosity (old_temperature_values, + old_old_temperature_values, + old_temperature_grads, + old_old_temperature_grads, + old_temperature_hessians, + old_old_temperature_hessians, + present_stokes_values, + gamma_values, + global_u_infty, + global_T_range.second - global_T_range.first, + global_Omega_diameter, cell->diameter(), + old_time_step); + + for (unsigned int q=0; q old_grad_T = old_temperature_grads[q]; + const Tensor<1,dim> old_old_grad_T = old_old_temperature_grads[q]; + + + Tensor<1,dim> present_u; + for (unsigned int d=0; dget_dof_indices (local_dof_indices); - temperature_constraints.distribute_local_to_global (local_rhs, - local_dof_indices, - temperature_rhs); - } + time_step * + present_u * + ((1+time_step/old_time_step) * old_grad_T + - + time_step / old_time_step * old_old_grad_T) * + phi_T[i] + - + time_step * + nu * + ((1+time_step/old_time_step) * old_grad_T + - + time_step / old_time_step * old_old_grad_T) * + grad_phi_T[i] + + + time_step * + gamma_values[q] * phi_T[i]) + * + temperature_fe_values.JxW(q); + } + else + { + for (unsigned int i=0; iget_dof_indices (local_dof_indices); + temperature_constraints.distribute_local_to_global (local_rhs, + local_dof_indices, + temperature_rhs); + } + + temperature_rhs.compress(); } @@ -1209,7 +1258,7 @@ void BoussinesqFlowProblem::assemble_temperature_system () template void BoussinesqFlowProblem::solve () { - std::cout << " Solving..." << std::endl; + pcout << " Solving..." << std::endl; { LinearSolvers::InverseMatrix::solve () SolverGMRES gmres(solver_control, SolverGMRES::AdditionalData(100)); - //stokes_solution = 0; gmres.solve(stokes_matrix, stokes_solution, stokes_rhs, preconditioner); - std::cout << " " - << solver_control.last_step() - << " GMRES iterations for Stokes subsystem." - << std::endl; - + pcout << " " + << solver_control.last_step() + << " GMRES iterations for Stokes subsystem." + << std::endl; + stokes_constraints.distribute (stokes_solution); } @@ -1261,10 +1309,10 @@ void BoussinesqFlowProblem::solve () temperature_constraints.distribute (temperature_solution); - std::cout << " " - << solver_control.last_step() - << " CG iterations for temperature." - << std::endl; + pcout << " " + << solver_control.last_step() + << " CG iterations for temperature." + << std::endl; double min_temperature = temperature_solution(0), max_temperature = temperature_solution(0); @@ -1276,9 +1324,9 @@ void BoussinesqFlowProblem::solve () temperature_solution(i)); } - std::cout << " Temperature range: " - << min_temperature << ' ' << max_temperature - << std::endl; + pcout << " Temperature range: " + << min_temperature << ' ' << max_temperature + << std::endl; } } @@ -1449,10 +1497,10 @@ void BoussinesqFlowProblem::run () do { - std::cout << "Timestep " << timestep_number - << ": t=" << time - << ", dt=" << time_step - << std::endl; + pcout << "Timestep " << timestep_number + << ": t=" << time + << ", dt=" << time_step + << std::endl; assemble_stokes_system (); build_stokes_preconditioner (); @@ -1462,7 +1510,7 @@ void BoussinesqFlowProblem::run () output_results (); - std::cout << std::endl; + pcout << std::endl; if ((timestep_number == 0) && (pre_refinement_step < n_pre_refinement_steps)) @@ -1495,7 +1543,7 @@ int main (int argc, char *argv[]) (void)argc; (void)argv; #endif - + try { deallog.depth_console (0); @@ -1527,7 +1575,7 @@ int main (int argc, char *argv[]) << std::endl; return 1; } - + #ifdef DEAL_II_COMPILER_SUPPORTS_MPI MPI_Finalize(); #endif