]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Clone the base elements of FESystem in parallel. 1165/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 19 Jul 2015 14:15:09 +0000 (09:15 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 19 Jul 2015 14:15:09 +0000 (09:15 -0500)
In the constructor of FESystem, we clone the base elements, which requires
constructing new finite element objects -- an expensive task. This patch
does this in parallel for the base elements.

source/fe/fe_system.cc

index e460f80e0b8a8a926f2e4e0d2299f9453f570648..3cffc168ef03a2eeebc56d320e531f14ddb6cf11 100644 (file)
@@ -1849,6 +1849,13 @@ void FESystem<dim,spacedim>::initialize (const std::vector<const FiniteElement<d
     if (multiplicities[i]>0)
       this->base_to_block_indices.push_back( multiplicities[i] );
 
+  std::vector<Threads::Task<FiniteElement<dim,spacedim>*> > clone_base_elements;
+
+  for (unsigned int i=0; i<fes.size(); i++)
+    if (multiplicities[i]>0)
+      clone_base_elements.push_back (Threads::new_task (&FiniteElement<dim,spacedim>::clone,
+                                                        *fes[i]));
+
   unsigned int ind=0;
   for (unsigned int i=0; i<fes.size(); i++)
     {
@@ -1856,7 +1863,7 @@ void FESystem<dim,spacedim>::initialize (const std::vector<const FiniteElement<d
         {
           base_elements[ind] =
             std::make_pair (std_cxx11::shared_ptr<const FiniteElement<dim,spacedim> >
-                            (fes[i]->clone()),
+                            (clone_base_elements[ind].return_value()),
                             multiplicities[i]);
           ++ind;
         }

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.