]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Merge Ralf's and my conditional ostream objects.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 29 Jul 2004 15:41:41 +0000 (15:41 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 29 Jul 2004 15:41:41 +0000 (15:41 +0000)
git-svn-id: https://svn.dealii.org/trunk@9535 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/conditional_ostream.h
deal.II/base/source/conditional_ostream.cc
deal.II/examples/step-17/step-17.cc

index 56171dcd2ed220cf75d770d74d102f0683fc08a7..632a2acb3625b8a349f97f04381e084240be32f1 100644 (file)
 
 
 /**
- * A class that allows printing to standard output, i.e. cout,
+ * A class that allows printing to an output stream, e.g. @p std::cout,
  * depending on the ConditionalOStream object being active (default)
  * or not. The condition of this object can be changed by
- * set_condition().
+ * set_condition() and in the constructor.
  *
- * The usual usage of this class is through the pregenerated object
- * <tt>pout</tt> which can be used within parallel computations as
- * follows:
+ * This class is mostly useful in parallel computations. Ordinarily, you would
+ * use @p std::cout to print messages like what the program is presently
+ * doing, or the number of degrees of freedom in each step. However, in
+ * parallel programs, this means that each of the MPI processes write to the
+ * screen, which yields many repetitions of the same text. To avoid it, one
+ * would have to have a designated process, say the one with MPI process
+ * number zero, do the output, and guard each write statement with an
+ * if-condition. This becomes cumbersome and clutters up the code. Rather than
+ * doing so, the present class can be used: objects of its type act just like
+ * a standard output stream, but they only print something based on a
+ * condition that can be set to, for example, <tt>mpi_process==0</tt>, so that
+ * only one process has a true condition and in all other processes writes to
+ * this object just disappear in nirvana.
+ *
+ * The usual usage of this class is as follows:
  *
  * @code
- * pout.set_condition(this_mpi_process==0);
+ * ConditionalOStream pout(this_mpi_process==0);
  *
  *                                  // all processes print following
  *                                  // information to standard output
- * cout << "Reading parameter file on process " << this_mpi_process << endl;
+ * std::cout << "Reading parameter file on process "
+ *           << this_mpi_process << std::endl;
  *
  *                                  // following is printed by
  *                                  // process 0 only
- * pout << "Solving ..." << endl;
+ * pout << "Solving ..." << std::endl;
  * solve();
- * pout << "done" << endl;
+ * pout << "done" << std::endl;
  * @endcode
  *
  * Here, `Reading parameter file on process xy' is printed by each
  *   system_matrix.print_formatted(cout);
  * @endcode
  *
- * @author Ralf Hartmann, 2004
+ * @author Ralf Hartmann, Wolfgang Bangerth, 2004
  */
 class ConditionalOStream
 {
   public:
                                     /**
-                                     * Constructor. Per default the
-                                     * condition of an object is
-                                     * active.
+                                     * Constructor. Set the stream to which
+                                     * we want to write, and the condition
+                                     * based on which writes are actually
+                                     * forwarded. Per default the condition
+                                     * of an object is active.
                                      */
-    ConditionalOStream();
+    ConditionalOStream (std::ostream &stream,
+                        const bool    active = true);
 
                                     /**
                                      * Depending on the
@@ -84,25 +100,41 @@ class ConditionalOStream
                                      * if and only if its condition
                                      * is active.
                                      */
-    void set_condition(bool active);
+    void set_condition (const bool active);
 
                                     /**
-                                     * Give read access to the
-                                     * condition of the object.
+                                     * Return the condition of the object.
                                      */
     bool is_active() const;
 
                                     /**
-                                     * Output a constant something
-                                     * through this stream.
+                                     * Output a constant something through
+                                     * this stream. This function must be @p
+                                     * const so that member objects of this
+                                     * type can also be used from @p const
+                                     * member functions of the surrounding
+                                     * class.
                                      */
     template <typename T>
-    ConditionalOStream & operator << (const T &t);
+    const ConditionalOStream &
+    operator << (const T &t) const;
 
                                     /**
-                                     * Treat ostream manipulators.
+                                     * Treat ostream manipulators. This
+                                     * function must be @p const so that
+                                     * member objects of this type can also
+                                     * be used from @p const member functions
+                                     * of the surrounding class.
+                                     *
+                                     * Note that compilers want to see this
+                                     * treated differently from the general
+                                     * template above since functions like @p
+                                     * std::endl are actually overloaded and
+                                     * can't be bound directly to a template
+                                     * type.
                                      */
-    ConditionalOStream & operator<< (std::ostream& (*p) (std::ostream&));
+    const ConditionalOStream &
+    operator<< (std::ostream& (*p) (std::ostream&)) const;
 
   private:
                                     /**
@@ -110,8 +142,14 @@ class ConditionalOStream
                                      * class could easily be extended
                                      * to treat streams different
                                      * to the standard output.
+                                     *
+                                     * This variable must be @p mutable so
+                                     * that we can write to it in above @p
+                                     * const @p operator<< functions. For the
+                                     * reason why they, in turn, need to be
+                                     * @p const, see there.
                                      */
-    std::ostream  *std_out;
+    mutable std::ostream  &output_stream;
 
                                     /**
                                      * Stores the actual condition
@@ -121,58 +159,29 @@ class ConditionalOStream
 };
 
 
+// --------------------------- inline and template functions -----------
+
 template <class T>
 inline
-ConditionalOStream &
-ConditionalOStream::operator<< (const T& t)
+const ConditionalOStream &
+ConditionalOStream::operator<< (const T& t) const
 {
-  if (active_flag)
-    *std_out << t;
+  if (active_flag == true)
+    output_stream << t;
 
   return *this;
 }
 
 
 inline
-ConditionalOStream &
-ConditionalOStream::operator<< (std::ostream& (*p) (std::ostream&))
+const ConditionalOStream &
+ConditionalOStream::operator<< (std::ostream& (*p) (std::ostream&)) const
 {
-  if (active_flag)
-    *std_out << p;
+  if (active_flag == true)
+    output_stream << p;
 
   return *this;
 }
 
 
-/**
- * Pregenerated object for conditional output to <tt>cout</tt>. Prints
- * to standard output depending on <tt>pout</tt> being active (default)
- * or not. The output can be disabled by set_condition(false).
- *
- * This object is particularly useful in the context of parallel
- * computations in order to avoid multiple but identical output by
- * several processes.
- *
- * @code
- * pout.set_condition(this_mpi_process==0);
- *
- *                                  // all processes print following
- *                                  // information to standard output
- * cout << "Reading parameter file on process " << this_mpi_process << endl;
- *
- *                                  // following is printed by
- *                                  // process 0 only
- * pout << "Solving ..." << endl;
- * solve();
- * pout << "done" << endl;
- * @endcode
- *
- * Here, `Reading parameter file on process xy' is printed by each
- * process separately. In contrast to that, `Solving ...' and `done'
- * is printed to standard output only once, namely by process 0.
- *
- * @author Ralf Hartmann, 2004
- */
-extern ConditionalOStream pout;
-
 #endif
index 8f6c41c253b8ab74dd242b2d31fc5cf81559f420..9fad723a5fc194e453c4fba67a03e294a6452186 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 2004 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
 #include <base/conditional_ostream.h>
 
 
-ConditionalOStream pout;
-
-
-ConditionalOStream::ConditionalOStream():
-               std_out(&std::cout),
-               active_flag(true)
+ConditionalOStream::ConditionalOStream(std::ostream &stream,
+                                       const bool    active)
+                :
+               output_stream (stream),
+               active_flag(active)
 {}
 
 
 void ConditionalOStream::set_condition(bool flag)
 {
-  active_flag=flag;
+  active_flag = flag;
 }
 
 
index 69e8f18ae82979e70af3541738a2748c99b7f9a5..4043f2e8eb4d3bc686f8d4465636c8a99ad48478 100644 (file)
@@ -43,7 +43,7 @@
                                  // program and that weren't in
                                  // step-8. First, we replace the
                                  // standard output ``std::cout'' by a
-                                 // new stream ``pout'' which is used
+                                 // new stream ``pcout'' which is used
                                  // in parallel computations for
                                  // generating output only on one of
                                  // the processes.
                                  // in that we let ``solve'' return a value,
                                  // namely the number of iterations it took to
                                  // converge, so that we can output this to
-                                 // the screen at the appropriate place.
+                                 // the screen at the appropriate place. In
+                                 // addition, we introduce a stream-like
+                                 // variable ``pcout'', explained below:
 template <int dim>
 class ElasticProblem 
 {
@@ -123,6 +125,33 @@ class ElasticProblem
     void refine_grid ();
     void output_results (const unsigned int cycle) const;
 
+                                     // The first variable is basically only
+                                     // for convenience: in parallel program,
+                                     // if each process outputs status
+                                     // information, then there quickly is a
+                                     // lot of clutter. Rather, we would want
+                                     // to only have one process output
+                                     // everything once, for example the one
+                                     // with process number
+                                     // zero. ``ConditionalOStream'' does
+                                     // exactly this: it acts as if it were a
+                                     // stream, but only forwards to a real,
+                                     // underlying stream if a flag is set. By
+                                     // setting this condition to
+                                     // ``this_mpi_process==0'', we make sure
+                                     // that output is only generated from the
+                                     // first process and that we don't get
+                                     // the same lines of output over and over
+                                     // again, once per process.
+                                     //
+                                     // With this simple trick, we make sure
+                                     // that we don't have to guard each and
+                                     // every write to ``std::cout'' by a
+                                     // prefixed ``if(this_mpi_process==0)''.
+    ConditionalOStream pcout;
+
+                                     // The next few variables are taken
+                                     // verbatim from step-8:
     Triangulation<dim>   triangulation;
     DoFHandler<dim>      dof_handler;
 
@@ -366,6 +395,8 @@ ElasticProblem<dim>::get_this_mpi_process (const MPI_Comm &mpi_communicator)
 template <int dim>
 ElasticProblem<dim>::ElasticProblem ()
                 :
+                pcout (std::cout,
+                       get_this_mpi_process(mpi_communicator) == 0),
                dof_handler (triangulation),
                fe (FE_Q<dim>(1), dim),
                 mpi_communicator (MPI_COMM_WORLD),
@@ -1227,31 +1258,20 @@ void ElasticProblem<dim>::refine_grid ()
 
 
 
-                                 // Lastly, here is the driver
-                                 // function. It is almost unchanged
-                                 // from step-8, with the exception
-                                 // that we replace ``std::cout'' by
-                                 // the ``pout'' stream. By setting
-                                 // its condition to
-                                 // ``this_mpi_process==0'', we make
-                                 // sure that output is only generated
-                                 // from the first process and that we
-                                 // don't get the same lines of output
-                                 // over and over again, once per
-                                 // process. Apart from this, the only
-                                 // other cosmetic change is that we
-                                 // output how many degrees of freedom
-                                 // there are per process, and how
-                                 // many iterations it took for the
-                                 // linear solver to converge:
+                                 // Lastly, here is the driver function. It is
+                                 // almost unchanged from step-8, with the
+                                 // exception that we replace ``std::cout'' by
+                                 // the ``pcout'' stream. Apart from this, the
+                                 // only other cosmetic change is that we
+                                 // output how many degrees of freedom there
+                                 // are per process, and how many iterations
+                                 // it took for the linear solver to converge:
 template <int dim>
 void ElasticProblem<dim>::run () 
 {
-  pout.set_condition(this_mpi_process == 0);
-
   for (unsigned int cycle=0; cycle<10; ++cycle)
     {
-      pout << "Cycle " << cycle << ':' << std::endl;
+      pcout << "Cycle " << cycle << ':' << std::endl;
 
       if (cycle == 0)
        {
@@ -1261,27 +1281,27 @@ void ElasticProblem<dim>::run ()
       else
        refine_grid ();
 
-      pout << "   Number of active cells:       "
-          << triangulation.n_active_cells()
-          << std::endl;
+      pcout << "   Number of active cells:       "
+            << triangulation.n_active_cells()
+            << std::endl;
 
       setup_system ();
 
-      pout << "   Number of degrees of freedom: "
-               << dof_handler.n_dofs()
-               << " (by partition:";
+      pcout << "   Number of degrees of freedom: "
+            << dof_handler.n_dofs()
+            << " (by partition:";
       for (unsigned int p=0; p<n_mpi_processes; ++p)
-       pout << (p==0 ? ' ' : '+')
-            << (DoFTools::
-                count_dofs_with_subdomain_association (dof_handler,
-                                                       p));
-      pout << ")" << std::endl;
+       pcout << (p==0 ? ' ' : '+')
+              << (DoFTools::
+                  count_dofs_with_subdomain_association (dof_handler,
+                                                         p));
+      pcout << ")" << std::endl;
       
       assemble_system ();
       const unsigned int n_iterations = solve ();
   
-      pout << "   Solver converged in " << n_iterations
-          << " iterations." << std::endl;
+      pcout << "   Solver converged in " << n_iterations
+            << " iterations." << std::endl;
       
       output_results (cycle);
     }

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.