phynexis-deps/mfem
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

added lexicographic full for the rest of the possible element and derivative configurations

Browse Source
This commit is contained in:
Victor DeCaria
2023-10-10 06:06:04 -06:00
parent 44b34f5e64
commit 21b0b8bcf6
2 changed files with 79 additions and 4 deletions
Download Patch File Download Diff File Expand all files Collapse all files
fem/fe/fe_base.cpp
+76 -4
View File
@@ -2535,11 +2535,16 @@ const DofToQuad &NodalTensorFiniteElement::GetDofToQuad(
DofToQuad::Mode mode,
const ElementDofOrdering ordering) const
{
MFEM_VERIFY(!(mode == DofToQuad::Mode::TENSOR &&
ordering == ElementDofOrdering::NATIVE),
"Invalide combination of DofToQuad::Mode and ElementDofOrdering.");
for (int i = 0; i < dof2quad_array.Size(); i++)
{
const DofToQuad &d2q = *dof2quad_array[i];
if (d2q.IntRule == &ir && d2q.mode == mode && d2q.ordering == ordering) { return d2q; }
}
//First create the DofToQuad map for ElementDofOrdering::NATIVE.
//Either return d2q, or create the FULL and LEXICOGRAPHIC map and return.
auto &d2q = GetDofToQuad(ir, mode);
if (mode == DofToQuad::Mode::FULL &&
ordering == ElementDofOrdering::LEXICOGRAPHIC)
@@ -2548,17 +2553,84 @@ const DofToQuad &NodalTensorFiniteElement::GetDofToQuad(
auto *d2q_new = new DofToQuad(d2q);
d2q_new->ordering = ElementDofOrdering::LEXICOGRAPHIC;
const int nqpt = ir.GetNPoints();
for (int i = 0; i < nqpt; i++)
if (range_type == SCALAR)
{
for (int d = 0; d < dim; d++)
for (int i = 0; i < nqpt; i++)
{
for (int j = 0; j < dof; j++)
{
d2q_new->G[i+nqpt*(d+dim*j)] = d2q_new->Gt[j+dof*(i+nqpt*d)] = d2q.G[i+nqpt*
(d+dim*dof_map[j])];
d2q_new->B[i+nqpt*j] = d2q_new->Bt[j+dof*i] = d2q.B[i+nqpt*dof_map[j]];
}
}
}
else if (range_type == VECTOR)
{
for (int i = 0; i < nqpt; i++)
{
for (int d = 0; d < dim; d++)
{
for (int j = 0; j < dof; j++)
{
d2q_new->B[i+nqpt*(d+dim*j)] = d2q_new->Bt[j+dof*(i+nqpt*d)] = d2q.B[i+nqpt*
(d+dim*dof_map[j])];
}
}
}
}
else
{
// Skip B and Bt for unknown range type
}
switch (deriv_type)
{
case GRAD:
{
for (int i = 0; i < nqpt; i++)
{
for (int d = 0; d < dim; d++)
{
for (int j = 0; j < dof; j++)
{
d2q_new->G[i+nqpt*(d+dim*j)] = d2q_new->Gt[j+dof*(i+nqpt*d)] = d2q.G[i+nqpt*
(d+dim*dof_map[j])];
}
}
}
break;
}
case DIV:
{
for (int i = 0; i < nqpt; i++)
{
for (int j = 0; j < dof; j++)
{
d2q_new->G[i+nqpt*j] = d2q_new->Gt[j+dof*i] = d2q.G[i+nqpt*dof_map[j]];
}
}
break;
}
case CURL:
{
for (int i = 0; i < nqpt; i++)
{
for (int d = 0; d < cdim; d++)
{
for (int j = 0; j < dof; j++)
{
d2q_new->G[i+nqpt*(d+cdim*j)] = d2q_new->Gt[j+dof*(i+nqpt*d)] = d2q.G[i+nqpt*
(d+cdim*dof_map[j])];
}
}
}
break;
}
case NONE:
default:
// Skip G and Gt for unknown derivative type
break;
}
dof2quad_array.Append(d2q_new);
return *d2q_new;
}