]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Update to auto-detection of the number of CPUs in a system.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 24 Feb 2000 13:13:32 +0000 (13:13 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 24 Feb 2000 13:13:32 +0000 (13:13 +0000)
git-svn-id: https://svn.dealii.org/trunk@2491 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/data_out.h
deal.II/deal.II/include/numerics/error_estimator.h
deal.II/deal.II/include/numerics/time_dependent.h
deal.II/deal.II/source/numerics/data_out.cc
deal.II/deal.II/source/numerics/error_estimator.cc

index 237c187f06dbf8417c7e71fb05403421cf39f51d..e3a9329b07dfafe65707d27b8d400263df8763fd 100644 (file)
@@ -17,6 +17,8 @@
 #include <lac/forward_declarations.h>
 #include <grid/forward_declarations.h>
 #include <base/data_out_base.h>
+#include <base/multithread_info.h>
+
 
 
 /**
@@ -391,13 +393,18 @@ class DataOut : public DataOut_DoFData<dim>
                                      * general documentation of this class
                                      * for further information.
                                      *
-                                     * The functions supports multithreading.
+                                     * The function supports
+                                     * multithreading, if deal.II is
+                                     * compiled in multithreading
+                                     * mode. The default number of
+                                     * threads to be used to build
+                                     * the patches is set to
+                                     * #multithread_info.n_default_threads#.
                                      */
     virtual void build_patches (const unsigned int n_subdivisions = 1,
-                               const unsigned int n_threads_     = 1);
-
+                               const unsigned int n_threads      = multithread_info.n_default_threads);
 
-/**
+                                    /**
                                      * Return the first cell which we
                                      * want output for. The default
                                      * implementation returns the first
@@ -453,14 +460,14 @@ class DataOut : public DataOut_DoFData<dim>
          {}
     };
                                     /**
-                                     * Builds every #n_threads's#
+                                     * Builds every #n_threads#'s
                                      * patch. This function may be
-                                     * called parallel.
+                                     * called in parallel.
                                      * If multithreading is not
                                      * used, the function is called
                                      * once and generates all patches.
                                      */
-    void DataOut<dim>::build_some_patches (Data data);
+    void DataOut<dim>::build_some_patches (Data data);
 };
 
 
index f2ab12e961f9980222b00622e9d4c2cf24861d81..fe7070207c0eb48ede2b2444d999011ed364dad3 100644 (file)
 
 #include <base/exceptions.h>
 #include <base/function.h>
+#include <base/multithread_info.h>
 #include <grid/forward_declarations.h>
 #include <dofs/dof_handler.h>
 #include <map>
 
-                                // if multithreaded set number
-                                // of threads to use.
-                                // The default number of threads can
-                                // be set during the compiling with
-                                // the flag N_THREADS.
-#ifdef DEAL_II_USE_MT
- #ifndef N_THREADS
-  #define N_THREADS 4
- #endif
-#else
- #ifdef N_THREADS
-  #undef N_THREADS
- #endif
- #define N_THREADS 1
-#endif
 
 
 /**
@@ -247,23 +233,25 @@ class KellyErrorEstimator
                                      *
                                      * The estimator supports
                                      * multithreading and splits the
-                                     * cells to #N_THREADS# (default)
-                                     * threads. The number of threads
-                                     * to be used in multithreaded
-                                     * mode can be set with the last
-                                     * parameter of the error
-                                     * estimator.  Multithreading is
-                                     * only implemented in two and
-                                     * three dimensions.
+                                     * cells to
+                                     * #multithread_info.n_default_threads#
+                                     * (default) threads. The number
+                                     * of threads to be used in
+                                     * multithreaded mode can be set
+                                     * with the last parameter of the
+                                     * error estimator.
+                                     * Multithreading is only
+                                     * implemented in two and three
+                                     * dimensions.
                                      */
     static void estimate (const DoFHandler<dim>   &dof,
                          const Quadrature<dim-1> &quadrature,
                          const FunctionMap       &neumann_bc,
                          const Vector<double>    &solution,
                          Vector<float>           &error,
-                         const vector<bool>      &component_mask_ = vector<bool>(),
+                         const vector<bool>      &component_mask = vector<bool>(),
                          const Function<dim>     *coefficients   = 0,
-                         unsigned int            n_threads=N_THREADS);
+                         unsigned int            n_threads = multithread_info.n_default_threads);
     
                                     /**
                                      * Exception
@@ -332,7 +320,7 @@ private:
                                      * is done by the C++ library and
                                      * has not to be handcoded, it
                                      * nevertheless seriously damages
-                                     * tha ability to efficiently run
+                                     * the ability to efficiently run
                                      * the functions of this class in
                                      * parallel, since they are quite
                                      * often blocked by these
@@ -422,7 +410,7 @@ private:
     };
 
 
-/**
+                                    /**
                                      * Computates the error on all cells
                                      * of the domain with the number n,
                                      * satisfying
@@ -458,9 +446,7 @@ private:
                                      * to avoid ending up with a function
                                      * of 500 lines of code.
                                      */
-
-
-static void integrate_over_regular_face (Data                       &data,
+    static void integrate_over_regular_face (Data                       &data,
                                             const unsigned int          this_thread,
                                             const active_cell_iterator &cell,
                                             const unsigned int          face_no,
@@ -468,7 +454,7 @@ static void integrate_over_regular_face (Data                       &data,
                                             FEFaceValues<dim>          &fe_face_values_neighbor);
 
 
-/**
+                                    /**
                                      * The same applies as for the function
                                      * above, except that integration is
                                      * over face #face_no# of #cell#, where
index a2f700ea5580f0f6b5e775a328ece3fb75669700..25a4f393340ffbba01dbe05eec8786ce86bfba2b 100644 (file)
@@ -663,9 +663,11 @@ class TimeDependent
                                      * The parameter denotes the
                                      * number of threads that shall
                                      * be spawned in parallel. It
-                                     * defaults to only one thread,
-                                     * and the value is ignored if
-                                     * not in multithread mode.
+                                     * 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);
 
@@ -728,6 +730,7 @@ class TimeDependent
 };
 
 
+
 /**
  * Base class for a time step in time dependent problems. This class provides
  * barely more than the basic framework, defining the necessary virtual
index 3974a81bed90d731da3de4e92c661ac51cccd3e1..55a14aeb3915df286a15a3b263365de9850ad14e 100644 (file)
@@ -146,6 +146,7 @@ vector<string> DataOut_DoFData<dim>::get_dataset_names () const
 };
 
 
+
 template <int dim>
 const vector<typename DataOutBase::Patch<dim> > &
 DataOut_DoFData<dim>::get_patches () const
@@ -154,8 +155,9 @@ DataOut_DoFData<dim>::get_patches () const
 };
 
 
+
 template <int dim>
-void DataOut<dim>::build_some_patches (Data data)
+void DataOut<dim>::build_some_patches (Data data)
 {
   QTrapez<1>     q_trapez;
   QIterated<dim> patch_points (q_trapez, data.n_subdivisions);
@@ -224,8 +226,6 @@ void * DataOut<dim>::build_some_patches (Data data)
           ++i,++patch,++cell_number,cell=next_cell(cell));
       
     };
-  
-  return 0;
 }
 
 
@@ -233,14 +233,12 @@ template <int dim>
 void DataOut<dim>::build_patches (const unsigned int n_subdivisions,
                                  const unsigned int n_threads_) 
 {
-  unsigned int n_threads = n_threads_;
-
   Assert (dofs != 0, ExcNoDoFHandlerSelected());
 
-                                  // if not in multithread mode
-                                  // set nr of threads to one
-#ifndef DEAL_II_USE_MT
-  n_threads = 1;
+#ifdef DEAL_II_USE_MT
+  const unsigned int n_threads = n_threads_;
+#else
+  const unsigned int n_threads = 1;
 #endif
 
 
@@ -301,7 +299,7 @@ void DataOut<dim>::build_patches (const unsigned int n_subdivisions,
   ThreadManager thread_manager;
   
   typedef ThreadManager::Mem_Fun_Data1
-    <DataOut<dim>, Data > MemFunData;
+    <DataOut<dim>, Data> MemFunData;
   
                                   // One struct of this type for every thread
   vector<MemFunData> mem_fun_data (n_threads,
index c8d794ab9ad21d8683c8d55cdbe421ac258d3db5..e1aebd488a8de8479265859b92e372659ea81e1c 100644 (file)
@@ -340,10 +340,9 @@ DoFHandler<dim>::active_cell_iterator cell=data.dof.begin_active();
                                   // be used, the threads would take widely
                                   // spread times to calculate their cells.
   for (unsigned int t=0;t<this_thread;++t,++cell);
-
                                   // loop over all cells for this thread
                                   // the iteration of cell is done at the end
-  for (;cell!=data.endc;)
+  for (; cell!=data.endc; )
     {
       
                                       // loop over all faces of this cell
@@ -381,7 +380,7 @@ DoFHandler<dim>::active_cell_iterator cell=data.dof.begin_active();
            };
 
 
-if (cell->face(face_no)->has_children() == false)
+         if (cell->face(face_no)->has_children() == false)
                                             // if the face is a regular one, i.e.
                                             // either on the other side there is
                                             // nirvana (face is at boundary), or
@@ -407,11 +406,13 @@ if (cell->face(face_no)->has_children() == false)
        };
 
                                       // next cell in this thread
-      for (unsigned int t=0;((t<data.n_threads)&&(cell!=data.endc));++t,++cell) {};
+      for (unsigned int t=0;((t<data.n_threads)&&(cell!=data.endc));++t,++cell)
+       {};
     };
 };
 
 
+
 template <int dim>
 void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>   &dof,
                                         const Quadrature<dim-1> &quadrature,
@@ -450,7 +451,7 @@ void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>   &dof,
                                   // of the cell.
                                   // the values for all faces are set to
                                   // -10e20. It would cost a lot of time
-                                  // to syncronisise the initialisation
+                                  // to synchronise the initialisation
                                   // of the map in multithreaded mode.
                                   // negative value indicates that the
                                   // face is not calculated.
@@ -460,7 +461,7 @@ void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>   &dof,
       data.face_integrals[cell->face(face_no)]=-10e20;
 
 
-// split all cells into threads
+                                  // split all cells into threads
                                   // if multithreading is used
 #ifdef DEAL_II_USE_MT
 
@@ -478,17 +479,15 @@ void KellyErrorEstimator<dim>::estimate (const DoFHandler<dim>   &dof,
              FunData (data,0,&KellyErrorEstimator::estimate_some));
 
 
-// get start cells for each thread
+                                  // get start cells for each thread
   for (unsigned int l=0;l<data.n_threads;++l)
-    {
-      fun_data[l].arg2=l;
-    };
+    fun_data[l].arg2=l;
+    
     
                                   // now spawn the threads
   for (unsigned int i=0;i<data.n_threads; ++i)
-    {
-      thread_manager.spawn(&fun_data[i],THR_SCOPE_SYSTEM | THR_DETACHED);
-    };
+    thread_manager.spawn(&fun_data[i],THR_SCOPE_SYSTEM | THR_DETACHED);
+    
                                   // wait for all threads to return
   thread_manager.wait();
   
@@ -665,7 +664,7 @@ integrate_over_regular_face (Data                       &data,
     };
 
 
-if (face->at_boundary() == true)
+  if (face->at_boundary() == true)
                                     // neumann boundary face. compute
                                     // difference between normal
                                     // derivative and boundary function
@@ -689,7 +688,7 @@ if (face->at_boundary() == true)
     };
 
 
-// now phi contains the following:
+                                  // now phi contains the following:
                                   // - for an internal face, phi=[a du/dn]
                                   // - for a neumann boundary face,
                                   //   phi=a du/dn-g
@@ -712,6 +711,7 @@ if (face->at_boundary() == true)
 };
 
 
+
 template <int dim>
 void KellyErrorEstimator<dim>::
 integrate_over_irregular_face (Data                       &data,
@@ -800,7 +800,7 @@ integrate_over_irregular_face (Data                       &data,
       data.normal_vectors[this_thread]=fe_face_values.get_normal_vectors();
 
 
-for (unsigned int component=0; component<data.n_components; ++component)
+      for (unsigned int component=0; component<data.n_components; ++component)
        for (unsigned int point=0; point<data.n_q_points; ++point)
          data.phi[this_thread][point][component] =
            data.psi[this_thread][point][component]*
@@ -847,7 +847,7 @@ for (unsigned int component=0; component<data.n_components; ++component)
     };
 
 
-// finally loop over all subfaces to
+                                  // finally loop over all subfaces to
                                   // collect the contributions of the
                                   // subfaces and store them with the
                                   // mother face
@@ -865,7 +865,6 @@ for (unsigned int component=0; component<data.n_components; ++component)
     };
 
   data.face_integrals[face] = sum;
-
 };
 
 

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.