]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Cleanup ::SolutionTransfer.
authorMarc Fehling <marc.fehling@gmx.net>
Mon, 8 Feb 2021 23:11:58 +0000 (16:11 -0700)
committerMarc Fehling <marc.fehling@gmx.net>
Thu, 11 Feb 2021 18:32:32 +0000 (11:32 -0700)
source/numerics/solution_transfer.cc

index 585c15f3dd7b2f94b06e629866661df50693aeba..80d4fc48570aef8da136fe6e514a8877336c5b5b 100644 (file)
@@ -100,12 +100,9 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::prepare_for_pure_refinement()
   std::vector<std::vector<types::global_dof_index>>(n_active_cells)
     .swap(indices_on_cell);
 
-  typename DoFHandlerType::active_cell_iterator cell =
-                                                  dof_handler->begin_active(),
-                                                endc = dof_handler->end();
-
-  for (unsigned int i = 0; cell != endc; ++cell, ++i)
+  for (const auto &cell : dof_handler->active_cell_iterators())
     {
+      const unsigned int i = cell->active_cell_index();
       indices_on_cell[i].resize(cell->get_fe().n_dofs_per_cell());
       // on each cell store the indices of the
       // dofs. after refining we get the values
@@ -138,14 +135,11 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::refine_interpolate(
 
   Vector<typename VectorType::value_type> local_values(0);
 
-  typename DoFHandlerType::cell_iterator cell = dof_handler->begin(),
-                                         endc = dof_handler->end();
-
   typename std::map<std::pair<unsigned int, unsigned int>,
                     Pointerstruct>::const_iterator pointerstruct,
     cell_map_end = cell_map.end();
 
-  for (; cell != endc; ++cell)
+  for (const auto &cell : dof_handler->cell_iterators())
     {
       pointerstruct =
         cell_map.find(std::make_pair(cell->level(), cell->index()));
@@ -265,23 +259,20 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::
   Assert(prepared_for != coarsening_and_refinement,
          ExcAlreadyPrepForCoarseAndRef());
 
+  clear();
+  n_dofs_old                 = dof_handler->n_dofs();
   const unsigned int in_size = all_in.size();
+
+#ifdef DEBUG
   Assert(in_size != 0,
          ExcMessage("The array of input vectors you pass to this "
                     "function has no elements. This is not useful."));
-
-  clear();
-
-  const unsigned int n_active_cells =
-    dof_handler->get_triangulation().n_active_cells();
-  (void)n_active_cells;
-  n_dofs_old = dof_handler->n_dofs();
-
   for (unsigned int i = 0; i < in_size; ++i)
     {
       Assert(all_in[i].size() == n_dofs_old,
              ExcDimensionMismatch(all_in[i].size(), n_dofs_old));
     }
+#endif
 
   // first count the number
   // of cells that will be coarsened
@@ -295,13 +286,12 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::
       else
         ++n_cells_to_stay_or_refine;
     }
-  Assert((n_cells_to_coarsen + n_cells_to_stay_or_refine) == n_active_cells,
+  Assert((n_cells_to_coarsen + n_cells_to_stay_or_refine) ==
+           dof_handler->get_triangulation().n_active_cells(),
          ExcInternalError());
 
   unsigned int n_coarsen_fathers = 0;
-  for (typename DoFHandlerType::cell_iterator cell = dof_handler->begin();
-       cell != dof_handler->end();
-       ++cell)
+  for (const auto &cell : dof_handler->cell_iterators())
     if (!cell->is_active() && cell->child(0)->coarsen_flag_set())
       ++n_coarsen_fathers;
   Assert(n_cells_to_coarsen >= 2 * n_coarsen_fathers, ExcInternalError());
@@ -328,9 +318,7 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::
   // the 'to_stay_or_refine' cells 'n_sr' and
   // the 'coarsen_fathers' cells 'n_cf',
   unsigned int n_sr = 0, n_cf = 0;
-  for (typename DoFHandlerType::cell_iterator cell = dof_handler->begin();
-       cell != dof_handler->end();
-       ++cell)
+  for (const auto &cell : dof_handler->cell_iterators())
     {
       // CASE 1: active cell that remains as it is
       if (cell->is_active() && !cell->coarsen_flag_set())
@@ -424,8 +412,9 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::interpolate(
   const std::vector<VectorType> &all_in,
   std::vector<VectorType> &      all_out) const
 {
-  Assert(prepared_for == coarsening_and_refinement, ExcNotPrepared());
   const unsigned int size = all_in.size();
+#ifdef DEBUG
+  Assert(prepared_for == coarsening_and_refinement, ExcNotPrepared());
   Assert(all_out.size() == size, ExcDimensionMismatch(all_out.size(), size));
   for (unsigned int i = 0; i < size; ++i)
     Assert(all_in[i].size() == n_dofs_old,
@@ -438,6 +427,7 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::interpolate(
       Assert(&all_in[i] != &all_out[j],
              ExcMessage("Vectors cannot be used as input and output"
                         " at the same time!"));
+#endif
 
   Vector<typename VectorType::value_type> local_values;
   std::vector<types::global_dof_index>    dofs;

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.