LR multipatch: special case corners, lines and faces for interior propegation

This commit is contained in:
Kjetil Andre Johannessen 2017-10-12 14:26:22 +02:00 committed by Knut Morten Okstad
parent 755da42e2c
commit 66bd5f7035
3 changed files with 119 additions and 27 deletions

View File

@ -145,11 +145,21 @@ public:
//! \brief Returns all functions whose support overlap with the input nodes.
//! \param[in] nodes List of (0-based) patch local node IDs
//! (typically requested by adaptive refinement)
//! \param[in] dir In which parametric direction the basis functions are
//! allowed to have bigger support (-1 for all)
//! \param[in] dir 3-bit binary mask for which parameter directions are allowed
//! to grow; i.e. bin(011)=dec(3) allows u-direction and v-direction to grow,
//! default is bin(111)=dec(7) all directions
//! \return 0-based node IDs for functions with overlapping support with
//! the ones in boundary
IntVec getOverlappingNodes(const IntSet& nodes, int dir = -1) const;
IntVec getOverlappingNodes(const IntSet& nodes, int dir = 7) const;
//! \brief Returns all functions whose support overlap with the input node.
//! \param[in] node 0-based patch local node ID
//! \param[in] dir 3-bit binary mask for which parameter directions can grow
//! \return 0-based node IDs for functions with overlapping support
IntVec getOverlappingNodes(int node, int dir = 7) const
{
return this->getOverlappingNodes(IntSet(&node,(&node)+1),dir);
}
//! \brief Remaps element-wise errors from geometry mesh to refinement mesh.
//! \param errors The remapped errors

View File

@ -331,20 +331,18 @@ IntVec ASMunstruct::getOverlappingNodes (const IntSet& nodes, int dir) const
LR::Basisfunction *b = geo->getBasisfunction(i);
for (auto el : b->support()) // for all elements where *b has support
for (auto basis : el->support()) // for all functions on this element
if (dir == -1)
result.insert(basis->getId());
else
{
bool support_only_bigger_in_allowed_direction = true;
for (int j = 0; j < b->nVariate() && support_only_bigger_in_allowed_direction; j++)
{
bool support_only_bigger_in_dir_direction = true;
for (int j = 0; j < b->nVariate() && support_only_bigger_in_dir_direction; j++)
if (j != dir)
if (b->getParmin(j) > basis->getParmin(j) ||
b->getParmax(j) < basis->getParmax(j))
support_only_bigger_in_dir_direction = false;
if (support_only_bigger_in_dir_direction)
result.insert(basis->getId());
if ((1<<j) & dir) continue; // the function is allowed to grow in the direction j
if (b->getParmin(j) > basis->getParmin(j) ||
b->getParmax(j) < basis->getParmax(j))
support_only_bigger_in_allowed_direction = false;
}
if (support_only_bigger_in_allowed_direction)
result.insert(basis->getId());
}
}
return IntVec(result.begin(), result.end());

View File

@ -15,7 +15,10 @@
#include "SIMoptions.h"
#include "ModelGenerator.h"
#include "ASMstruct.h"
#include "ASMunstruct.h"
#ifdef HAS_LRSPLINE
#include "ASMu2D.h"
#include "ASMu3D.h"
#endif
#include "GlbL2projector.h"
#include "LinSolParams.h"
#include "Functions.h"
@ -1137,32 +1140,113 @@ bool SIMinput::refine (const LR::RefineData& prm,
}
}
#ifdef HAS_LRSPLINE
for (size_t i = 0; i < myModel.size(); i++)
{
pch = dynamic_cast<ASMunstruct*>(myModel[i]);
int nedg = (pch->getNoParamDim() == 2) ? 4 : 6;
std::vector<IntVec> bndry(nedg);
// OPTIMIZATION NOTE: if we by some clever datastructures already knew which edge
// each node in conformingIndices[i] was on, then we don't have to brute-force search
// for it like we do here. pch->getBoundaryNodes() above seem to compute this,
// but it is only for the sending patch boundary index, not the recieving patch boundary
// index.
for (int j = 1; j <= nedg; j++)
pch->getBoundaryNodes(j, bndry[j-1], 1, 1, 0, true);
for (int j : conformingIndices[i]) // add refinement from neighbours
{
bool done_with_this_node = false;
for (int edge = 0; edge < nedg && !done_with_this_node; edge++)
for (int edgeNode : bndry[edge])
if (myModel[i]->getNoParamDim() == 2) {
ASMu2D* lr = dynamic_cast<ASMu2D*>(myModel[i]);
int nedg0d = 4;
int nedg1d = 4;
IntVec bndry0;
std::vector<IntVec> bndry1(nedg1d);
for (int j = 0; j < nedg0d; j++)
bndry0.push_back(lr->getCorner((j%2)*2-1, (j/2)*2-1, 1));
for (int j = 1; j <= nedg1d; j++)
lr->getBoundaryNodes(j, bndry1[j-1], 1, 1, 0, true);
for (int j : conformingIndices[i]) // add refinement from neighbours
{
bool done_with_this_node = false;
// check if node is a corner node, compute large extended domain (all directions)
for (int edgeNode : bndry0)
if (edgeNode-1 == j)
{
IntVec secondary = pch->getOverlappingNodes(conformingIndices[i], edge/2);
IntVec secondary = lr->getOverlappingNodes(j);
refineIndices[i].insert(secondary.begin(),secondary.end());
done_with_this_node = true;
break;
}
// check if node is an edge node, compute small extended domain (one direction)
for (int edge1d = 0; edge1d < nedg1d && !done_with_this_node; edge1d++)
for (int edgeNode : bndry1[edge1d])
if (edgeNode-1 == j)
{
IntVec secondary = lr->getOverlappingNodes(j, edge1d/2+1);
refineIndices[i].insert(secondary.begin(),secondary.end());
done_with_this_node = true;
break;
}
}
}
else //volumetric models
{
ASMu3D* lr = dynamic_cast<ASMu3D*>(myModel[i]);
int nedg1d = 12;
int nedg2d = 6;
IntVec bndry0;
std::vector<IntVec> bndry1;
std::vector<IntVec> bndry2(nedg2d);
for (int K = -1; K < 2; K += 2)
for (int J = -1; J < 2; J += 2)
for (int I = -1; I < 2; I += 2)
bndry0.push_back(lr->getCorner(I,J,K,1));
for (int j = 1; j <= nedg1d; j++)
bndry1.push_back(lr->getEdge(j, true, 1, 0));
for (int j = 1; j <= nedg2d; j++)
lr->getBoundaryNodes(j, bndry2[j-1], 1, 1, 0, true);
for (int j : conformingIndices[i]) // add refinement from neighbours
{
bool done_with_this_node = false;
// check if node is a corner node, compute large extended domain (all directions)
for (int edgeNode : bndry0)
if (edgeNode-1 == j)
{
IntVec secondary = lr->getOverlappingNodes(j);
refineIndices[i].insert(secondary.begin(),secondary.end());
done_with_this_node = true;
break;
}
// check if node is an edge node, compute moderate extended domain (2 directions)
for (int edge1d = 0; edge1d < nedg1d && !done_with_this_node; edge1d++)
for (int edgeNode : bndry1[edge1d])
if (edgeNode-1 == j)
{
int allowed_direction;
if (edge1d < 4)
allowed_direction = 6; // bin(110), allowed to grow in v- and w-direction
else if (edge1d < 8)
allowed_direction = 5; // bin(101), allowed to grow in u- and w-direction
else
allowed_direction = 3; // bin(011), allowed to grow in u- and v-direction
IntVec secondary = lr->getOverlappingNodes(j, allowed_direction);
refineIndices[i].insert(secondary.begin(),secondary.end());
done_with_this_node = true;
break;
}
// check if node is a face node, compute small extended domain (1 direction)
for (int edge2d = 0; edge2d < nedg2d && !done_with_this_node; edge2d++)
for (int edgeNode : bndry2[edge2d])
if (edgeNode-1 == j)
{
IntVec secondary = lr->getOverlappingNodes(j, (1<<edge2d/2));
refineIndices[i].insert(secondary.begin(), secondary.end());
done_with_this_node = true;
break;
}
}
} // end volumetric
}
#else
return false;
#endif
Vectors lsols;
lsols.reserve(sol.size()*myModel.size());