]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Fix the same bug in PolynomialsABF.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 28 Jan 2010 17:04:05 +0000 (17:04 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 28 Jan 2010 17:04:05 +0000 (17:04 +0000)
git-svn-id: https://svn.dealii.org/trunk@20468 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/polynomials_abf.h
deal.II/base/source/polynomials_abf.cc
deal.II/doc/news/changes.h

index 0abd61ca1f79f704298da4879dd7fb7119f48191..583bf75ca8dbc917f23fffce6b8f1e84793ae85d 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2004, 2005, 2006, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -22,6 +22,7 @@
 #include <base/polynomial_space.h>
 #include <base/tensor_product_polynomials.h>
 #include <base/table.h>
+#include <base/thread_management.h>
 
 #include <vector>
 
@@ -70,7 +71,7 @@ class PolynomialsABF
                                      * Destructor deleting the polynomials.
                                      */
     ~PolynomialsABF ();
-    
+
                                     /**
                                      * Computes the value and the
                                      * first and second derivatives
@@ -98,20 +99,20 @@ class PolynomialsABF
     void compute (const Point<dim>            &unit_point,
                   std::vector<Tensor<1,dim> > &values,
                   std::vector<Tensor<2,dim> > &grads,
-                  std::vector<Tensor<3,dim> > &grad_grads) const;    
-    
+                  std::vector<Tensor<3,dim> > &grad_grads) const;
+
                                     /**
                                      * Returns the number of ABF polynomials.
                                      */
     unsigned int n () const;
-    
+
                                     /**
                                      * Returns the degree of the ABF
                                      * space, which is two less than
                                      * the highest polynomial degree.
                                      */
     unsigned int degree () const;
-    
+
                                     /**
                                      * Return the number of
                                      * polynomials in the space
@@ -122,14 +123,14 @@ class PolynomialsABF
                                      * classes.
                                      */
     static unsigned int compute_n_pols(unsigned int degree);
-    
+
   private:
                                     /**
                                      * The degree of this object as
                                      * given to the constructor.
                                      */
     const unsigned int my_degree;
-    
+
                                     /**
                                      * An object representing the
                                      * polynomial space for a single
@@ -144,17 +145,23 @@ class PolynomialsABF
                                      * polynomials.
                                      */
     unsigned int n_pols;
-    
+
+                                    /**
+                                     * A mutex that guards the
+                                     * following scratch arrays.
+                                     */
+    mutable Threads::Mutex mutex;
+
                                     /**
                                      * Auxiliary memory.
                                      */
-    mutable std::vector<double> p_values;    
-    
+    mutable std::vector<double> p_values;
+
                                     /**
                                      * Auxiliary memory.
                                      */
     mutable std::vector<Tensor<1,dim> > p_grads;
-    
+
                                     /**
                                      * Auxiliary memory.
                                      */
index 5bc0a39f5813c4e4c183bf59794ed0e528f9e1e3..b67f325a1d996da7367f1f5f3467182c9706af60 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2004, 2005, 2006, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -60,10 +60,17 @@ PolynomialsABF<dim>::compute (const Point<dim>            &unit_point,
         ExcDimensionMismatch(grad_grads.size(), n_pols));
 
   const unsigned int n_sub = polynomial_space->n();
+                                  // guard access to the scratch
+                                  // arrays in the following block
+                                  // using a mutex to make sure they
+                                  // are not used by multiple threads
+                                  // at once
+  Threads::Mutex::ScopedLock lock(mutex);
+
   p_values.resize((values.size() == 0) ? 0 : n_sub);
   p_grads.resize((grads.size() == 0) ? 0 : n_sub);
   p_grad_grads.resize((grad_grads.size() == 0) ? 0 : n_sub);
-  
+
   for (unsigned int d=0;d<dim;++d)
     {
                                       // First we copy the point. The
@@ -80,16 +87,16 @@ PolynomialsABF<dim>::compute (const Point<dim>            &unit_point,
       Point<dim> p;
       for (unsigned int c=0;c<dim;++c)
        p(c) = unit_point((c+d)%dim);
-      
+
       polynomial_space->compute (p, p_values, p_grads, p_grad_grads);
-      
+
       for (unsigned int i=0;i<p_values.size();++i)
-         values[i+d*n_sub][d] = p_values[i];
-      
+       values[i+d*n_sub][d] = p_values[i];
+
       for (unsigned int i=0;i<p_grads.size();++i)
        for (unsigned int d1=0;d1<dim;++d1)
          grads[i+d*n_sub][d][(d1+d)%dim] = p_grads[i][d1];
-      
+
       for (unsigned int i=0;i<p_grad_grads.size();++i)
        for (unsigned int d1=0;d1<dim;++d1)
          for (unsigned int d2=0;d2<dim;++d2)
@@ -107,7 +114,7 @@ PolynomialsABF<dim>::compute_n_pols(unsigned int k)
   if (dim == 2) return 2*(k+1)*(k+3);
   //TODO:Check what are the correct numbers ...
   if (dim == 3) return 3*(k+1)*(k+1)*(k+2);
-  
+
   Assert(false, ExcNotImplemented());
   return 0;
 }
index f3842f5a72edaf4c65783707cb57380f65bb43f8..48c2580f5e3362cb66169ce1e7671fedccb88236 100644 (file)
@@ -246,8 +246,8 @@ inconvenience this causes.
 <h3>base</h3>
 
 <ol>
-  <li><p>Fixed: The PolynomialsBDM had a race condition in multithreaded
-  programs. This now fixed.
+  <li><p>Fixed: The PolynomialsBDM and PolynomialsABF classes had a race
+  condition in multithreaded programs. This now fixed.
   <br>
   (WB 2010/01/27)
   </p></li>

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.