From e42ba419aca7370bd0892074d69a336aa71c05e5 Mon Sep 17 00:00:00 2001 From: kronbichler Date: Sat, 25 May 2013 07:49:44 +0000 Subject: [PATCH] Use locks only when modifying entries, not whenever accessing through get_prolongation_matrix git-svn-id: https://svn.dealii.org/trunk@29582 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/source/fe/fe_q_base.cc | 32 ++++++++++++++----- deal.II/source/fe/fe_system.cc | 56 +++++++++++++++------------------- 2 files changed, 49 insertions(+), 39 deletions(-) diff --git a/deal.II/source/fe/fe_q_base.cc b/deal.II/source/fe/fe_q_base.cc index 0e106d7c2e..43d46d3dc9 100644 --- a/deal.II/source/fe/fe_q_base.cc +++ b/deal.II/source/fe/fe_q_base.cc @@ -1020,9 +1020,15 @@ FE_Q_Base ExcIndexRange(child,0,GeometryInfo::n_children(refinement_case))); // initialization upon first request - Threads::Mutex::ScopedLock lock(this->mutex); if (this->prolongation[refinement_case-1][child].n() == 0) { + Threads::Mutex::ScopedLock lock(this->mutex); + + // if matrix got updated while waiting for the lock + if (this->prolongation[refinement_case-1][child].n() == + this->dofs_per_cell) + return this->prolongation[refinement_case-1][child]; + // distinguish q/q_dg0 case: only treat Q dofs first const unsigned int q_dofs_per_cell = Utilities::fixed_power(this->degree+1); @@ -1072,9 +1078,7 @@ FE_Q_Base } } - FullMatrix& prolongate = const_cast&> - (this->prolongation[refinement_case-1][child]); - prolongate.reinit(this->dofs_per_cell, this->dofs_per_cell); + FullMatrix prolongate (this->dofs_per_cell, this->dofs_per_cell); // go through the points in diagonal to capture variation in all // directions simultaneously @@ -1167,6 +1171,10 @@ FE_Q_Base Assert (std::fabs(sum-1.) < eps, ExcInternalError()); } #endif + + // swap matrices + std::swap(const_cast &> + (this->prolongation[refinement_case-1][child]), prolongate); } // finally return the matrix @@ -1189,12 +1197,16 @@ FE_Q_Base ExcIndexRange(child,0,GeometryInfo::n_children(refinement_case))); // initialization upon first request - Threads::Mutex::ScopedLock lock(this->mutex); if (this->restriction[refinement_case-1][child].n() == 0) { - FullMatrix &restriction = - const_cast &>(this->restriction[refinement_case-1][child]); + Threads::Mutex::ScopedLock lock(this->mutex); + + // if matrix got updated while waiting for the lock... + if (this->restriction[refinement_case-1][child].n() == + this->dofs_per_cell) + return this->restriction[refinement_case-1][child]; + FullMatrix restriction(this->dofs_per_cell, this->dofs_per_cell); // distinguish q/q_dg0 case const unsigned int q_dofs_per_cell = Utilities::fixed_power(this->degree+1); @@ -1291,7 +1303,13 @@ FE_Q_Base restriction(this->dofs_per_cell-1,this->dofs_per_cell-1) = 1./GeometryInfo::n_children(RefinementCase(refinement_case)); } + + // swap matrices + std::swap(const_cast &> + (this->restriction[refinement_case-1][child]), restriction); + } + return this->restriction[refinement_case-1][child]; } diff --git a/deal.II/source/fe/fe_system.cc b/deal.II/source/fe/fe_system.cc index 76c79d2b8e..f20f55657b 100644 --- a/deal.II/source/fe/fe_system.cc +++ b/deal.II/source/fe/fe_system.cc @@ -607,6 +607,13 @@ FESystem Threads::Mutex::ScopedLock lock(this->mutex); if (this->restriction[refinement_case-1][child].n() == 0) { + Threads::Mutex::ScopedLock lock(this->mutex); + + // check if updated while waiting for lock + if (this->restriction[refinement_case-1][child].n() == + this->dofs_per_cell) + return this->restriction[refinement_case-1][child]; + // Check if some of the matrices of the base elements are void. bool do_restriction = true; @@ -627,10 +634,8 @@ FESystem // if we did not encounter void matrices, initialize the matrix sizes if (do_restriction) { - FullMatrix &restriction = - const_cast&>(this->restriction[refinement_case-1] - [child]); - restriction.reinit(this->dofs_per_cell, this->dofs_per_cell); + FullMatrix restriction(this->dofs_per_cell, + this->dofs_per_cell); // distribute the matrices of the base finite elements to the // matrices of this object. for this, loop over all degrees of @@ -664,6 +669,9 @@ FESystem // entries of the matrices: restriction(i,j) = (*base_matrices[base])(base_index_i,base_index_j); } + + std::swap(restriction, const_cast &> + (this->restriction[refinement_case-1][child])); } } @@ -685,17 +693,19 @@ FESystem Assert (child::n_children(refinement_case), ExcIndexRange(child,0,GeometryInfo::n_children(refinement_case))); - // initialization upon first request - Threads::Mutex::ScopedLock lock(this->mutex); + // initialization upon first request, construction completely analogous to + // restriction matrix if (this->prolongation[refinement_case-1][child].n() == 0) { - // Check if some of the matrices of the base elements are void. - bool do_prolongation = true; + Threads::Mutex::ScopedLock lock(this->mutex); + + if (this->prolongation[refinement_case-1][child].n() == + this->dofs_per_cell) + return this->prolongation[refinement_case-1][child]; - // shortcut for accessing local prolongations further down + bool do_prolongation = true; std::vector *> base_matrices(this->n_base_elements()); - for (unsigned int i=0; in_base_elements(); ++i) { base_matrices[i] = @@ -706,45 +716,27 @@ FESystem Assert(do_prolongation, (typename FiniteElement::ExcEmbeddingVoid())); - if (do_prolongation) { - FullMatrix &prolongate = - const_cast &>(this->prolongation[refinement_case-1][child]); - prolongate.reinit(this->dofs_per_cell, this->dofs_per_cell); + FullMatrix prolongate (this->dofs_per_cell, + this->dofs_per_cell); - // distribute the matrices of the base finite elements to the - // matrices of this object. for this, loop over all degrees of - // freedom and take the respective entry of the underlying base - // element. - // - // note that we by definition of a base element, they are - // independent, i.e. do not couple. only DoFs that belong to the - // same instance of a base element may couple for (unsigned int i=0; idofs_per_cell; ++i) for (unsigned int j=0; jdofs_per_cell; ++j) { - // first find out to which base element indices i and j - // belong, and which instance thereof in case the base element - // has a multiplicity greater than one. if they should not - // happen to belong to the same instance of a base element, - // then they cannot couple, so go on with the next index if (this->system_to_base_table[i].first != this->system_to_base_table[j].first) continue; - - // so get the common base element and the indices therein: const unsigned int base = this->system_to_base_table[i].first.first; const unsigned int base_index_i = this->system_to_base_table[i].second, base_index_j = this->system_to_base_table[j].second; - - // if we are sure that DoFs i and j may couple, then copy - // entries of the matrices: prolongate(i,j) = (*base_matrices[base])(base_index_i,base_index_j); } + std::swap(prolongate, const_cast &> + (this->prolongation[refinement_case-1][child])); } } -- 2.39.5