// $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
#include <base/polynomial_space.h>
#include <base/tensor_product_polynomials.h>
#include <base/table.h>
+#include <base/thread_management.h>
#include <vector>
* Destructor deleting the polynomials.
*/
~PolynomialsABF ();
-
+
/**
* Computes the value and the
* first and second derivatives
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
* 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
* 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.
*/
// $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
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
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)
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;
}
<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>