From: bangerth Date: Thu, 9 Feb 2012 16:30:49 +0000 (+0000) Subject: Some steps to clean up the program. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=6cdeae8c69d6f9956f037378d6017b53182acdc7;p=dealii-svn.git Some steps to clean up the program. git-svn-id: https://svn.dealii.org/trunk@25022 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/examples/step-41/step-41.cc b/deal.II/examples/step-41/step-41.cc index 911404cf08..df1c1b7d96 100644 --- a/deal.II/examples/step-41/step-41.cc +++ b/deal.II/examples/step-41/step-41.cc @@ -68,22 +68,18 @@ namespace Step41 void assemble_mass_matrix (); void projection_active_set (); void solve (); - void output_results (const std::string& title) const; + void output_results (const unsigned int iteration) const; Triangulation triangulation; FE_Q fe; DoFHandler dof_handler; - unsigned int n_refinements; - ConstraintMatrix constraints; - SparsityPattern sparsity_pattern; TrilinosWrappers::SparseMatrix system_matrix; TrilinosWrappers::SparseMatrix system_matrix_complete; TrilinosWrappers::SparseMatrix mass_matrix; TrilinosWrappers::Vector solution; - TrilinosWrappers::Vector tmp_solution; TrilinosWrappers::Vector system_rhs; TrilinosWrappers::Vector system_rhs_complete; TrilinosWrappers::Vector resid_vector; @@ -193,15 +189,14 @@ namespace Step41 void ObstacleProblem::make_grid () { GridGenerator::hyper_cube (triangulation, -1, 1); - n_refinements = 8; - triangulation.refine_global (n_refinements); + triangulation.refine_global (7); - std::cout << " Number of active cells: " + std::cout << "Number of active cells: " << triangulation.n_active_cells() << std::endl - << " Total number of cells: " + << "Total number of cells: " << triangulation.n_cells() - << std::endl; + << std::endl; } // @sect4{ObstacleProblem::setup_system} @@ -211,20 +206,19 @@ namespace Step41 { dof_handler.distribute_dofs (fe); - std::cout << " Number of degrees of freedom: " + std::cout << "Number of degrees of freedom: " << dof_handler.n_dofs() + << std::endl << std::endl; CompressedSparsityPattern c_sparsity(dof_handler.n_dofs()); DoFTools::make_sparsity_pattern (dof_handler, c_sparsity, constraints, false); - sparsity_pattern.copy_from(c_sparsity); - system_matrix.reinit (sparsity_pattern); - system_matrix_complete.reinit (sparsity_pattern); - mass_matrix.reinit (sparsity_pattern); + system_matrix.reinit (c_sparsity); + system_matrix_complete.reinit (c_sparsity); + mass_matrix.reinit (c_sparsity); solution.reinit (dof_handler.n_dofs()); - tmp_solution.reinit (dof_handler.n_dofs()); system_rhs.reinit (dof_handler.n_dofs()); system_rhs_complete.reinit (dof_handler.n_dofs()); resid_vector.reinit (dof_handler.n_dofs()); @@ -246,6 +240,8 @@ namespace Step41 template void ObstacleProblem::assemble_system () { + std::cout << " Assembling system..." << std::endl; + QGauss quadrature_formula(2); const RightHandSide right_hand_side; @@ -357,6 +353,8 @@ namespace Step41 template void ObstacleProblem::projection_active_set () { + std::cout << " Updating active set..." << std::endl; + const Obstacle obstacle; std::vector vertex_touched (triangulation.n_vertices(), false); @@ -403,7 +401,8 @@ namespace Step41 } } } - std::cout<< "Number of Contact-Constaints: " << counter_contact_constraints < void ObstacleProblem::solve () { + std::cout << " Solving system..." << std::endl; + ReductionControl reduction_control (100, 1e-12, 1e-3); SolverCG solver (reduction_control); TrilinosWrappers::PreconditionAMG precondition; precondition.initialize (system_matrix); solver.solve (system_matrix, solution, system_rhs, precondition); + constraints.distribute (solution); - std::cout << "Initial error: " << reduction_control.initial_value() < " << reduction_control.last_value() + << " in " + << reduction_control.last_step() + << " CG iterations." << std::endl; - - constraints.distribute (solution); } // @sect4{ObstacleProblem::output_results} @@ -441,18 +442,22 @@ namespace Step41 // The file contains the displacement, // the residual and active set vectors. template - void ObstacleProblem::output_results (const std::string& title) const + void ObstacleProblem::output_results (const unsigned int iteration) const { + std::cout << " Writing graphical output..." << std::endl; + DataOut data_out; data_out.attach_dof_handler (dof_handler); - data_out.add_data_vector (tmp_solution, "Displacement"); - data_out.add_data_vector (resid_vector, "Residual"); - data_out.add_data_vector (active_set, "ActiveSet"); + data_out.add_data_vector (solution, "displacement"); + data_out.add_data_vector (resid_vector, "residual"); + data_out.add_data_vector (active_set, "active_set"); data_out.build_patches (); - std::ofstream output_vtk ((title + ".vtk").c_str ()); + std::ofstream output_vtk ((std::string("output_") + + Utilities::int_to_string (iteration) + + ".vtk").c_str ()); data_out.write_vtk (output_vtk); } @@ -469,11 +474,14 @@ namespace Step41 template void ObstacleProblem::run () { - std::cout << "Solving problem in " << dim << " space dimensions." << std::endl; - make_grid(); setup_system (); + // TODO: can't some of this be + // merged with the first Newton + // iteration? + std::cout << "Initial start-up step" << std::endl; + constraints.clear (); VectorTools::interpolate_boundary_values (dof_handler, 0, @@ -499,48 +507,47 @@ namespace Step41 for (unsigned int j=0; j