]> https://gitweb.dealii.org/ - dealii.git/commitdiff
step-32: correct comment and cleanup calculation of get_extrapolated_temperature_range().
authorTimo Heister <timo.heister@gmail.com>
Wed, 26 Jan 2011 08:02:09 +0000 (08:02 +0000)
committerTimo Heister <timo.heister@gmail.com>
Wed, 26 Jan 2011 08:02:09 +0000 (08:02 +0000)
git-svn-id: https://svn.dealii.org/trunk@23263 0785d39b-7218-0410-832d-ea1e28bc413d

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

index ad34bce157b24ed445452b34baf8c8466da230b7..78c0229f2f1577bb3bc284b6170c99582f28acf0 100644 (file)
@@ -69,6 +69,7 @@
 
 #include <fstream>
 #include <iostream>
+#include <limits>
 
                                 // This is the only include file that is new:
                                 // We use an IndexSet to describe the
@@ -1133,15 +1134,17 @@ BoussinesqFlowProblem<dim>::get_extrapolated_temperature_range () const
   std::vector<double> old_temperature_values(n_q_points);
   std::vector<double> old_old_temperature_values(n_q_points);
 
+                                  // This presets the minimum with a bigger
+                                  // and the maximum with a smaller number
+                                  // than one that is going to appear. Will
+                                  // be overwritten in the cell loop or in
+                                  // the communication step at the
+                                  // latest.
+  double min_local_temperature = std::numeric_limits<double>::max(),
+        max_local_temperature = -std::numeric_limits<double>::max();
+
   if (timestep_number != 0)
     {
-                                      // use different formula here than in
-                                      // step-31 to avoid computing the
-                                      // linfty_norm of the temperature
-                                      // (which requires a communication)
-      double min_local_temperature = 1e30,
-            max_local_temperature = -1e30;
-
       typename DoFHandler<dim>::active_cell_iterator
        cell = temperature_dof_handler.begin_active(),
        endc = temperature_dof_handler.end();
@@ -1167,25 +1170,9 @@ BoussinesqFlowProblem<dim>::get_extrapolated_temperature_range () const
                                                  temperature);
              }
          }
-
-      double min_temperature, max_temperature;
-#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
-      MPI_Allreduce (&max_local_temperature, &max_temperature, 1, MPI_DOUBLE,
-                    MPI_MAX, MPI_COMM_WORLD);
-      MPI_Allreduce (&min_local_temperature, &min_temperature, 1, MPI_DOUBLE,
-                    MPI_MIN, MPI_COMM_WORLD);
-#else
-      min_temperature = min_local_temperature;
-      max_temperature = max_local_temperature;
-#endif
-
-      return std::make_pair(min_temperature, max_temperature);
     }
   else
     {
-      double min_local_temperature = 1e30,
-            max_local_temperature = -1e30;      
-
       typename DoFHandler<dim>::active_cell_iterator
        cell = temperature_dof_handler.begin_active(),
        endc = temperature_dof_handler.end();
@@ -1207,20 +1194,20 @@ BoussinesqFlowProblem<dim>::get_extrapolated_temperature_range () const
                                                  temperature);
              }
          }
-
-      double min_temperature, max_temperature;
+    }
+  
+  double min_temperature, max_temperature;
 #ifdef DEAL_II_COMPILER_SUPPORTS_MPI
-      MPI_Allreduce (&max_local_temperature, &max_temperature, 1, MPI_DOUBLE,
-                    MPI_MAX, MPI_COMM_WORLD);
-      MPI_Allreduce (&min_local_temperature, &min_temperature, 1, MPI_DOUBLE,
-                    MPI_MIN, MPI_COMM_WORLD);
+  MPI_Allreduce (&max_local_temperature, &max_temperature, 1, MPI_DOUBLE,
+                MPI_MAX, MPI_COMM_WORLD);
+  MPI_Allreduce (&min_local_temperature, &min_temperature, 1, MPI_DOUBLE,
+                MPI_MIN, MPI_COMM_WORLD);
 #else
-      min_temperature = min_local_temperature;
-      max_temperature = max_local_temperature;
+  min_temperature = min_local_temperature;
+  max_temperature = max_local_temperature;
 #endif
 
-      return std::make_pair(min_temperature, max_temperature);
-    }
+  return std::make_pair(min_temperature, max_temperature);
 }
 
 

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.