From: bangerth Date: Thu, 28 Jan 2010 16:28:25 +0000 (+0000) Subject: Fix a race condition in PolynomialsBDM. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=abac0a8c8cf8aeccee54507981f15df21eb08f5f;p=dealii-svn.git Fix a race condition in PolynomialsBDM. git-svn-id: https://svn.dealii.org/trunk@20467 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/polynomials_bdm.h b/deal.II/base/include/base/polynomials_bdm.h index 89ed5a7124..71df841afb 100644 --- a/deal.II/base/include/base/polynomials_bdm.h +++ b/deal.II/base/include/base/polynomials_bdm.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2004, 2005, 2006, 2009 by the deal.II authors +// Copyright (C) 2004, 2005, 2006, 2009, 2010 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -21,6 +21,7 @@ #include #include #include +#include #include @@ -46,7 +47,7 @@ DEAL_II_NAMESPACE_OPEN *
In 2D: *
The 2D-curl of the functions xk+1y * and xyk+1. - *
In 3D: + *
In 3D: *
For any i=0,...,k the curls of * (0,0,xyi+1zk-i), * (xk-iyzi+1,0,0) and @@ -102,20 +103,20 @@ class PolynomialsBDM void compute (const Point &unit_point, std::vector > &values, std::vector > &grads, - std::vector > &grad_grads) const; - + std::vector > &grad_grads) const; + /** * Returns the number of BDM polynomials. */ unsigned int n () const; - + /** * Returns the degree of the BDM * space, which is one less than * the highest polynomial degree. */ unsigned int degree () const; - + /** * Return the number of * polynomials in the space @@ -126,7 +127,7 @@ class PolynomialsBDM * classes. */ static unsigned int compute_n_pols(unsigned int degree); - + private: /** * An object representing the @@ -144,23 +145,29 @@ class PolynomialsBDM * degree zero to k. */ std::vector > monomials; - + /** * Number of BDM * polynomials. */ unsigned int n_pols; - + + /** + * A mutex that guards the + * following scratch arrays. + */ + mutable Threads::Mutex mutex; + /** * Auxiliary memory. */ mutable std::vector p_values; - + /** * Auxiliary memory. */ mutable std::vector > p_grads; - + /** * Auxiliary memory. */ diff --git a/deal.II/base/source/polynomials_bdm.cc b/deal.II/base/source/polynomials_bdm.cc index 782348ee7d..c4ae72a235 100644 --- a/deal.II/base/source/polynomials_bdm.cc +++ b/deal.II/base/source/polynomials_bdm.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors +// Copyright (C) 2004, 2005, 2006, 2007, 2008, 2009, 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,40 +60,41 @@ PolynomialsBDM::compute (const Point &unit_point, ExcDimensionMismatch(grad_grads.size(), n_pols)); const unsigned int n_sub = polynomial_space.n(); - 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); - - // Compute values of complete space - // and insert into tensors. Result - // will have first all polynomials - // in the x-component, then y and - // z. - polynomial_space.compute (unit_point, p_values, p_grads, p_grad_grads); - - std::fill(values.begin(), values.end(), Tensor<1,dim>()); - for (unsigned int i=0;i()); + for (unsigned int i=0;i()); - for (unsigned int i=0;i()); + for (unsigned int i=0;i()); - for (unsigned int i=0;i()); + for (unsigned int i=0;i::compute (const Point &unit_point, // derivatives std::vector > monovali(dim, std::vector(4)); std::vector > monovalk(dim, std::vector(4)); - + if (dim == 2) { for (unsigned int d=0;d::compute_node_matrix (Table<2,double>& A) const double orientation = 1.; if ((face==0) || (face==3)) orientation = -1.; - + for (unsigned int k=0;k::compute_node_matrix (Table<2,double>& A) const p(0) = 1.; case 3: p(1) = x; - break; + break; } // std::cerr << p // << '\t' << moment_weight[0].value(x) // << '\t' << moment_weight[1].value(x) // ; - + compute (p, values, grads, grad_grads); - + for (unsigned int i=0;i::compute_node_matrix (Table<2,double>& A) const // std::cerr << std::endl; } } - + // Volume integrals are missing // // This degree is one larger diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index a2bfc6fb4e..f3842f5a72 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -246,6 +246,12 @@ inconvenience this causes.

base

    +
  1. Fixed: The PolynomialsBDM had a race condition in multithreaded + programs. This now fixed. +
    + (WB 2010/01/27) +

  2. +
  3. New: The new class IndexSet can represent sets and ranges of indices.
    (WB 2009/10/09)