diff --git a/fem/fe/fe_base.cpp b/fem/fe/fe_base.cpp index 645968b708..3afe7f2f65 100644 --- a/fem/fe/fe_base.cpp +++ b/fem/fe/fe_base.cpp @@ -665,56 +665,56 @@ const { // If the new Dof2Quad is already present, e.g. added in a previous call // or added by another omp thread, return. - if (DofToQuad::SearchArray(dof2quad_array, ir, - DofToQuad::LEXICOGRAPHIC_FULL)) - { return; } - - // Undo the native ordering which is what FiniteElement::GetDofToQuad - // returns. - auto *d2q_new = new DofToQuad(d2q); - d2q_new->mode = DofToQuad::LEXICOGRAPHIC_FULL; - const int nqpt = ir.GetNPoints(); - - const int b_dim = (range_type == VECTOR) ? dim : 1; - - for (int i = 0; i < nqpt; i++) + if (DofToQuad::SearchArray(dof2quad_array, ir, DofToQuad::LEXICOGRAPHIC_FULL) == nullptr) { - for (int d = 0; d < b_dim; d++) + // Undo the native ordering which is what FiniteElement::GetDofToQuad + // returns. + auto *d2q_new = new DofToQuad(d2q); + d2q_new->mode = DofToQuad::LEXICOGRAPHIC_FULL; + const int nqpt = ir.GetNPoints(); + + const int b_dim = (range_type == VECTOR) ? dim : 1; + + for (int i = 0; i < nqpt; i++) { - for (int j = 0; j < dof; j++) + for (int d = 0; d < b_dim; d++) { - const double val = d2q.B[i + nqpt*(d+b_dim*lex_ordering[j])]; - d2q_new->B[i+nqpt*(d+b_dim*j)] = val; - d2q_new->Bt[j+dof*(i+nqpt*d)] = val; + for (int j = 0; j < dof; j++) + { + const double val = d2q.B[i + nqpt*(d+b_dim*lex_ordering[j])]; + d2q_new->B[i+nqpt*(d+b_dim*j)] = val; + d2q_new->Bt[j+dof*(i+nqpt*d)] = val; + } } } - } - const int g_dim = [this]() - { - switch (deriv_type) + const int g_dim = [this]() { - case GRAD: return dim; - case DIV: return 1; - case CURL: return cdim; - default: return 0; - } - }(); - - for (int i = 0; i < nqpt; i++) - { - for (int d = 0; d < g_dim; d++) - { - for (int j = 0; j < dof; j++) + switch (deriv_type) { - const double val = d2q.G[i + nqpt*(d+g_dim*lex_ordering[j])]; - d2q_new->G[i+nqpt*(d+g_dim*j)] = val; - d2q_new->Gt[j+dof*(i+nqpt*d)] = val; + case GRAD: return dim; + case DIV: return 1; + case CURL: return cdim; + default: return 0; + } + }(); + + for (int i = 0; i < nqpt; i++) + { + for (int d = 0; d < g_dim; d++) + { + for (int j = 0; j < dof; j++) + { + const double val = d2q.G[i + nqpt*(d+g_dim*lex_ordering[j])]; + d2q_new->G[i+nqpt*(d+g_dim*j)] = val; + d2q_new->Gt[j+dof*(i+nqpt*d)] = val; + } } } + + dof2quad_array.Append(d2q_new); } - dof2quad_array.Append(d2q_new); } }