]> https://gitweb.dealii.org/ - dealii.git/commitdiff
step-32: calculate and output l2 divergence error for testing
authorTimo Heister <timo.heister@gmail.com>
Tue, 22 Feb 2011 10:32:28 +0000 (10:32 +0000)
committerTimo Heister <timo.heister@gmail.com>
Tue, 22 Feb 2011 10:32:28 +0000 (10:32 +0000)
git-svn-id: https://svn.dealii.org/trunk@23411 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 228afda708794135e10bceca51f2d58f62bf4bb0..f05dbd6b6905ec12bae699cb44f2a1a78b683f36 100644 (file)
@@ -3207,6 +3207,57 @@ void BoussinesqFlowProblem<dim>::output_results ()
 {
   computing_timer.enter_section ("Postprocessing");
 
+                                  //calculate l2 divergence error
+  {
+    double my_cells_error = 0.0;
+    QGauss<1>      q_base(stokes_degree+1);
+    QIterated<dim> err_quadrature(q_base, 2);
+    
+    const unsigned int n_q_points =  err_quadrature.size();
+    FEValues<dim> fe_values (stokes_fe,  err_quadrature, update_JxW_values | update_gradients);
+    const unsigned int dofs_per_cell = fe_values.get_fe().dofs_per_cell;
+    const FEValuesExtractors::Vector velocities (0);
+
+    std::vector<unsigned int> local_dof_indices (fe_values.dofs_per_cell);
+    
+    
+    typename DoFHandler<dim>::active_cell_iterator
+      cell = stokes_dof_handler.begin_active(),
+      endc = stokes_dof_handler.end();
+    for (; cell!=endc; ++cell)
+      if (cell->subdomain_id() == 
+         Utilities::System::get_this_mpi_process(MPI_COMM_WORLD))
+       {
+         fe_values.reinit (cell);
+         cell->get_dof_indices(local_dof_indices);
+         
+         double cell_error = 0.0;
+         for (unsigned int q = 0; q < n_q_points; ++q)
+           {
+             double tmp = 0;
+             for (unsigned int i=0; i < dofs_per_cell; ++i)
+               {
+                 tmp +=
+                   fe_values[velocities].divergence(i, q)
+                   * stokes_solution(local_dof_indices[i]);
+               }
+             cell_error += tmp * tmp * fe_values.JxW(q);
+           }
+         my_cells_error += cell_error; // = sqrt(cell_error)^2
+       }
+
+  double div_error = 0.0;
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+  MPI_Allreduce (&my_cells_error, &div_error, 1, MPI_DOUBLE,
+                MPI_SUM, MPI_COMM_WORLD);
+#else
+  div_error = my_cells_error;
+#endif
+
+  div_error = std::sqrt(div_error);
+  pcout << "> divergence error:" << div_error << std::endl;
+  }
+
   const FESystem<dim> joint_fe (stokes_fe, 1,
                                 temperature_fe, 1);
 

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.