]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
changing the functions for computation time
authorfrohne <frohne@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 21 Feb 2013 00:38:41 +0000 (00:38 +0000)
committerfrohne <frohne@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 21 Feb 2013 00:38:41 +0000 (00:38 +0000)
git-svn-id: https://svn.dealii.org/trunk@28502 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-42/step-42.cc

index a24556b2c2db3e0bfbe6d1ba235ce6385dc0ee4e..065acfebeb6260f6c82b63716e92a73b7035eb28 100644 (file)
@@ -137,7 +137,6 @@ namespace Step42
     IndexSet             locally_relevant_dofs;
 
     unsigned int         number_iterations;
-    std::vector<double>  run_time;
 
     ConstraintMatrix     constraints;
     ConstraintMatrix     constraints_hanging_nodes;
@@ -157,7 +156,6 @@ namespace Step42
 
     TrilinosWrappers::PreconditionAMG::AdditionalData additional_data;
     TrilinosWrappers::PreconditionAMG preconditioner_u;
-    TrilinosWrappers::PreconditionAMG preconditioner_t;
 
     std::unique_ptr<Input<dim> >               input_obstacle;
     std::unique_ptr<ConstitutiveLaw<dim> >     plast_lin_hard;
@@ -166,6 +164,8 @@ namespace Step42
     double gamma;      // Parameter for the linear isotropic hardening
     double e_modul;    // E-Modul
     double nu;         // Poisson ratio
+
+    TimerOutput          computing_timer;
   };
 
   // As explained above this class
@@ -414,10 +414,8 @@ namespace Step42
         yield = 0;
         stress_strain_tensor = stress_strain_tensor_mu;
         double beta = 1.0;
-        if (deviator_stress_tensor_norm >
-            sigma_0)
+        if (deviator_stress_tensor_norm > sigma_0)
           {
-//            beta = (sigma_0 + gamma/(2*mu + gamma)*(deviator_stress_tensor_norm - sigma_0))/deviator_stress_tensor_norm;
             beta = sigma_0/deviator_stress_tensor_norm;
             stress_strain_tensor *= (gamma + (1 - gamma)*beta);
             yield = 1;
@@ -447,19 +445,13 @@ namespace Step42
         stress_strain_tensor = stress_strain_tensor_mu;
         stress_strain_tensor_linearized = stress_strain_tensor_mu;
         double beta = 1.0;
-        if (deviator_stress_tensor_norm >
-            sigma_0)
+        if (deviator_stress_tensor_norm > sigma_0)
           {
-//            beta = (sigma_0 + gamma/(2*mu + gamma)*(deviator_stress_tensor_norm - sigma_0))/deviator_stress_tensor_norm;
             beta = sigma_0/deviator_stress_tensor_norm;
             stress_strain_tensor *= (gamma + (1 - gamma)*beta);
             stress_strain_tensor_linearized *= (gamma + (1 - gamma)*beta);
             deviator_stress_tensor /= deviator_stress_tensor_norm;
             stress_strain_tensor_linearized -= (1 - gamma)*beta*2*mu*outer_product(deviator_stress_tensor, deviator_stress_tensor);
-
-//            std::cout << "Plastisch mit beta = " << beta
-//                << ", sigma_0 = " << sigma_0
-//                << std::endl;
           }
 
         stress_strain_tensor += stress_strain_tensor_kappa;
@@ -615,7 +607,10 @@ namespace Step42
     sigma_0 (400),
     gamma (0.01),
     e_modul (2.0e+5),
-    nu (0.3)
+    nu (0.3),
+    computing_timer (pcout,
+                     TimerOutput::summary,
+                     TimerOutput::wall_times)
   {
     // double _E, double _nu, double _sigma_0, double _gamma
     plast_lin_hard.reset (new ConstitutiveLaw<dim> (e_modul, nu, sigma_0, gamma, mpi_communicator, pcout));
@@ -739,6 +734,8 @@ namespace Step42
   template <int dim>
   void PlasticityContactProblem<dim>::assemble_nl_system (TrilinosWrappers::MPI::Vector &u)
   {
+    computing_timer.enter_section("Assembling");
+
     QGauss<dim>  quadrature_formula(2);
     QGauss<dim-1>  face_quadrature_formula(2);
 
@@ -856,6 +853,8 @@ namespace Step42
 
     system_matrix_newton.compress ();
     system_rhs_newton.compress (VectorOperation::add);
+
+    computing_timer.exit_section("Assembling");
   }
 
   template <int dim>
@@ -1034,7 +1033,7 @@ namespace Step42
   template <int dim>
   void PlasticityContactProblem<dim>::update_solution_and_constraints ()
   {
-    clock_t                        start_proj, end_proj;
+    computing_timer.enter_section("Update solution and constraints");
 
     const EquationData::Obstacle<dim>     obstacle (input_obstacle);
     std::vector<bool>                     vertex_touched (dof_handler.n_dofs (), false);
@@ -1115,6 +1114,8 @@ namespace Step42
     constraints.close ();
 
     constraints.merge (constraints_dirichlet_hanging_nodes);
+
+    computing_timer.exit_section("Update solution and constraints");
   }
 
   // @sect4{PlasticityContactProblem::dirichlet_constraints}
@@ -1190,7 +1191,7 @@ namespace Step42
   template <int dim>
   void PlasticityContactProblem<dim>::solve ()
   {
-    Timer t;
+    computing_timer.enter_section ("Solving and Preconditioning");
 
     TrilinosWrappers::MPI::Vector    distributed_solution (system_rhs_newton);
     distributed_solution = solution;
@@ -1200,18 +1201,13 @@ namespace Step42
     distributed_solution.compress(VectorOperation::insert);
     system_rhs_newton.compress(VectorOperation::insert);
 
-    MPI_Barrier (mpi_communicator);
-    t.restart();
+    computing_timer.enter_section("Preconditioning");
 
     preconditioner_u.initialize (system_matrix_newton, additional_data);
 
-    MPI_Barrier (mpi_communicator);
-    t.stop();
-    if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
-      run_time[6] += t.wall_time();
+    computing_timer.exit_section("Preconditioning");
 
-    MPI_Barrier (mpi_communicator);
-    t.restart();
+    computing_timer.enter_section("Solving");
 
     PrimitiveVectorMemory<TrilinosWrappers::MPI::Vector> mem;
     TrilinosWrappers::MPI::Vector    tmp (system_rhs_newton);
@@ -1231,16 +1227,15 @@ namespace Step42
           << " FGMRES iterations."
           << std::endl;
 
-    MPI_Barrier (mpi_communicator);
-    t.stop();
-    if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
-      run_time[7] += t.wall_time();
+    computing_timer.exit_section("Solving");
 
     number_iterations += solver_control.last_step();
 
     constraints.distribute (distributed_solution);
 
     solution = distributed_solution;
+
+    computing_timer.exit_section("Solving and Preconditioning");
   }
 
 
@@ -1252,7 +1247,6 @@ namespace Step42
     double                         resid_old=100000;
     TrilinosWrappers::MPI::Vector  res (system_rhs_newton);
     TrilinosWrappers::MPI::Vector  tmp_vector (system_rhs_newton);
-    Timer                          t;
 
     std::vector<std::vector<bool> > constant_modes;
     DoFTools::extract_constant_modes (dof_handler,
@@ -1283,37 +1277,17 @@ namespace Step42
         pcout<< "   Newton iteration " << j <<std::endl;
         pcout<< "      Updating active set..." <<std::endl;
 
-        MPI_Barrier (mpi_communicator);
-        t.restart();
-
         update_solution_and_constraints ();
 
-        MPI_Barrier (mpi_communicator);
-        t.stop();
-        if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
-          run_time[5] += t.wall_time();
-
         pcout<< "      Assembling system... " <<std::endl;
-        MPI_Barrier (mpi_communicator);
-        t.restart();
         system_matrix_newton = 0;
         system_rhs_newton = 0;
         assemble_nl_system (solution);  //compute Newton-Matrix
-        MPI_Barrier (mpi_communicator);
-        t.stop();
-        if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
-          run_time[1] += t.wall_time();
 
         number_assemble_system += 1;
 
-        MPI_Barrier (mpi_communicator);
-        t.restart();
         pcout<< "      Solving system... " <<std::endl;
         solve ();
-        MPI_Barrier (mpi_communicator);
-        t.stop();
-        if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
-          run_time[2] += t.wall_time();
 
         TrilinosWrappers::MPI::Vector    distributed_solution (system_rhs_newton);
         distributed_solution = solution;
@@ -1330,8 +1304,8 @@ namespace Step42
 //            constraints.distribute (old_solution);
             old_solution.compress (VectorOperation::add);
 
-            MPI_Barrier (mpi_communicator);
-            t.restart();
+            computing_timer.enter_section("Residual and lambda");
+
             system_rhs_newton = 0;
 
             solution = old_solution;
@@ -1352,10 +1326,7 @@ namespace Step42
             if (resid<resid_old)
               damped=1;
 
-            MPI_Barrier (mpi_communicator);
-            t.stop();
-            if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
-              run_time[3] += t.wall_time();
+            computing_timer.exit_section("Residual and lambda");
 
             pcout << "      Residual of the non-contact part of the system: " << resid
                               << std::endl
@@ -1509,18 +1480,14 @@ namespace Step42
     input_obstacle.reset (new Input<dim>("obstacle_file.pbm"));
     pcout << "Ostacle is available now." << std::endl;
 
-    Timer             t;
-    run_time.resize (8);
-
     const unsigned int n_cycles = 6;
     for (unsigned int cycle=0; cycle<n_cycles; ++cycle)
       {
+        computing_timer.enter_section("Mesh refinement and setup system");
+
         pcout << "" <<std::endl;
         pcout << "Cycle " << cycle << ':' << std::endl;
 
-        MPI_Barrier (mpi_communicator);
-        t.restart();
-
         if (cycle == 0)
           {
             make_grid();
@@ -1530,33 +1497,21 @@ namespace Step42
 
         setup_system ();
 
-        MPI_Barrier (mpi_communicator);
-        t.stop();
-        run_time[0] += t.wall_time();;
+        computing_timer.exit_section("Mesh refinement and setup system");
 
         solve_newton ();
 
         pcout<< "      Writing graphical output..." <<std::endl;
-        MPI_Barrier (mpi_communicator);
-        t.restart();
+        computing_timer.enter_section("Graphical output");
+
         std::ostringstream filename_solution;
         filename_solution << "solution-";
         filename_solution << cycle;
         output_results (filename_solution.str ());
-        MPI_Barrier (mpi_communicator);
-        t.stop();
-        if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
-          run_time[4] += t.wall_time();
-
-        pcout << "      Computing time for:" << std::endl
-              << "         making grid and setup = " << run_time[0] << std::endl
-              << "         updating active set = " << run_time[5] <<std::endl
-              << "         assembling system = " << run_time[1] <<std::endl
-              << "         solving system = " << run_time[2] <<std::endl
-              << "         preconditioning = " << run_time[6] <<std::endl
-              << "         solving with FGMRES = " << run_time[7] <<std::endl
-              << "         computing error and lambda = " << run_time[3] <<std::endl
-              << "         writing graphical output = " << run_time[4] <<std::endl;
+
+        computing_timer.exit_section("Graphical output");
+
+        computing_timer.print_summary();
       }
   }
 }

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.