From: Wolfgang Bangerth Date: Mon, 10 Apr 2017 03:58:36 +0000 (-0600) Subject: Use lambdas and tasks to parallelize construction of FESystem. X-Git-Tag: v9.0.0-rc1~1702^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F4241%2Fhead;p=dealii.git Use lambdas and tasks to parallelize construction of FESystem. --- diff --git a/include/deal.II/fe/fe_system.h b/include/deal.II/fe/fe_system.h index 8af24fe5b7..73029e2c91 100644 --- a/include/deal.II/fe/fe_system.h +++ b/include/deal.II/fe/fe_system.h @@ -972,23 +972,6 @@ private: unsigned int> > base_elements; - /** - * Initialize the @p unit_support_points field of the FiniteElement class. - * Called from the constructor. - */ - void initialize_unit_support_points (); - - /** - * Initialize the @p unit_face_support_points field of the FiniteElement - * class. Called from the constructor. - */ - void initialize_unit_face_support_points (); - - /** - * Initialize the @p adjust_quad_dof_index_for_face_orientation_table field - * of the FiniteElement class. Called from the constructor. - */ - void initialize_quad_dof_index_permutation (); /** * This function is simply singled out of the constructors since there are * several of them. It sets up the index table for the system as well as @p diff --git a/source/fe/fe_system.cc b/source/fe/fe_system.cc index adcc12048a..a24c48b84a 100644 --- a/source/fe/fe_system.cc +++ b/source/fe/fe_system.cc @@ -1507,153 +1507,134 @@ void FESystem::initialize (const std::vector init_tasks; + init_tasks += Threads::new_task ([&]() + { + this->build_interface_constraints (); + }); + init_tasks += Threads::new_task ([&]() + { + // if one of the base elements has no support points, then it makes no sense + // to define support points for the composed element, so return an empty + // array to demonstrate that fact. Note that we ignore FE_Nothing in this logic. + for (unsigned int base_el=0; base_eln_base_elements(); ++base_el) + if (!base_element(base_el).has_support_points() && base_element(base_el).dofs_per_cell!=0) + { + this->unit_support_points.resize(0); + return; + } + // generate unit support points from unit support points of sub elements + this->unit_support_points.resize(this->dofs_per_cell); -template -void -FESystem:: -initialize_unit_support_points () -{ - // if one of the base elements has no support points, then it makes no sense - // to define support points for the composed element, so return an empty - // array to demonstrate that fact. Note that we ignore FE_Nothing in this logic. - for (unsigned int base_el=0; base_eln_base_elements(); ++base_el) - if (!base_element(base_el).has_support_points() && base_element(base_el).dofs_per_cell!=0) + for (unsigned int i=0; idofs_per_cell; ++i) { - this->unit_support_points.resize(0); - return; + const unsigned int + base = this->system_to_base_table[i].first.first, + base_index = this->system_to_base_table[i].second; + Assert (basen_base_elements(), ExcInternalError()); + Assert (base_indexunit_support_points[i] = base_element(base).unit_support_points[base_index]; } + }); - // generate unit support points from unit support points of sub elements - this->unit_support_points.resize(this->dofs_per_cell); - - for (unsigned int i=0; idofs_per_cell; ++i) + // initialize face support points (for dim==2,3). same procedure as above + if (dim > 1) + init_tasks += Threads::new_task ([&]() { - const unsigned int - base = this->system_to_base_table[i].first.first, - base_index = this->system_to_base_table[i].second; - Assert (basen_base_elements(), ExcInternalError()); - Assert (base_indexunit_support_points[i] = base_element(base).unit_support_points[base_index]; - } -} - - - -template -void -FESystem:: -initialize_unit_face_support_points () -{ - // Nothing to do in 1D - if (dim == 1) - return; - - // if one of the base elements has no support points, then it makes no sense - // to define support points for the composed element. In that case, return - // an empty array to demonstrate that fact (note that we ask whether the - // base element has no support points at all, not only none on the face!) - // - // on the other hand, if there is an element that simply has no degrees of - // freedom on the face at all, then we don't care whether it has support - // points or not. this is, for example, the case for the stable Stokes - // element Q(p)^dim \times DGP(p-1). - for (unsigned int base_el=0; base_eln_base_elements(); ++base_el) - if (!base_element(base_el).has_support_points() - && - (base_element(base_el).dofs_per_face > 0)) - { - this->unit_face_support_points.resize(0); - return; - } - + // if one of the base elements has no support points, then it makes no sense + // to define support points for the composed element. In that case, return + // an empty array to demonstrate that fact (note that we ask whether the + // base element has no support points at all, not only none on the face!) + // + // on the other hand, if there is an element that simply has no degrees of + // freedom on the face at all, then we don't care whether it has support + // points or not. this is, for example, the case for the stable Stokes + // element Q(p)^dim \times DGP(p-1). + for (unsigned int base_el=0; base_eln_base_elements(); ++base_el) + if (!base_element(base_el).has_support_points() + && + (base_element(base_el).dofs_per_face > 0)) + { + this->unit_face_support_points.resize(0); + return; + } - // generate unit face support points from unit support points of sub - // elements - this->unit_face_support_points.resize(this->dofs_per_face); - for (unsigned int i=0; idofs_per_face; ++i) - { - const unsigned int base_i = this->face_system_to_base_table[i].first.first; - const unsigned int index_in_base = this->face_system_to_base_table[i].second; + // generate unit face support points from unit support points of sub + // elements + this->unit_face_support_points.resize(this->dofs_per_face); - Assert (index_in_base < base_element(base_i).unit_face_support_points.size(), - ExcInternalError()); - - this->unit_face_support_points[i] - = base_element(base_i).unit_face_support_points[index_in_base]; - } -} + for (unsigned int i=0; idofs_per_face; ++i) + { + const unsigned int base_i = this->face_system_to_base_table[i].first.first; + const unsigned int index_in_base = this->face_system_to_base_table[i].second; + Assert (index_in_base < base_element(base_i).unit_face_support_points.size(), + ExcInternalError()); + this->unit_face_support_points[i] + = base_element(base_i).unit_face_support_points[index_in_base]; + } + }); -template -void -FESystem::initialize_quad_dof_index_permutation () -{ - // nothing to do in other dimensions than 3 - if (dim < 3) - return; - - // the array into which we want to write should have the correct size - // already. - Assert (this->adjust_quad_dof_index_for_face_orientation_table.n_elements()== - 8*this->dofs_per_quad, ExcInternalError()); - - // to obtain the shifts for this composed element, copy the shift - // information of the base elements - unsigned int index = 0; - for (unsigned int b=0; bn_base_elements(); ++b) + // initialize quad dof index permutation in 3d and higher + if (dim >= 3) + init_tasks += Threads::new_task ([&]() { - const Table<2,int> &temp - = this->base_element(b).adjust_quad_dof_index_for_face_orientation_table; - for (unsigned int c=0; celement_multiplicity(b); ++c) + // the array into which we want to write should have the correct size + // already. + Assert (this->adjust_quad_dof_index_for_face_orientation_table.n_elements()== + 8*this->dofs_per_quad, ExcInternalError()); + + // to obtain the shifts for this composed element, copy the shift + // information of the base elements + unsigned int index = 0; + for (unsigned int b=0; bn_base_elements(); ++b) { - for (unsigned int i=0; iadjust_quad_dof_index_for_face_orientation_table(index+i,j)= - temp(i,j); - index += temp.size(0); + const Table<2,int> &temp + = this->base_element(b).adjust_quad_dof_index_for_face_orientation_table; + for (unsigned int c=0; celement_multiplicity(b); ++c) + { + for (unsigned int i=0; iadjust_quad_dof_index_for_face_orientation_table(index+i,j)= + temp(i,j); + index += temp.size(0); + } } - } - Assert (index == this->dofs_per_quad, - ExcInternalError()); + Assert (index == this->dofs_per_quad, + ExcInternalError()); - // additionally compose the permutation information for lines - Assert (this->adjust_line_dof_index_for_line_orientation_table.size()== - this->dofs_per_line, ExcInternalError()); - index = 0; - for (unsigned int b=0; bn_base_elements(); ++b) - { - const std::vector &temp2 - = this->base_element(b).adjust_line_dof_index_for_line_orientation_table; - for (unsigned int c=0; celement_multiplicity(b); ++c) + // additionally compose the permutation information for lines + Assert (this->adjust_line_dof_index_for_line_orientation_table.size()== + this->dofs_per_line, ExcInternalError()); + index = 0; + for (unsigned int b=0; bn_base_elements(); ++b) { - std::copy(temp2.begin(), temp2.end(), - this->adjust_line_dof_index_for_line_orientation_table.begin() - +index); - index += temp2.size(); + const std::vector &temp2 + = this->base_element(b).adjust_line_dof_index_for_line_orientation_table; + for (unsigned int c=0; celement_multiplicity(b); ++c) + { + std::copy(temp2.begin(), temp2.end(), + this->adjust_line_dof_index_for_line_orientation_table.begin() + +index); + index += temp2.size(); + } } - } - Assert (index == this->dofs_per_line, - ExcInternalError()); + Assert (index == this->dofs_per_line, + ExcInternalError()); + }); + + // wait for all of this to finish + init_tasks.join_all(); }