]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Deprecate the function TimeDependent::end_sweep with an argument indicating the numbe...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 17 Mar 2013 17:55:59 +0000 (17:55 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 17 Mar 2013 17:55:59 +0000 (17:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@28928 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/numerics/time_dependent.h
deal.II/source/numerics/time_dependent.cc
tests/deal.II/time_dependent_01.cc [new file with mode: 0644]
tests/deal.II/time_dependent_01/cmp/generic [new file with mode: 0644]

index 275a9b86296b570446c1d0a7fd1386368659a7a1..0a0277e0aca0b9b15ae1166a69ecb13424ec4495 100644 (file)
@@ -24,11 +24,13 @@ inconvenience this causes.
 </p>
 
 <ol>
-<!--
-<li>
+
+<li> Changed: The TimeDependent::end_sweep function with an argument indicating
+the number of threads has been removed. Use the corresponding function without
+an argument. Since the argument had a default value, few users will have used
+this function.
 <br>
-(NAME, 2013/02/16)
--->
+(Wolfgang Bangerth, 2013/03/17)
 
 </ol>
 
index cf7a73cec9a3c5719bff1068bc30c14dc4d31534..5b34cc9ac12104c99024e2f0d0b5c925c9c55e52 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2009, 2010, 2012 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2009, 2010, 2012, 2013 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -649,11 +649,11 @@ public:
    * function as for the previous
    * one.
    *
-   * This function does not
+   * @note This function does not
    * guarantee that @p end_sweep is
    * called for successive time
-   * steps, rather the order of
-   * time steps for which the
+   * steps successively, rather the order of
+   * time step objects for which the
    * function is called is
    * arbitrary. You should
    * therefore not assume that that
@@ -661,22 +661,19 @@ public:
    * previous time steps
    * already. If in multithread
    * mode, the @p end_sweep function
-   * of several time steps is
+   * of several time steps may be
    * called at once, so you should
    * use synchronization
-   * mechanisms, if your program
+   * mechanisms if your program
    * requires so.
-   *
-   * The parameter denotes the
-   * number of threads that shall
-   * be spawned in parallel. It
-   * defaults to only one thread to
-   * avoid hidden synchronisation
-   * problems, and the value is
-   * ignored if not in multithread
-   * mode.
-   */
-  virtual void end_sweep (const unsigned int n_threads = 1);
+   */
+  virtual void end_sweep ();
+
+  /**
+   * @deprecated Use the function without an argument.
+   * @arg n_threads This argument is ignored.
+   */
+  virtual void end_sweep (const unsigned int n_threads) DEAL_II_DEPRECATED;
 
   /**
    * Determine an estimate for the
index 52d8d4f726ed5e935918e93bc3bf79acb092c58f..8944ad436cadd1804214a72ef5575f0b93bfce1c 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2005, 2006, 2008, 2009, 2010, 2011, 2012 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2005, 2006, 2008, 2009, 2010, 2011, 2012, 2013 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -16,6 +16,7 @@
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/thread_management.h>
 #include <deal.II/base/utilities.h>
+#include <deal.II/base/parallel.h>
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_iterator.h>
@@ -203,30 +204,19 @@ void TimeDependent::start_sweep (const unsigned int s)
 
 
 
-void TimeDependent::end_sweep (const unsigned int n_threads)
+void TimeDependent::end_sweep (const unsigned int)
 {
-#ifdef DEAL_II_WITH_THREADS
-  if (n_threads > 1)
-    {
-      const unsigned int stride = timesteps.size() / n_threads;
-      Threads::ThreadGroup<> threads;
-      void (TimeDependent::*p) (const unsigned int, const unsigned int)
-        = &TimeDependent::end_sweep;
-      for (unsigned int i=0; i<n_threads; ++i)
-        threads += Threads::new_thread (p, *this, i*stride,
-                                        (i == n_threads-1 ?
-                                         timesteps.size() :
-                                         (i+1)*stride));
-      threads.join_all();
-    }
-  else
-    {
-      // now do the work
-      end_sweep (0, timesteps.size());
-    }
-#else
-  end_sweep (0, timesteps.size());
-#endif
+  end_sweep ();
+}
+
+
+void TimeDependent::end_sweep ()
+{
+  void (TimeDependent::*p) (const unsigned int, const unsigned int)
+    = &TimeDependent::end_sweep;
+  parallel::apply_to_subranges (0U, timesteps.size(),
+                               std_cxx1x::bind (p, this, std_cxx1x::_1, std_cxx1x::_2),
+                               1);
 }
 
 
diff --git a/tests/deal.II/time_dependent_01.cc b/tests/deal.II/time_dependent_01.cc
new file mode 100644 (file)
index 0000000..ed54954
--- /dev/null
@@ -0,0 +1,86 @@
+//----------------------------  time_dependent_01.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2013 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  time_dependent_01.cc  ---------------------------
+
+// verify that a change to TimeDependent::end_sweep works as expected:
+// TimeStep::end_sweep must be called once for every time step object
+
+
+#include "../tests.h"
+#include <deal.II/numerics/time_dependent.h>
+
+#include <fstream>
+#include <algorithm>
+#include <cmath>
+
+
+std::ofstream logfile("time_dependent_01/output");
+
+
+std::vector<bool> end_sweep_flags;
+
+
+class TimeStep : public TimeStepBase
+{
+public:
+  TimeStep (const unsigned int time_step_number)
+  :
+  TimeStepBase(0),
+  time_step_number (time_step_number)
+    {}
+
+  virtual void end_sweep ()
+    {
+      static Threads::Mutex mutex;
+      Threads::Mutex::ScopedLock lock(mutex);
+      end_sweep_flags[time_step_number] = true;
+    }
+
+  virtual void solve_primal_problem () {}
+
+private:
+  const unsigned int time_step_number;
+};
+
+
+void test ()
+{
+  // create time steps, more than there are likely threads on current machines
+  TimeDependent td (TimeDependent::TimeSteppingData(0,0),
+                   TimeDependent::TimeSteppingData(0,0),
+                   TimeDependent::TimeSteppingData(0,0));
+  const unsigned int n_time_steps = 10000;
+  for (unsigned int i=0; i<n_time_steps; ++i)
+    td.add_timestep (new TimeStep(i));
+
+  end_sweep_flags = std::vector<bool> (n_time_steps, false);
+  td.end_sweep ();
+
+  // make sure we have called TimeStep::end_sweep once for every time step object
+  Assert (end_sweep_flags == std::vector<bool> (n_time_steps, true),
+         ExcInternalError());
+
+  deallog << "OK" << std::endl;
+}
+
+
+int main ()
+{
+  deallog << std::setprecision(4);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test ();
+
+  return 0;
+}
diff --git a/tests/deal.II/time_dependent_01/cmp/generic b/tests/deal.II/time_dependent_01/cmp/generic
new file mode 100644 (file)
index 0000000..0fd8fc1
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK

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.