OilfieldToolsNet/IFEM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

extend: addConnection to thick overlaps

Browse Source 
this adds the possibility to establish connections
with non-trivial overlap.

in the process make the serial connectPatch and the parallel
DomainDecomposition share more code.
This commit is contained in:
timovanopstal
2016-12-09 10:43:53 +01:00
parent 38858395b0
commit 9930e34b29
22 changed files with 273 additions and 320 deletions
Download Patch File Download Diff File Expand all files Collapse all files

5
src/ASM/ASMbase.h
View File

@@ -210,8 +210,11 @@ public:
//! \param[in] lIndex Local index of the boundary face/edge
//! \param nodes Array of node numbers
//! \param basis Which basis to grab nodes for (0 for all)
//! \param thick Thickness of connection
//! \param local If true, return patch-local numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes,
int basis = 0, bool local = false) const = 0;
int basis = 0, int thick = 1,
bool local = false) const = 0;
//! \brief Finds the node that is closest to the given point.
virtual std::pair<size_t,double> findClosestNode(const Vec3&) const

96
src/ASM/ASMs1D.C
View File

@@ -309,9 +309,11 @@ bool ASMs1D::generateTwistedFEModel (const RealFunc& twist, const Vec3& Zaxis)
}
bool ASMs1D::connectPatch (int vertex, ASMs1D& neighbor, int nvertex)
bool ASMs1D::connectPatch (int vertex, ASMs1D& neighbor, int nvertex, int thick)
{
if (!this->connectBasis(vertex,neighbor,nvertex))
int slave = vertex == 1 ? 0 : this->getSize(1)-thick;
int master = nvertex == 1 ? 0 : neighbor.getSize(1)-thick;
if (!this->connectBasis(vertex,neighbor,nvertex,1,slave,master,thick))
return false;
this->addNeighbor(&neighbor);
@@ -320,64 +322,40 @@ bool ASMs1D::connectPatch (int vertex, ASMs1D& neighbor, int nvertex)
bool ASMs1D::connectBasis (int vertex, ASMs1D& neighbor, int nvertex,
int basis, int slave, int master)
int basis, int slave, int master, int thick)
{
if (shareFE && neighbor.shareFE)
return true;
else if (shareFE || neighbor.shareFE)
{
std::cerr <<" *** ASMs1D::connectPatch: Logic error, cannot"
std::cerr <<" *** ASMs1D::connectBasis: Logic error, cannot"
<<" connect a sharedFE patch with an unshared one"<< std::endl;
return false;
}
// Set up the slave node number for this curve patch
int n1 = this->getSize(basis);
switch (vertex)
{
case 2: // Positive I-direction
slave += n1;
case 1: // Negative I-direction
slave += 1;
break;
default:
std::cerr <<" *** ASMs1D::connectPatch: Invalid slave vertex "
<< vertex << std::endl;
return false;
}
// Set up the master node number for the neighboring patch
n1 = neighbor.getSize(basis);
switch (nvertex)
{
case 2: // Positive I-direction
master += n1;
case 1: // Negative I-direction
master += 1;
break;
default:
std::cerr <<" *** ASMs1D::connectPatch: Invalid master vertex "
<< nvertex << std::endl;
return false;
}
const double xtol = 1.0e-4;
if (!neighbor.getCoord(master).equal(this->getCoord(slave),xtol))
{
std::cerr <<" *** ASMs1D::connectPatch: Non-matching nodes "
<< master <<": "<< neighbor.getCoord(master)
<<"\n and "
<< slave <<": "<< this->getCoord(slave) << std::endl;
} else if (vertex < 1 || vertex > 2) {
std::cerr <<" *** ASMs1D::connectBasis: Invalid slave vertex "
<< vertex << std::endl;
return false;
} else if (nvertex < 1 || nvertex > 2) {
std::cerr <<" *** ASMs1D::connectBasis: Invalid master vertex "
<< nvertex << std::endl;
return false;
}
else
ASMbase::collapseNodes(neighbor,master,*this,slave);
const double xtol = 1.0e-4;
for (int i = 0; i < thick; ++i) {
slave += 1; // add first, getCoord is 1-based
master += 1;
if (!neighbor.getCoord(master).equal(this->getCoord(slave),xtol))
{
std::cerr <<" *** ASMs1D::connectBasis: Non-matching nodes "
<< master <<": "<< neighbor.getCoord(master)
<<"\n and "
<< slave <<": "<< this->getCoord(slave) << std::endl;
return false;
}
else
ASMbase::collapseNodes(neighbor,master,*this,slave);
}
return true;
}
@@ -585,19 +563,21 @@ bool ASMs1D::updateRotations (const Vector& displ, bool reInit)
void ASMs1D::getBoundaryNodes (int lIndex, IntVec& nodes,
int, bool local) const
int, int thick, bool local) const
{
if (!curv) return; // silently ignore empty patches
size_t iel = lIndex == 1 ? 0 : nel-1;
if (MLGE[iel] > 0)
{
int node;
if (lIndex == 1)
node = MNPC[iel].front();
else if (lIndex == 2)
node = MNPC[iel][curv->order()-1];
nodes.push_back(local ? node+1 : MLGN[node]);
for (int i = 0; i < thick; ++i) {
int node;
if (lIndex == 1)
node = MNPC[iel][i];
else if (lIndex == 2)
node = MNPC[iel][curv->order()-thick+i];
nodes.push_back(local ? node+1 : MLGN[node]);
}
}
}

10
src/ASM/ASMs1D.h
View File

@@ -108,9 +108,10 @@ public:
//! \brief Finds the global (or patch-local) node number on a patch end.
//! \param[in] lIndex Local index of the end point
//! \param nodes Array of global boundary node numbers
//! \param thick Thickness of connection
//! \param local If \e true, return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes,
int, bool local) const;
int, int thick, bool local) const;
//! \brief Finds the node that is closest to the given point.
//! \param[in] X Global coordinates of point to search for
@@ -149,7 +150,9 @@ public:
//! \param[in] vertex Local vertex index of this patch, in range [1,2]
//! \param neighbor The neighbor patch
//! \param[in] nvertex Local vertex index of neighbor patch, in range [1,2]
virtual bool connectPatch(int vertex, ASMs1D& neighbor, int nvertex);
//! \param[in] thick Thickness of connection
virtual bool connectPatch(int vertex, ASMs1D& neighbor, int nvertex,
int thick = 1);
//! \brief Makes the two end vertices of the curve periodic.
//! \param[in] basis Which basis to connect (mixed methods), 0 means both
@@ -296,8 +299,9 @@ protected:
//! \param[in] basis Which basis to connect the nodes for (mixed methods)
//! \param[in] slave 0-based index of the first slave node in this basis
//! \param[in] master 0-based index of the first master node in this basis
//! \param[in] thick Thickness of connection
bool connectBasis(int vertex, ASMs1D& neighbor, int nvertex,
int basis = 1, int slave = 0, int master = 0);
int basis = 1, int slave = 0, int master = 0, int thick = 1);
//! \brief Extracts parameter values of the Gauss points.
//! \param[out] uGP Parameter values for all points

105
src/ASM/ASMs2D.C
View File

@@ -591,7 +591,7 @@ bool ASMs2D::assignNodeNumbers (BlockNodes& nodes, int basis)
bool ASMs2D::connectPatch (int edge, ASMs2D& neighbor, int nedge,
bool revers, int, bool coordCheck)
bool revers, int, bool coordCheck, int thick)
{
if (swapV && edge > 2) // Account for swapped parameter direction
edge = 7-edge;
@@ -599,7 +599,7 @@ bool ASMs2D::connectPatch (int edge, ASMs2D& neighbor, int nedge,
if (neighbor.swapV && nedge > 2) // Account for swapped parameter direction
nedge = 7-nedge;
if (!this->connectBasis(edge,neighbor,nedge,revers,1,0,0,coordCheck))
if (!this->connectBasis(edge,neighbor,nedge,revers,1,0,0,coordCheck,thick))
return false;
this->addNeighbor(&neighbor);
@@ -608,96 +608,52 @@ bool ASMs2D::connectPatch (int edge, ASMs2D& neighbor, int nedge,
bool ASMs2D::connectBasis (int edge, ASMs2D& neighbor, int nedge, bool revers,
int basis, int slave, int master, bool coordCheck)
int basis, int slave, int master,
bool coordCheck, int thick)
{
if (this->isShared() && neighbor.isShared())
return true;
else if (this->isShared() || neighbor.isShared())
{
std::cerr <<" *** ASMs2D::connectPatch: Logic error, cannot"
std::cerr <<" *** ASMs2D::connectBasis: Logic error, cannot"
<<" connect a shared patch with an unshared one"<< std::endl;
return false;
}
// Set up the slave node numbers for this surface patch
int n1, n2;
if (!this->getSize(n1,n2,basis)) return false;
int node = slave+1, i1 = 0;
switch (edge)
{
case 2: // Positive I-direction
node += n1-1;
case 1: // Negative I-direction
i1 = n1;
n1 = n2;
break;
case 4: // Positive J-direction
node += n1*(n2-1);
case 3: // Negative J-direction
i1 = 1;
break;
default:
std::cerr <<" *** ASMs2D::connectPatch: Invalid slave edge "
<< edge << std::endl;
return false;
}
int i;
IntVec slaveNodes(n1,0);
for (i = 0; i < n1; i++, node += i1)
slaveNodes[i] = node;
IntVec slaveNodes;
this->getBoundaryNodes(edge, slaveNodes, basis, thick, true);
for (int& it : slaveNodes)
it += slave;
// Set up the master node numbers for the neighboring surface patch
IntVec masterNodes;
neighbor.getBoundaryNodes(nedge, masterNodes, basis, thick, true);
for (int& it : masterNodes)
it += master;
if (!neighbor.getSize(n1,n2,basis)) return false;
node = master+1; i1 = 0;
switch (nedge)
{
case 2: // Positive I-direction
node += n1-1;
case 1: // Negative I-direction
i1 = n1;
n1 = n2;
break;
case 4: // Positive J-direction
node += n1*(n2-1);
case 3: // Negative J-direction
i1 = 1;
break;
default:
std::cerr <<" *** ASMs2D::connectPatch: Invalid master edge "
<< nedge << std::endl;
return false;
}
if (n1 != (int)slaveNodes.size())
if (masterNodes.size() != slaveNodes.size())
{
std::cerr <<" *** ASMs2D::connectPatch: Non-matching edges, sizes "
<< n1 <<" and "<< slaveNodes.size() << std::endl;
std::cerr <<" *** ASMs2D::connectBasis: Non-matching edges, sizes "
<< masterNodes.size() <<" and "<< slaveNodes.size() << std::endl;
return false;
}
const double xtol = 1.0e-4;
for (i = 0; i < n1; i++, node += i1)
for (size_t i = 0; i < masterNodes.size(); ++i)
{
int slave = slaveNodes[revers ? n1-i-1 : i];
int node = masterNodes[i];
int slave = slaveNodes[revers ? slaveNodes.size()-i-1 : i];
if (!coordCheck)
ASMbase::collapseNodes(neighbor,node,*this,slave);
else if (neighbor.getCoord(node).equal(this->getCoord(slave),xtol))
ASMbase::collapseNodes(neighbor,node,*this,slave);
else
{
std::cerr <<" *** ASMs2D::connectPatch: Non-matching nodes "
<< node <<": "<< neighbor.getCoord(node)
<<"\n and "
<< slave <<": "<< this->getCoord(slave) << std::endl;
std::cerr <<" *** ASMs2D::connectBasis: Non-matching nodes "
<< node <<": "<< neighbor.getCoord(node)
<<"\n and "
<< slave <<": "<< this->getCoord(slave) << std::endl;
return false;
}
}
@@ -1292,7 +1248,8 @@ bool ASMs2D::updateCoords (const Vector& displ)
}
void ASMs2D::getBoundaryNodes (int lIndex, IntVec& nodes, int basis, bool local) const
void ASMs2D::getBoundaryNodes (int lIndex, IntVec& nodes, int basis,
int thick, bool local) const
{
if (basis == 0)
basis = 1;
@@ -1307,23 +1264,25 @@ void ASMs2D::getBoundaryNodes (int lIndex, IntVec& nodes, int basis, bool local)
for (char i = 1; i <= basis; i++)
if (!this->getSize(n1,n2,i))
return;
else if (i < basis)
else if (i < basis && !local)
node += n1*n2;
switch (lIndex)
{
case 2: // Right edge (positive I-direction)
node += n1-1;
node += n1-thick;
case 1: // Left edge (negative I-direction)
for (int i2 = 1; i2 <= n2; i2++, node += n1)
nodes.push_back(local ? node : this->getNodeID(node));
for (int t = 0; t < thick; ++t)
nodes.push_back(local ? node+t : this->getNodeID(node+t));
break;
case 4: // Back edge (positive J-direction)
node += n1*(n2-1);
node += n1*(n2-thick);
case 3: // Front edge (negative J-direction)
for (int i1 = 1; i1 <= n1; i1++, node++)
nodes.push_back(local ? node : this->getNodeID(node));
for (int t = 0; t < thick; ++t)
nodes.push_back(local ? node + t*n1 : this->getNodeID(node + t*n1));
break;
}

12
src/ASM/ASMs2D.h
View File

@@ -183,9 +183,11 @@ public:
//! \param[in] lIndex Local index of the boundary edge
//! \param nodes Array of node numbers
//! \param basis Which basis to grab nodes for (0 for all)
//! \param local If \e true, return patch-local node numbers
//! \param thick Thickness of connection
//! \param local If \e true return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes,
int basis, bool local) const;
int basis, int thick = 1,
bool local = false) const;
//! \brief Returns the node index for a given corner.
virtual int getCorner(int I, int J, int basis = 1) const;
@@ -273,8 +275,9 @@ public:
//! \param[in] nedge Local edge index of neighbor patch, in range [1,4]
//! \param[in] revers Indicates whether the two edges have opposite directions
//! \param[in] coordCheck False to disable coordinate checks (periodic connections)
//! \param[in] thick Thickness of connection
virtual bool connectPatch(int edge, ASMs2D& neighbor, int nedge, bool revers,
int = 0, bool coordCheck = true);
int = 0, bool coordCheck = true, int thick = 1);
//! \brief Makes two opposite boundary edges periodic.
//! \param[in] dir Parameter direction defining the periodic edges
@@ -480,9 +483,10 @@ protected:
//! \param[in] slave 0-based index of the first slave node in this basis
//! \param[in] master 0-based index of the first master node in this basis
//! \param[in] coordCheck False to disable coordinate checks (periodic connections)
//! \param[in] thick Thickness of connection
bool connectBasis(int edge, ASMs2D& neighbor, int nedge, bool revers,
int basis = 1, int slave = 0, int master = 0,
bool coordCheck = true);
bool coordCheck = true, int thick = 1);
//! \brief Extracts parameter values of the Gauss points in one direction.
//! \param[out] uGP Parameter values in given direction for all points

12
src/ASM/ASMs2Dmx.C
View File

@@ -293,7 +293,7 @@ bool ASMs2Dmx::generateFEMTopology ()
bool ASMs2Dmx::connectPatch (int edge, ASMs2D& neighbor, int nedge, bool revers,
int basis, bool coordCheck)
int basis, bool coordCheck, int thick)
{
ASMs2Dmx* neighMx = dynamic_cast<ASMs2Dmx*>(&neighbor);
if (!neighMx) return false;
@@ -301,7 +301,7 @@ bool ASMs2Dmx::connectPatch (int edge, ASMs2D& neighbor, int nedge, bool revers,
size_t nb1 = 0, nb2 = 0;
for (size_t i = 1; i <= m_basis.size(); ++i) {
if (basis == 0 || i == (size_t)basis)
if (!this->connectBasis(edge,neighbor,nedge,revers,i,nb1,nb2,coordCheck))
if (!this->connectBasis(edge,neighbor,nedge,revers,i,nb1,nb2,coordCheck,thick))
return false;
nb1 += nb[i-1];
@@ -1014,12 +1014,12 @@ void ASMs2Dmx::generateThreadGroups (const Integrand& integrand, bool silence,
}
void ASMs2Dmx::getBoundaryNodes (int lIndex, IntVec& nodes,
int basis, bool local) const
void ASMs2Dmx::getBoundaryNodes (int lIndex, IntVec& nodes, int basis,
int thick, bool local) const
{
if (basis > 0)
this->ASMs2D::getBoundaryNodes(lIndex, nodes, basis, local);
this->ASMs2D::getBoundaryNodes(lIndex, nodes, basis, thick, local);
else
for (size_t b = 1; b <= this->getNoBasis(); ++b)
this->ASMs2D::getBoundaryNodes(lIndex, nodes, b, local);
this->ASMs2D::getBoundaryNodes(lIndex, nodes, b, thick, local);
}

9
src/ASM/ASMs2Dmx.h
View File

@@ -103,8 +103,9 @@ public:
//! \param[in] revers Indicates whether the two edges have opposite directions
//! \param[in] basis The basis to connect (for mixed problems)
//! \param[in] coordCheck False to disable coordinate checks (periodic connections)
//! \param[in] thick Thickness of connection
virtual bool connectPatch(int edge, ASMs2D& neighbor, int nedge, bool revers,
int basis = 0, bool coordCheck = true);
int basis = 0, bool coordCheck = true, int thick = 1);
//! \brief Makes two opposite boundary edges periodic.
//! \param[in] dir Parameter direction defining the periodic edges
@@ -212,8 +213,10 @@ public:
//! \param[in] lIndex Local index of the boundary edge
//! \param nodes Array of node numbers
//! \param basis Which basis to grab nodes for (0 for all)
virtual void getBoundaryNodes(int lIndex, IntVec& nodes,
int basis, bool local = false) const;
//! \param thick Thickness of connection
//! \param local If \e true return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes, int basis,
int thick = 1, bool local = false) const;
protected:
std::vector<std::shared_ptr<Go::SplineSurface>> m_basis; //!< Vector of bases

185
src/ASM/ASMs3D.C
View File

@@ -617,7 +617,7 @@ bool ASMs3D::assignNodeNumbers (BlockNodes& nodes, int basis)
bool ASMs3D::connectPatch (int face, ASMs3D& neighbor, int nface,
int norient, int, bool coordCheck)
int norient, int, bool coordCheck, int thick)
{
if (swapW && face > 4) // Account for swapped parameter direction
face = 11-face;
@@ -625,7 +625,7 @@ bool ASMs3D::connectPatch (int face, ASMs3D& neighbor, int nface,
if (neighbor.swapW && nface > 4) // Account for swapped parameter direction
nface = 11-nface;
if (!this->connectBasis(face,neighbor,nface,norient,1,0,0,coordCheck))
if (!this->connectBasis(face,neighbor,nface,norient,1,0,0,coordCheck,thick))
return false;
this->addNeighbor(&neighbor);
@@ -634,7 +634,8 @@ bool ASMs3D::connectPatch (int face, ASMs3D& neighbor, int nface,
bool ASMs3D::connectBasis (int face, ASMs3D& neighbor, int nface, int norient,
int basis, int slave, int master, bool coordCheck)
int basis, int slave, int master,
bool coordCheck, int thick)
{
if (this->isShared() && neighbor.isShared())
return true;
@@ -643,94 +644,38 @@ bool ASMs3D::connectBasis (int face, ASMs3D& neighbor, int nface, int norient,
std::cerr <<" *** ASMs3D::connectPatch: Logic error, cannot"
<<" connect a shared patch with an unshared one"<< std::endl;
return false;
}
// Set up the slave node numbers for this volume patch
int n1, n2, n3;
if (!this->getSize(n1,n2,n3,basis)) return false;
int node = slave+1, i1 = 0, i2 = 0;
switch (face)
{
case 2: // Positive I-direction
node += n1-1;
case 1: // Negative I-direction
i1 = n1;
n1 = n2;
n2 = n3;
break;
case 4: // Positive J-direction
node += n1*(n2-1);
case 3: // Negative J-direction
i2 = n1*(n2-1);
i1 = 1;
n2 = n3;
break;
case 6: // Positive K-direction
node += n1*n2*(n3-1);
case 5: // Negative K-direction
i1 = 1;
break;
default:
std::cerr <<" *** ASMs3D::connectPatch: Invalid slave face "
<< face << std::endl;
return false;
}
int i, j;
IntMat slaveNodes(n1,IntVec(n2,0));
for (j = 0; j < n2; j++, node += i2)
for (i = 0; i < n1; i++, node += i1)
slaveNodes[i][j] = node;
// Set up the master node numbers for the neighboring volume patch
if (!neighbor.getSize(n1,n2,n3,basis)) return false;
node = master+1; i1 = i2 = 0;
switch (nface)
{
case 2: // Positive I-direction
node += n1-1;
case 1: // Negative I-direction
i1 = n1;
n1 = n2;
n2 = n3;
break;
case 4: // Positive J-direction
node += n1*(n2-1);
case 3: // Negative J-direction
i2 = n1*(n2-1);
i1 = 1;
n2 = n3;
break;
case 6: // Positive K-direction
node += n1*n2*(n3-1);
case 5: // Negative K-direction
i1 = 1;
break;
default:
std::cerr <<" *** ASMs3D::connectPatch: Invalid master face "
<< nface << std::endl;
return false;
}
if (norient < 0 || norient > 7)
{
} else if (face < 1 || face > 6) {
std::cerr <<" *** ASMs3D::connectPatch: Invalid slave face "
<< face << std::endl;
return false;
} else if (nface < 1 || nface > 6) {
std::cerr <<" *** ASMs3D::connectPatch: Invalid master face "
<< nface << std::endl;
return false;
} else if (norient < 0 || norient > 7) {
std::cerr <<" *** ASMs3D::connectPatch: Orientation flag "
<< norient <<" is out of range [0,7]"<< std::endl;
<< norient <<" is out of range [0,7]"<< std::endl;
return false;
}
int m1 = slaveNodes.size();
int m2 = slaveNodes.front().size();
int m1, m2;
if (!this->getFaceSize(m1,m2,basis,face)) return false;
int n1, n2;
if (!neighbor.getFaceSize(n1,n2,basis,nface)) return false;
// Set up the slave node numbers for this volume patch
IntVec slaveNodes;
this->getBoundaryNodes(face, slaveNodes, basis, thick, true);
for (int& it : slaveNodes)
it += slave;
// Set up the master node numbers for the neighboring volume patch
IntVec masterNodes;
neighbor.getBoundaryNodes(nface, masterNodes, basis, thick, true);
for (int& it : masterNodes)
it += master;
if (norient < 4 ? (n1 != m1 || n2 != m2) : (n2 != m1 || n1 != m2))
{
std::cerr <<" *** ASMs3D::connectPatch: Non-matching faces, sizes "
@@ -739,8 +684,9 @@ bool ASMs3D::connectBasis (int face, ASMs3D& neighbor, int nface, int norient,
}
const double xtol = 1.0e-4;
for (j = 0; j < n2; j++, node += i2)
for (i = 0; i < n1; i++, node += i1)
int node = 1;
for (int j = 0; j < n2; j++)
for (int i = 0; i < n1; i++, ++node)
{
int k = i, l = j;
switch (norient)
@@ -755,18 +701,22 @@ bool ASMs3D::connectBasis (int face, ASMs3D& neighbor, int nface, int norient,
default: k = i ; l = j ;
}
int slave = slaveNodes[k][l];
if (!coordCheck)
ASMbase::collapseNodes(neighbor,node,*this,slave);
else if (neighbor.getCoord(node).equal(this->getCoord(slave),xtol))
ASMbase::collapseNodes(neighbor,node,*this,slave);
else
for (int t = 0; t < thick; ++t)
{
std::cerr <<" *** ASMs3D::connectPatch: Non-matching nodes "
<< node <<": "<< neighbor.getCoord(node)
<<"\n and "
<< slave <<": "<< this->getCoord(slave) << std::endl;
return false;
int slave = slaveNodes[(l*n1+k)*thick+t];
int node2 = masterNodes[(node-1)*thick+t];
if (!coordCheck)
ASMbase::collapseNodes(neighbor,node2,*this,slave);
else if (neighbor.getCoord(node2).equal(this->getCoord(slave),xtol))
ASMbase::collapseNodes(neighbor,node2,*this,slave);
else
{
std::cerr <<" *** ASMs3D::connectPatch: Non-matching nodes "
<< node2 <<": "<< neighbor.getCoord(node2)
<<"\n and "
<< slave <<": "<< this->getCoord(slave) << std::endl;
return false;
}
}
}
@@ -1543,7 +1493,7 @@ bool ASMs3D::updateCoords (const Vector& displ)
void ASMs3D::getBoundaryNodes (int lIndex, IntVec& nodes,
int basis, bool local) const
int basis, int thick, bool local) const
{
if (basis == 0)
basis = 1;
@@ -1556,30 +1506,36 @@ void ASMs3D::getBoundaryNodes (int lIndex, IntVec& nodes,
int n1, n2, n3;
int node = this->findStartNode(n1,n2,n3,basis);
if (local)
node = 1;
switch (lIndex)
{
case 2: // Right face (positive I-direction)
node += n1-1;
node += n1-thick;
case 1: // Left face (negative I-direction)
for (int i3 = 1; i3 <= n3; i3++)
for (int i2 = 1; i2 <= n2; i2++, node += n1)
nodes.push_back(local ? node : this->getNodeID(node));
for (int t = 0; t < thick; ++t)
nodes.push_back(local ? node+t : this->getNodeID(node+t));
break;
case 4: // Back face (positive J-direction)
node += n1*(n2-1);
node += n1*(n2-thick);
case 3: // Front face (negative J-direction)
for (int i3 = 1; i3 <= n3; i3++, node += n1*(n2-1))
for (int i1 = 1; i1 <= n1; i1++, node++)
nodes.push_back(local ? node : this->getNodeID(node));
for (int t = 0; t < thick; ++t)
nodes.push_back(local ? node+t*n1 : this->getNodeID(node+t*n1));
break;
case 6: // Top face (positive K-direction)
node += n1*n2*(n3-1);
node += n1*n2*(n3-thick);
case 5: // Bottom face (negative K-direction)
for (int i2 = 1; i2 <= n2; i2++)
for (int i1 = 1; i1 <= n1; i1++, node++)
nodes.push_back(local ? node : this->getNodeID(node));
for (int t = 0; t < thick; ++t)
nodes.push_back(local ? node+t*n1*n2 : this->getNodeID(node+t*n1*n2));
break;
}
@@ -3240,3 +3196,18 @@ int ASMs3D::getCorner (int I, int J, int K, int basis) const
return node;
}
bool ASMs3D::getFaceSize(int& n1, int& n2, int basis, int face) const
{
int n3;
if (!this->getSize(n1,n2,n3,basis))
return false;
if (face == 1 || face == 2)
n1 = n2, n2 = n3;
else if (face == 3 || face == 4)
n2 = n3;
return true;
}

21
src/ASM/ASMs3D.h
View File

@@ -199,9 +199,11 @@ public:
//! \param[in] lIndex Local index of the boundary face
//! \param nodes Array of node numbers
//! \param basis Which basis to grab nodes for (0 for all)
//! \param local If \e true, return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& glbNodes, int basis,
bool local) const;
//! \param thick Thickness of connection
//! \param local If true return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes,
int basis, int thick = 1,
bool local = false) const;
//! \brief Returns the node index for a given corner.
virtual int getCorner(int I, int J, int K, int basis = 1) const;
@@ -326,6 +328,7 @@ public:
//! \param[in] nface Local face index of neighbor patch, in range [1,6]
//! \param[in] norient Relative face orientation flag (see below)
//! \param[in] coordCheck False to disable coordinate checks (periodic connections)
//! \param[in] thick Thickness of connection
//!
//! \details The face orientation flag \a norient must be in range [0,7].
//! When interpreted as a binary number, its 3 digits are decoded as follows:
@@ -333,7 +336,7 @@ public:
//! - middle digit = 1: Parameter \a u in neighbor patch face is reversed
//! - right digit = 1: Parameter \a v in neighbor patch face is reversed
virtual bool connectPatch(int face, ASMs3D& neighbor, int nface, int norient,
int = 0, bool coordCheck = true);
int = 0, bool coordCheck = true, int thick = 1);
//! \brief Makes two opposite boundary faces periodic.
//! \param[in] dir Parameter direction defining the periodic faces
@@ -544,9 +547,10 @@ protected:
//! \param[in] slave 0-based index of the first slave node in this basis
//! \param[in] master 0-based index of the first master node in this basis
//! \param[in] coordCheck False to turn off coordinate checks
//! \param[in] thick Thickness of connection
bool connectBasis(int face, ASMs3D& neighbor, int nface, int norient,
int basis = 1, int slave = 0, int master = 0,
bool coordCheck = true);
bool coordCheck = true, int thick = 1);
//! \brief Extracts parameter values of the Gauss points in one direction.
//! \param[out] uGP Parameter values in given direction for all points
@@ -683,6 +687,13 @@ private:
//! \return 1-based index of start node for basis
int findStartNode(int& n1, int& n2, int& n3, char basis) const;
//! \brief Find local sizes for a given face.
//! \param[out] n1 Number of nodes in first local parameter direction on face
//! \param[out] n2 Number of nodes in second local parameter direction face
//! \param[in] basis Basis to obtain sizes for
//! \param[in] face Face to obtain sizes for
bool getFaceSize(int& n1, int& n2, int basis, int face) const;
protected:
Go::SplineVolume* svol; //!< Pointer to the actual spline volume object
bool swapW; //!< Has the w-parameter direction been swapped?

12
src/ASM/ASMs3Dmx.C
View File

@@ -299,7 +299,7 @@ bool ASMs3Dmx::generateFEMTopology ()
bool ASMs3Dmx::connectPatch (int face, ASMs3D& neighbor, int nface, int norient,
int basis, bool coordCheck)
int basis, bool coordCheck, int thick)
{
ASMs3Dmx* neighMx = dynamic_cast<ASMs3Dmx*>(&neighbor);
if (!neighMx) return false;
@@ -313,7 +313,7 @@ bool ASMs3Dmx::connectPatch (int face, ASMs3D& neighbor, int nface, int norient,
size_t nb1 = 0, nb2 = 0;
for (size_t i = 1; i <= m_basis.size(); ++i) {
if (basis == 0 || i == (size_t)basis)
if (!this->connectBasis(face,neighbor,nface,norient,i,nb1,nb2,coordCheck))
if (!this->connectBasis(face,neighbor,nface,norient,i,nb1,nb2,coordCheck,thick))
return false;
nb1 += nb[i-1];
@@ -1109,12 +1109,12 @@ double ASMs3Dmx::getParametricArea (int iel, int dir) const
}
void ASMs3Dmx::getBoundaryNodes (int lIndex, IntVec& nodes,
int basis, bool local) const
void ASMs3Dmx::getBoundaryNodes (int lIndex, IntVec& nodes, int basis,
int thick, bool local) const
{
if (basis > 0)
this->ASMs3D::getBoundaryNodes(lIndex, nodes, basis, local);
this->ASMs3D::getBoundaryNodes(lIndex, nodes, basis, thick, local);
else
for (size_t b = 1; b <= this->getNoBasis(); ++b)
this->ASMs3D::getBoundaryNodes(lIndex, nodes, b, local);
this->ASMs3D::getBoundaryNodes(lIndex, nodes, b, thick, local);
}

12
src/ASM/ASMs3Dmx.h
View File

@@ -96,8 +96,9 @@ public:
//! \param[in] norient Relative face orientation flag (see class ASMs3D)
//! \param[in] basis Which basis to connect, or 0 for all.
//! \param[in] coordCheck False to disable coordinate checks (periodic connections)
//! \param[in] thick Thickness of connection
virtual bool connectPatch(int face, ASMs3D& neighbor, int nface, int norient,
int basis = 0, bool coordCheck = true);
int basis = 0, bool coordCheck = true, int thick = 1);
//! \brief Makes two opposite boundary faces periodic.
//! \param[in] dir Parameter direction defining the periodic faces
@@ -216,11 +217,12 @@ protected:
//! \brief Finds the global (or patch-local) node numbers on a patch boundary.
//! \param[in] lIndex Local index of the boundary edge
//! \param nodes Array of global boundary node numbers
//! \param nodes Array of node numbers
//! \param basis Which basis to grab nodes for (0 for all)
//! \param local If \e true, return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes,
int basis, bool local) const;
//! \param thick Thickness of connection
//! \param local If \e true return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes, int basis = 0,
int thick = 1, bool local = false) const;
private:
std::vector<std::shared_ptr<Go::SplineVolume>> m_basis; //!< Vector of bases

65
src/ASM/DomainDecomposition.C
View File

@@ -295,7 +295,7 @@ std::vector<IntVec> DomainDecomposition::calcSubdomains3D(size_t nel1, size_t ne
void DomainDecomposition::setupNodeNumbers(int basis, IntVec& lNodes,
std::set<int>& cbasis,
const ASMbase* pch,
int dim, int lidx)
int dim, int lidx, int thick)
{
if (basis != 0) // specified base
cbasis = utl::getDigits(basis);
@@ -304,7 +304,7 @@ void DomainDecomposition::setupNodeNumbers(int basis, IntVec& lNodes,
for (size_t b = 1; b <= pch->getNoBasis(); ++b)
cbasis.insert(b);
else // directly add nodes, cbasis remains empty
pch->getBoundaryNodes(lidx, lNodes);
pch->getBoundaryNodes(lidx, lNodes, 0, thick, false);
const ASM2D* pch2D = dynamic_cast<const ASM2D*>(pch);
const ASM3D* pch3D = dynamic_cast<const ASM3D*>(pch);
@@ -335,7 +335,7 @@ void DomainDecomposition::setupNodeNumbers(int basis, IntVec& lNodes,
for (const int& it : eNodes)
lNodes.push_back(pch->getNodeID(it));
} else
pch->getBoundaryNodes(lidx, lNodes, it2);
pch->getBoundaryNodes(lidx, lNodes, it2, thick, false);
}
@@ -393,7 +393,8 @@ bool DomainDecomposition::calcGlobalNodeNumbers(const ProcessAdm& adm,
std::set<int> cbasis;
IntVec lNodes;
setupNodeNumbers(it.basis, lNodes, cbasis, sim.getPatch(sidx), it.dim, it.sidx);
setupNodeNumbers(it.basis, lNodes, cbasis, sim.getPatch(sidx),
it.dim, it.sidx, it.thick);
int nRecv;
adm.receive(nRecv, getPatchOwner(it.master));
@@ -412,10 +413,12 @@ bool DomainDecomposition::calcGlobalNodeNumbers(const ProcessAdm& adm,
it.orient, it.sidx, b, it.dim);
auto it_n = iter.begin();
for (size_t i = 0; i < iter.size(); ++i, ++it_n) {
int node = MLGN[lNodes[i+ofs]-1];
old2new[node] = glbNodes[*it_n + ofs];
for (int t = 0; t < it.thick; ++t) {
int node = MLGN[lNodes[i*it.thick+t+ofs]-1];
old2new[node] = glbNodes[*it_n*it.thick+t + ofs];
}
}
ofs += iter.size();
ofs += iter.size()*it.thick;
}
}
}
@@ -455,7 +458,8 @@ bool DomainDecomposition::calcGlobalNodeNumbers(const ProcessAdm& adm,
std::set<int> cbasis;
std::vector<int> glbNodes;
setupNodeNumbers(it.basis, glbNodes, cbasis, sim.getPatch(midx), it.dim, it.midx);
setupNodeNumbers(it.basis, glbNodes, cbasis, sim.getPatch(midx),
it.dim, it.midx, it.thick);
for (size_t i = 0; i < glbNodes.size(); ++i)
glbNodes[i] = MLGN[glbNodes[i]-1];
@@ -472,7 +476,7 @@ bool DomainDecomposition::calcGlobalNodeNumbers(const ProcessAdm& adm,
std::vector<int> DomainDecomposition::setupEquationNumbers(const SIMbase& sim,
int pidx, int lidx,
const std::set<int>& cbasis,
int dim)
int dim, int thick)
{
std::vector<IntVec> lNodes(sim.getPatch(pidx)->getNoBasis());
std::vector<int> result;
@@ -491,7 +495,8 @@ std::vector<int> DomainDecomposition::setupEquationNumbers(const SIMbase& sim,
if (lNodes[basis-1].empty()) {
IntSet dummy;
setupNodeNumbers(basis, lNodes[basis-1], dummy, sim.getPatch(pidx), dim, lidx);
setupNodeNumbers(basis, lNodes[basis-1], dummy, sim.getPatch(pidx),
dim, lidx, thick);
}
std::set<int> components;
@@ -616,7 +621,8 @@ bool DomainDecomposition::calcGlobalEqNumbers(const ProcessAdm& adm,
if (it.basis != 0)
cbasis = utl::getDigits(it.basis);
IntVec locEqs = setupEquationNumbers(sim, sidx, it.sidx, cbasis, it.dim);
IntVec locEqs = setupEquationNumbers(sim, sidx, it.sidx,
cbasis, it.dim, it.thick);
int nRecv;
adm.receive(nRecv, getPatchOwner(it.master));
@@ -655,25 +661,27 @@ bool DomainDecomposition::calcGlobalEqNumbers(const ProcessAdm& adm,
auto it_n = iter.begin();
int nodeDofs = components.size();
for (size_t i = 0; i < iter.size(); ++i, ++it_n) {
for (size_t comp = 0; comp < components.size(); ++comp) {
int leq = locEqs[ofs+i*nodeDofs+comp];
int geq = glbEqs[ofs+(*it_n)*nodeDofs+comp];
if ((leq < 1 && geq > 0) || (leq > 0 && geq < 1)) {
std::cerr <<"\n *** DomainDecomposition::calcGlobalEqNumbers(): "
<< "Topology error on process " << adm.getProcId()
<< ", dof constraint mismatch "
<< "local: " << leq << " global " << geq << std::endl;
return false;
for (int t = 0; t < it.thick; ++t) {
for (size_t comp = 0; comp < components.size(); ++comp) {
int leq = locEqs[ofs+(i*it.thick+t)*nodeDofs+comp];
int geq = glbEqs[ofs+(*it_n*it.thick+t)*nodeDofs+comp];
if ((leq < 1 && geq > 0) || (leq > 0 && geq < 1)) {
std::cerr <<"\n *** DomainDecomposition::calcGlobalEqNumbers(): "
<< "Topology error on process " << adm.getProcId()
<< ", dof constraint mismatch "
<< "local: " << leq << " global " << geq << std::endl;
return false;
}
if (leq < 1)
continue;
old2new[block][leq] = geq;
o2nu[block][leq-nEqs[block]-1] = true;
}
if (leq < 1)
continue;
old2new[block][leq] = geq;
o2nu[block][leq-nEqs[block]-1] = true;
}
}
ofs += iter.size()*nodeDofs;
ofs += iter.size()*nodeDofs*it.thick;
}
}
}
@@ -732,7 +740,8 @@ bool DomainDecomposition::calcGlobalEqNumbers(const ProcessAdm& adm,
if (it.basis != 0)
cbasis = utl::getDigits(it.basis);
IntVec glbEqs = setupEquationNumbers(sim, midx, it.midx, cbasis, it.dim);
IntVec glbEqs = setupEquationNumbers(sim, midx, it.midx,
cbasis, it.dim, it.thick);
adm.send(int(glbEqs.size()), getPatchOwner(it.slave));
adm.send(glbEqs, getPatchOwner(it.slave));
}

7
src/ASM/DomainDecomposition.h
View File

@@ -41,6 +41,7 @@ public:
int orient; //!< Orientation.
int dim; //!< Dimension of boundary.
int basis; //!< Basis of boundary.
int thick; //!< Thickness of connection.
};
//! \brief Functor to order ghost connections.
@@ -184,10 +185,11 @@ private:
//! \param pidx Patch index
//! \param lidx Boundary index on patch
//! \param cbasis If non-empty, bases to connect
//! \param thick Thickness of connection (subdivisions)
std::vector<int> setupEquationNumbers(const SIMbase& sim,
int pidx, int lidx,
const std::set<int>& cbasis,
int dim);
int dim, int thick);
//! \brief Setup node numbers for all bases on a boundary.
//! \param basis Bases to grab nodes for
@@ -196,10 +198,11 @@ private:
//! \param pch Patch to obtain nodes for
//! \param dim Dimension of boundary to obtain nodes for
//! \param lidx Local index of boundary to obtain nodes for
//! \param thick Thickness of connection (subdivisions)
void setupNodeNumbers(int basis, std::vector<int>& lNodes,
std::set<int>& cbasis,
const ASMbase* pch,
int dim, int lidx);
int dim, int lidx, int thick);
//! \brief Calculate the global node numbers for given finite element model.
bool calcGlobalNodeNumbers(const ProcessAdm& adm, const SIMbase& sim);

2
src/ASM/LR/ASMu2D.C
View File

@@ -1813,7 +1813,7 @@ bool ASMu2D::evalSolution (Matrix& sField, const IntegrandBase& integrand,
void ASMu2D::getBoundaryNodes (int lIndex, IntVec& nodes, int basis,
bool local) const
int, bool local) const
{
if (basis == 0)
basis = 1;

5
src/ASM/LR/ASMu2D.h
View File

@@ -94,9 +94,10 @@ public:
//! \brief Finds the global (or patch-local) node numbers on a patch boundary.
//! \param[in] lIndex Local index of the boundary edge
//! \param nodes Array of node numbers
//! \param local If true return patch-local node numbers
//! \param basis Which basis to grab nodes for (0 for all)
//! \param local If \e true return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes,
int basis, bool local = false) const;
int basis, int, bool local = false) const;
//! \brief Returns the polynomial order in each parameter direction.
//! \param[out] p1 Order in first (u) direction

2
src/ASM/LR/ASMu3D.C
View File

@@ -1966,7 +1966,7 @@ std::vector<int> ASMu3D::getFaceNodes (int face, int basis) const
}
void ASMu3D::getBoundaryNodes (int lIndex, IntVec& nodes, int basis, bool) const
void ASMu3D::getBoundaryNodes (int lIndex, IntVec& nodes, int basis, int, bool) const
{
if (basis == 0)
basis = 1;

4
src/ASM/LR/ASMu3D.h
View File

@@ -89,10 +89,10 @@ public:
//! \brief Finds the global (or patch-local) node numbers on a patch boundary.
//! \param[in] lIndex Local index of the boundary edge
//! \param nodes Array of global boundary node numbers
//! \param glbNodes Array of global boundary node numbers
//! \param local If \e true return patch-local node numbers
virtual void getBoundaryNodes(int lIndex, IntVec& nodes,
int basis, bool local) const;
int basis, int, bool local) const;
//! \brief Returns the node index for a given corner.
virtual int getCorner(int I, int J, int K, int basis = 1) const;

8
src/SIM/SIM2D.C
View File

@@ -69,7 +69,7 @@ SIM2D::SIM2D (IntegrandBase* itg, unsigned char n, bool check) : SIMgeneric(itg)
bool SIM2D::addConnection (int master, int slave, int mIdx,
int sIdx, int orient, int basis,
bool coordCheck, int dim)
bool coordCheck, int dim, int thick)
{
if (orient < 0 || orient > 1)
{
@@ -99,13 +99,13 @@ bool SIM2D::addConnection (int master, int slave, int mIdx,
bases = utl::getDigits(basis);
for (const int& b : bases)
if (!spch->connectPatch(sIdx,*mpch,mIdx,orient,b,coordCheck))
if (!spch->connectPatch(sIdx,*mpch,mIdx,orient,b,coordCheck,thick))
return false;
}
else
adm.dd.ghostConnections.insert(DomainDecomposition::Interface{master, slave,
mIdx, sIdx,
orient, dim, basis});
mIdx, sIdx, orient,
dim, basis, thick});
return true;
}

7
src/SIM/SIM2D.h
View File

@@ -80,9 +80,10 @@ public:
//! \param[in] basis Which bases to connect (0 for all)
//! \param[in] coordCheck If \e false, do not check for matching coordinates
//! \param[in] dim Dimensionality of connection
virtual bool addConnection(int master, int slave, int mEdge, int sEdge,
int orient, int basis = 0, bool coordCheck = true,
int dim = 1);
//! \param thick Thickness of connection
bool addConnection(int master, int slave, int mIdx, int sIdx,
int orient, int basis = 0, bool coordCheck = true,
int dim = 1, int thick = 1);
//! \brief Evaluates the primary solution at the given point.
//! \param[in] psol Primary solution vector

8
src/SIM/SIM3D.C
View File

@@ -47,7 +47,7 @@ SIM3D::SIM3D (IntegrandBase* itg, unsigned char n, bool check) : SIMgeneric(itg)
bool SIM3D::addConnection (int master, int slave, int mIdx,
int sIdx, int orient, int basis,
bool coordCheck, int dim)
bool coordCheck, int dim, int thick)
{
if (orient < 0 || orient > 7)
{
@@ -77,13 +77,13 @@ bool SIM3D::addConnection (int master, int slave, int mIdx,
bases = utl::getDigits(basis);
for (const int& b : bases)
if (!spch->connectPatch(sIdx,*mpch,mIdx,orient,b,coordCheck))
if (!spch->connectPatch(sIdx,*mpch,mIdx,orient,b,coordCheck,thick))
return false;
}
else
adm.dd.ghostConnections.insert(DomainDecomposition::Interface{master, slave,
mIdx, sIdx,
orient, dim, basis});
mIdx, sIdx, orient,
dim, basis, thick});
return true;
}

3
src/SIM/SIM3D.h
View File

@@ -74,9 +74,10 @@ public:
//! \param[in] basis Which bases to connect (0 for all)
//! \param[in] coordCheck If \e false, do not check for matching coordinates
//! \param[in] dim Dimensionality of connection
//! \param thick Thickness of connection
virtual bool addConnection(int master, int slave, int mFace, int sFace,
int orient, int basis = 0, bool coordCheck = true,
int dim = 2);
int dim = 2, int thick = 1);
//! \brief Evaluates the primary solution at the given point.
//! \param[in] psol Primary solution vector

3
src/SIM/SIMinput.h
View File

@@ -184,9 +184,10 @@ public:
//! \param[in] basis Which bases to connect (0 for all)
//! \param[in] coordCheck If \e false, do not check for matching coordinates
//! \param[in] dim Dimensionality of connection
//! \param[in] thick Thickness of connection
virtual bool addConnection(int master, int slave, int mIdx, int sIdx,
int orient, int basis = 0, bool coordCheck = true,
int dim = 1) { return false; }
int dim = 1, int thick = 1) { return false; }
protected:
//! \brief Reads global node data for a patch from given input stream.