]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Fix by Patrick Sodre: FEFieldFunction was not thread-safe due to the cell_hint cache...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Sep 2011 23:21:05 +0000 (23:21 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Sep 2011 23:21:05 +0000 (23:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@24281 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/numerics/fe_field_function.h
deal.II/include/deal.II/numerics/fe_field_function.templates.h

index c23460535577073ff112e07c2edc7d25e078df58..6adb15686de02a8b583f40553a89bf43237b136a 100644 (file)
@@ -283,6 +283,13 @@ and DoF handlers embedded in higher dimensional space have been added.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> Fixed: FEFieldFunction class was not thread-safe because it keeps a
+cache on the side that was invalidated when different
+threads kept pouncing on it. This is now fixed.
+<br>
+(Patrick Sodr&eacute;, 2011/09/07)
+
+
 <li> New: There is now a function GridTools::volume computing the volume
 of a triangulation.
 <br>
index c57fca493c98eeee6b6ce137ab8f067ebc2db77e..151c5255e9163df3f8bc3a6547b5fe0a789e56bb 100644 (file)
@@ -19,6 +19,7 @@
 #include <deal.II/base/function.h>
 #include <deal.II/base/point.h>
 #include <deal.II/base/tensor.h>
+#include <deal.II/base/thread_local_storage.h>
 
 #include <deal.II/lac/vector.h>
 
@@ -150,10 +151,10 @@ namespace Functions
                                        * points lie, you can tell
                                        * this object by calling this
                                        * function. This will speed
-                                       * things up a little.
+                                       * things up a little. 
                                        */
       void set_active_cell(typename DH::active_cell_iterator &newcell);
-  
+
                                       /**
                                        * Get ONE vector value at the
                                        * given point. It is
@@ -329,6 +330,12 @@ namespace Functions
                              std::vector<std::vector<unsigned int> > &maps) const;
     
     private:
+                                      /**
+                                       * Typedef holding the local cell_hint.
+                                       */
+      typedef Threads::ThreadLocalStorage < 
+       typename DH::active_cell_iterator > cell_hint_t; 
+
                                       /**
                                        * Pointer to the dof handler.
                                        */
@@ -347,10 +354,9 @@ namespace Functions
       const Mapping<dim> & mapping;
   
                                       /**
-                                       * The current cell in which we
-                                       * are evaluating.
+                                       * The latest cell hint.
                                        */
-      mutable typename DH::active_cell_iterator cell;
+      mutable cell_hint_t cell_hint;
   
                                       /**
                                        * Store the number of
index 2a220a7e78680b7fd36d8e842b5627f08eb16e05..70dce9ed4417cf7efc9c41d5f56cbcf2303de667 100644 (file)
@@ -31,9 +31,9 @@ namespace Functions
                  dh(&mydh, "FEFieldFunction"),
                  data_vector(myv),
                  mapping(mymapping),
+                 cell_hint(dh->end()),
                  n_components(mydh.get_fe().n_components())
   {
-    cell = dh->begin_active();
   }
 
 
@@ -43,17 +43,20 @@ namespace Functions
   FEFieldFunction<dim, DH, VECTOR>::
   set_active_cell(typename DH::active_cell_iterator &newcell)
   {
-    cell = newcell;
+    cell_hint.get() = newcell;
   }
 
 
-  
+
   template <int dim, typename DH, typename VECTOR>
   void FEFieldFunction<dim, DH, VECTOR>::vector_value (const Point<dim> &p,
                                                       Vector<double>   &values) const 
   { 
     Assert (values.size() == n_components, 
            ExcDimensionMismatch(values.size(), n_components));
+    typename DH::active_cell_iterator cell = cell_hint.get();
+    if (cell == dh->end()) 
+       cell = dh->begin_active();
     Point<dim> qp = mapping.transform_real_to_unit_cell(cell, p);
   
                                     // Check if we already have all we need
@@ -64,6 +67,8 @@ namespace Functions
        cell = my_pair.first;
        qp = my_pair.second;
       }                                    
+
+    cell_hint.get() = cell;
   
                                     // Now we can find out about the point
     Quadrature<dim> quad(qp);
@@ -96,6 +101,9 @@ namespace Functions
   { 
     Assert (gradients.size() == n_components, 
            ExcDimensionMismatch(gradients.size(), n_components));
+    typename DH::active_cell_iterator cell = cell_hint.get();
+    if(cell == dh->end())
+       cell = dh->begin_active();
     Point<dim> qp = mapping.transform_real_to_unit_cell(cell, p);
   
                                     // Check if we already have all we need
@@ -106,7 +114,9 @@ namespace Functions
        cell = my_pair.first;
        qp = my_pair.second;
       }                                    
-  
+
+    cell_hint.get() = cell;
+   
                                     // Now we can find out about the point
     Quadrature<dim> quad(qp);
     FEValues<dim> fe_v(mapping, dh->get_fe(), quad, 
@@ -268,6 +278,9 @@ namespace Functions
     bool left_over = true;
   
                                     // Current quadrature point
+    typename DH::active_cell_iterator cell = cell_hint.get();
+    if (cell == dh->end())
+       cell = dh->begin_active();
     Point<dim> qp = mapping.transform_real_to_unit_cell(cell, points[0]);
   
                                     // Check if we already have a

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.