[cig-commits] r18424 - in short/3D/PyLith/trunk: . libsrc/pylith/topology
knepley at geodynamics.org
knepley at geodynamics.org
Sun May 22 16:16:42 PDT 2011
Author: knepley
Date: 2011-05-22 16:16:42 -0700 (Sun, 22 May 2011)
New Revision: 18424
Modified:
short/3D/PyLith/trunk/TODO
short/3D/PyLith/trunk/libsrc/pylith/topology/MeshRefiner.cc
short/3D/PyLith/trunk/libsrc/pylith/topology/RefineFace4Edges2.cc
short/3D/PyLith/trunk/libsrc/pylith/topology/RefineVol8Face4Edges2.cc
Log:
Fixed bugs with vertex offsets in refinement
Modified: short/3D/PyLith/trunk/TODO
===================================================================
--- short/3D/PyLith/trunk/TODO 2011-05-22 19:05:20 UTC (rev 18423)
+++ short/3D/PyLith/trunk/TODO 2011-05-22 23:16:42 UTC (rev 18424)
@@ -132,6 +132,13 @@
* Scalable distribution [MATT]
+ + It appears that topology adjustment is not currently limiting the runs
+ + Need a more memory scalable Distribution
+ + Consider simple distribution, followed by cleanup
+ + I think the Overlap structure is probably taking up all the memory
+ + send/recvMeshOverlap is only used in completeConesV() to call SimpleCopy::copy()
+ + Can we use a lighter-weight structure to sit in copy()?
+
----------------------------------------------------------------------
SECONDARY PRIORITIES
----------------------------------------------------------------------
Modified: short/3D/PyLith/trunk/libsrc/pylith/topology/MeshRefiner.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/topology/MeshRefiner.cc 2011-05-22 19:05:20 UTC (rev 18423)
+++ short/3D/PyLith/trunk/libsrc/pylith/topology/MeshRefiner.cc 2011-05-22 23:16:42 UTC (rev 18424)
@@ -673,21 +673,28 @@
// Create the parallel overlap
// Get offsets for points across processors and add points in overlap from original mesh to new overlap
- int* oldVerticesStartP = (mesh->commSize() > 0) ? new int[mesh->commSize()] : 0;
- int* newVerticesStartP = (newMesh->commSize() > 0) ? new int[newMesh->commSize()] : 0;
+ int* oldVerticesStartNormalP = (mesh->commSize() > 0) ? new int[mesh->commSize()] : 0;
+ int* oldVerticesStartFaultP = (mesh->commSize() > 0) ? new int[mesh->commSize()] : 0;
+ int* newVerticesStartNormalP = (newMesh->commSize() > 0) ? new int[newMesh->commSize()] : 0;
+ int* newVerticesStartFaultP = (newMesh->commSize() > 0) ? new int[newMesh->commSize()] : 0;
int ierr = 0;
- const int oldVerticesStart = _orderOldMesh->verticesNormal().min();
- const int newVerticesStart = _orderNewMesh->verticesNormal().min();
+ const int oldVerticesStartNormal = _orderOldMesh->verticesNormal().min();
+ const int oldVerticesStartFault = _orderOldMesh->verticesCensored().min();
+ const int newVerticesStartNormal = _orderNewMesh->verticesNormal().min();
+ const int newVerticesStartFault = _orderNewMesh->verticesCensored().min();
- ierr = MPI_Allgather((void *) &oldVerticesStart, 1, MPI_INT, oldVerticesStartP, 1, MPI_INT, mesh->comm());CHKERRXX(ierr);
- ierr = MPI_Allgather((void *) &newVerticesStart, 1, MPI_INT, newVerticesStartP, 1, MPI_INT, newMesh->comm());CHKERRXX(ierr);
+ ierr = MPI_Allgather((void *) &oldVerticesStartNormal, 1, MPI_INT, oldVerticesStartNormalP, 1, MPI_INT, mesh->comm());CHKERRXX(ierr);
+ ierr = MPI_Allgather((void *) &oldVerticesStartFault, 1, MPI_INT, oldVerticesStartFaultP, 1, MPI_INT, mesh->comm());CHKERRXX(ierr);
+ ierr = MPI_Allgather((void *) &newVerticesStartNormal, 1, MPI_INT, newVerticesStartNormalP, 1, MPI_INT, newMesh->comm());CHKERRXX(ierr);
+ ierr = MPI_Allgather((void *) &newVerticesStartFault, 1, MPI_INT, newVerticesStartFaultP, 1, MPI_INT, newMesh->comm());CHKERRXX(ierr);
Obj<mesh_type::send_overlap_type> newSendOverlap = newMesh->getSendOverlap();
Obj<mesh_type::recv_overlap_type> newRecvOverlap = newMesh->getRecvOverlap();
const Obj<mesh_type::send_overlap_type>& sendOverlap = mesh->getSendOverlap();
const Obj<mesh_type::recv_overlap_type>& recvOverlap = mesh->getRecvOverlap();
Obj<mesh_type::send_overlap_type::traits::capSequence> sendPoints = sendOverlap->cap();
- const mesh_type::send_overlap_type::source_type localOffset = newVerticesStart - oldVerticesStart;
+ const mesh_type::send_overlap_type::source_type localOffsetNormal = newVerticesStartNormal - oldVerticesStartNormal;
+ const mesh_type::send_overlap_type::source_type localOffsetFault = newVerticesStartFault - oldVerticesStartFault;
for(mesh_type::send_overlap_type::traits::capSequence::iterator p_iter = sendPoints->begin(); p_iter != sendPoints->end(); ++p_iter) {
const Obj<mesh_type::send_overlap_type::traits::supportSequence>& ranks = sendOverlap->support(*p_iter);
@@ -696,9 +703,17 @@
for(mesh_type::send_overlap_type::traits::supportSequence::iterator r_iter = ranks->begin(); r_iter != ranks->end(); ++r_iter) {
const int rank = *r_iter;
const mesh_type::send_overlap_type::source_type& remotePoint = r_iter.color();
- const mesh_type::send_overlap_type::source_type remoteOffset = newVerticesStartP[rank] - oldVerticesStartP[rank];
+ const mesh_type::send_overlap_type::source_type remoteOffsetNormal = newVerticesStartNormalP[rank] - oldVerticesStartNormalP[rank];
+ const mesh_type::send_overlap_type::source_type remoteOffsetFault = newVerticesStartFaultP[rank] - oldVerticesStartFaultP[rank];
- newSendOverlap->addArrow(localPoint+localOffset, rank, remotePoint+remoteOffset);
+ if (localPoint >= oldVerticesStartNormal && localPoint < oldVerticesStartFault) {
+ assert(remotePoint >= oldVerticesStartNormalP[rank] && remotePoint < oldVerticesStartFaultP[rank]);
+ newSendOverlap->addArrow(localPoint+localOffsetNormal, rank, remotePoint+remoteOffsetNormal);
+ } else {
+ assert(localPoint >= oldVerticesStartFault);
+ assert(remotePoint >= oldVerticesStartFaultP[rank]);
+ newSendOverlap->addArrow(localPoint+localOffsetFault, rank, remotePoint+remoteOffsetFault);
+ }
} // for
} // for
Obj<mesh_type::recv_overlap_type::traits::baseSequence> recvPoints = recvOverlap->base();
@@ -710,14 +725,24 @@
for(mesh_type::recv_overlap_type::traits::coneSequence::iterator r_iter = ranks->begin(); r_iter != ranks->end(); ++r_iter) {
const int rank = *r_iter;
const mesh_type::recv_overlap_type::target_type& remotePoint = r_iter.color();
- const mesh_type::recv_overlap_type::target_type remoteOffset = newVerticesStartP[rank] - oldVerticesStartP[rank];
+ const mesh_type::send_overlap_type::source_type remoteOffsetNormal = newVerticesStartNormalP[rank] - oldVerticesStartNormalP[rank];
+ const mesh_type::send_overlap_type::source_type remoteOffsetFault = newVerticesStartFaultP[rank] - oldVerticesStartFaultP[rank];
- newRecvOverlap->addArrow(rank, localPoint+localOffset, remotePoint+remoteOffset);
+ if (localPoint >= oldVerticesStartNormal && localPoint < oldVerticesStartFault) {
+ assert(remotePoint >= oldVerticesStartNormalP[rank] && remotePoint < oldVerticesStartFaultP[rank]);
+ newRecvOverlap->addArrow(rank, localPoint+localOffsetNormal, remotePoint+remoteOffsetNormal);
+ } else {
+ assert(localPoint >= oldVerticesStartFault);
+ assert(remotePoint >= oldVerticesStartFaultP[rank]);
+ newRecvOverlap->addArrow(rank, localPoint+localOffsetFault, remotePoint+remoteOffsetFault);
+ }
} // for
} // for
newMesh->setCalculatedOverlap(true);
- delete [] oldVerticesStartP; oldVerticesStartP = PETSC_NULL;
- delete [] newVerticesStartP; newVerticesStartP = PETSC_NULL;
+ delete [] oldVerticesStartNormalP; oldVerticesStartNormalP = PETSC_NULL;
+ delete [] oldVerticesStartFaultP; oldVerticesStartFaultP = PETSC_NULL;
+ delete [] newVerticesStartNormalP; newVerticesStartNormalP = PETSC_NULL;
+ delete [] newVerticesStartFaultP; newVerticesStartFaultP = PETSC_NULL;
refiner.overlapAddNewVertices(newMesh, *_orderNewMesh, mesh, *_orderOldMesh);
} // _calcNewOverlap
Modified: short/3D/PyLith/trunk/libsrc/pylith/topology/RefineFace4Edges2.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/topology/RefineFace4Edges2.cc 2011-05-22 19:05:20 UTC (rev 18423)
+++ short/3D/PyLith/trunk/libsrc/pylith/topology/RefineFace4Edges2.cc 2011-05-22 23:16:42 UTC (rev 18424)
@@ -241,7 +241,8 @@
assert(!oldRecvOverlap.isNull());
// Offset used when converting from old mesh numbering to new
- const int localOffset = orderNewMesh.verticesNormal().min() - orderOldMesh.verticesNormal().min();
+ const int localNormalOffset = orderNewMesh.verticesNormal().min() - orderOldMesh.verticesNormal().min();
+ const int localCensoredOffset = orderNewMesh.verticesCensored().min() - orderOldMesh.verticesCensored().min();
// Check edges in edgeToVertex for both endpoints sent to same process
// Put it in section with point being the lowest numbered vertex and value (other endpoint, new vertex)
@@ -265,6 +266,9 @@
std::insert_iterator<std::set<int> >(ranks, ranks.begin()));
if(ranks.size()) {
+ const int localOffset = orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) ? localNormalOffset : localCensoredOffset;
+
+ assert(orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) == orderOldMesh.verticesNormal().hasPoint(e_iter->first.second));
newVerticesSection->addFiberDimension(std::min(e_iter->first.first, e_iter->first.second)+localOffset, 1);
for(std::set<int>::const_iterator r_iter = ranks.begin(); r_iter != ranks.end(); ++r_iter) {
bndryEdgeToRank[e_iter->first].push_back(*r_iter);
@@ -283,6 +287,9 @@
value_type* values = (dim > 0) ? new value_type[dim] : 0;
for(std::map<EdgeType, std::vector<int> >::const_iterator e_iter = bndryEdgeToRank.begin(); e_iter != bndryEdgeToRank.end() && v < dim; ++e_iter) {
+ const int localOffset = orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) ? localNormalOffset : localCensoredOffset;
+
+ assert(orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) == orderOldMesh.verticesNormal().hasPoint(e_iter->first.second));
if (std::min(e_iter->first.first, e_iter->first.second)+localOffset == p) {
values[v++] = EdgeType(std::max(e_iter->first.first, e_iter->first.second)+localOffset, _edgeToVertex[e_iter->first]);
} // if
@@ -300,17 +307,19 @@
const point_type localPoint = _edgeToVertex[e_iter->first];
for(std::vector<int>::const_iterator r_iter = e_iter->second.begin(); r_iter != e_iter->second.end(); ++r_iter) {
- point_type remoteLeft = -1, remoteRight = -1;
- const int rank = *r_iter;
+ point_type remoteLeft = -1, remoteRight = -1;
+ const int rank = *r_iter;
+ const int localOffsetL = orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) ? localNormalOffset : localCensoredOffset;
- const Obj<mesh_type::send_overlap_type::traits::supportSequence>& leftRanks = newSendOverlap->support(e_iter->first.first+localOffset);
+ const Obj<mesh_type::send_overlap_type::traits::supportSequence>& leftRanks = newSendOverlap->support(e_iter->first.first+localOffsetL);
for(mesh_type::send_overlap_type::traits::supportSequence::iterator lr_iter = leftRanks->begin(); lr_iter != leftRanks->end(); ++lr_iter) {
if (rank == *lr_iter) {
remoteLeft = lr_iter.color();
break;
} // if
} // for
- const Obj<mesh_type::send_overlap_type::traits::supportSequence>& rightRanks = newSendOverlap->support(e_iter->first.second+localOffset);
+ const int localOffsetR = orderOldMesh.verticesNormal().hasPoint(e_iter->first.second) ? localNormalOffset : localCensoredOffset;
+ const Obj<mesh_type::send_overlap_type::traits::supportSequence>& rightRanks = newSendOverlap->support(e_iter->first.second+localOffsetR);
for(mesh_type::send_overlap_type::traits::supportSequence::iterator rr_iter = rightRanks->begin(); rr_iter != rightRanks->end(); ++rr_iter) {
if (rank == *rr_iter) {
remoteRight = rr_iter.color();
Modified: short/3D/PyLith/trunk/libsrc/pylith/topology/RefineVol8Face4Edges2.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/topology/RefineVol8Face4Edges2.cc 2011-05-22 19:05:20 UTC (rev 18423)
+++ short/3D/PyLith/trunk/libsrc/pylith/topology/RefineVol8Face4Edges2.cc 2011-05-22 23:16:42 UTC (rev 18424)
@@ -315,7 +315,8 @@
assert(!newVerticesSection.isNull());
std::map<EdgeType, std::vector<int> > bndryEdgeToRank;
- const int localOffset = orderNewMesh.verticesNormal().min() - orderOldMesh.verticesNormal().min();
+ const int localNormalOffset = orderNewMesh.verticesNormal().min() - orderOldMesh.verticesNormal().min();
+ const int localCensoredOffset = orderNewMesh.verticesCensored().min() - orderOldMesh.verticesCensored().min();
for(std::map<EdgeType, point_type>::const_iterator e_iter = _edgeToVertex.begin(); e_iter != _edgeToVertex.end(); ++e_iter) {
const point_type left = e_iter->first.first;
@@ -333,6 +334,9 @@
std::insert_iterator<std::set<int> >(ranks, ranks.begin()));
if(ranks.size()) {
+ const int localOffset = orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) ? localNormalOffset : localCensoredOffset;
+
+ assert(orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) == orderOldMesh.verticesNormal().hasPoint(e_iter->first.second));
newVerticesSection->addFiberDimension(std::min(e_iter->first.first, e_iter->first.second)+localOffset, 1);
for(std::set<int>::const_iterator r_iter = ranks.begin(); r_iter != ranks.end(); ++r_iter) {
bndryEdgeToRank[e_iter->first].push_back(*r_iter);
@@ -351,6 +355,9 @@
value_type* values = (dim > 0) ? new value_type[dim] : 0;
for(std::map<EdgeType, std::vector<int> >::const_iterator e_iter = bndryEdgeToRank.begin(); e_iter != bndryEdgeToRank.end() && v < dim; ++e_iter) {
+ const int localOffset = orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) ? localNormalOffset : localCensoredOffset;
+
+ assert(orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) == orderOldMesh.verticesNormal().hasPoint(e_iter->first.second));
if (std::min(e_iter->first.first, e_iter->first.second)+localOffset == p) {
values[v++] = EdgeType(std::max(e_iter->first.first, e_iter->first.second)+localOffset, _edgeToVertex[e_iter->first]);
} // if
@@ -399,7 +406,9 @@
minVertex = std::min(nextVertex, minVertex);
}
if (ranks.size()) {
- assert(f_iter->second >= orderNewMesh.verticesNormal().min() && f_iter->second < orderNewMesh.verticesNormal().max());
+ const int localOffset = orderOldMesh.verticesNormal().hasPoint(minVertex) ? localNormalOffset : localCensoredOffset;
+
+ assert(orderNewMesh.verticesNormal().hasPoint(f_iter->second));
newFaceVerticesSection->addFiberDimension(minVertex+localOffset, 1);
for(std::set<int>::const_iterator r_iter = ranks.begin(); r_iter != ranks.end(); ++r_iter) {
bndryFaceToRank[f_iter->first].push_back(*r_iter);
@@ -420,11 +429,13 @@
for(std::map<FaceType, std::vector<int> >::const_iterator f_iter = bndryFaceToRank.begin(); f_iter != bndryFaceToRank.end() && v < dim; ++f_iter) {
const point_type first = f_iter->first.points[0];
point_type minVertex = first;
+ int localOffset;
for(int i = 1; i < 4; ++i) {
const point_type nextVertex = f_iter->first.points[i];
minVertex = std::min(nextVertex, minVertex);
}
+ localOffset = orderOldMesh.verticesNormal().hasPoint(minVertex) ? localNormalOffset : localCensoredOffset;
if (minVertex+localOffset == p) {
FaceType face;
@@ -432,6 +443,8 @@
for(int i = 0; i < 4; ++i) {
if (f_iter->first.points[i] != minVertex) {
+ localOffset = orderOldMesh.verticesNormal().hasPoint(f_iter->first.points[i]) ? localNormalOffset : localCensoredOffset;
+
face.points[k++] = f_iter->first.points[i]+localOffset;
}
}
@@ -458,17 +471,19 @@
const point_type localPoint = _edgeToVertex[e_iter->first];
for(std::vector<int>::const_iterator r_iter = e_iter->second.begin(); r_iter != e_iter->second.end(); ++r_iter) {
- point_type remoteLeft = -1, remoteRight = -1;
- const int rank = *r_iter;
+ point_type remoteLeft = -1, remoteRight = -1;
+ const int rank = *r_iter;
+ const int localOffsetL = orderOldMesh.verticesNormal().hasPoint(e_iter->first.first) ? localNormalOffset : localCensoredOffset;
- const Obj<mesh_type::send_overlap_type::traits::supportSequence>& leftRanks = newSendOverlap->support(e_iter->first.first+localOffset);
+ const Obj<mesh_type::send_overlap_type::traits::supportSequence>& leftRanks = newSendOverlap->support(e_iter->first.first+localOffsetL);
for(mesh_type::send_overlap_type::traits::supportSequence::iterator lr_iter = leftRanks->begin(); lr_iter != leftRanks->end(); ++lr_iter) {
if (rank == *lr_iter) {
remoteLeft = lr_iter.color();
break;
} // if
} // for
- const Obj<mesh_type::send_overlap_type::traits::supportSequence>& rightRanks = newSendOverlap->support(e_iter->first.second+localOffset);
+ const int localOffsetR = orderOldMesh.verticesNormal().hasPoint(e_iter->first.second) ? localNormalOffset : localCensoredOffset;
+ const Obj<mesh_type::send_overlap_type::traits::supportSequence>& rightRanks = newSendOverlap->support(e_iter->first.second+localOffsetR);
for(mesh_type::send_overlap_type::traits::supportSequence::iterator rr_iter = rightRanks->begin(); rr_iter != rightRanks->end(); ++rr_iter) {
if (rank == *rr_iter) {
remoteRight = rr_iter.color();
@@ -487,6 +502,7 @@
break;
} // if
} // for
+ assert(remotePoint >= 0);
newSendOverlap->addArrow(localPoint, rank, remotePoint);
newRecvOverlap->addArrow(rank, localPoint, remotePoint);
} // for
@@ -500,6 +516,7 @@
const int rank = *r_iter;
for(int i = 0; i < 4; ++i) {
+ const int localOffset = orderOldMesh.verticesNormal().hasPoint(f_iter->first.points[i]) ? localNormalOffset : localCensoredOffset;
const Obj<mesh_type::send_overlap_type::traits::supportSequence>& faceRanks = newSendOverlap->support(f_iter->first.points[i]+localOffset);
for(mesh_type::send_overlap_type::traits::supportSequence::iterator fr_iter = faceRanks->begin(); fr_iter != faceRanks->end(); ++fr_iter) {
if (rank == *fr_iter) {
@@ -534,6 +551,7 @@
} // if
} // for
assert(localPoint >= orderNewMesh.verticesNormal().min() && localPoint < orderNewMesh.verticesNormal().max());
+ assert(remotePoint >= 0);
newSendOverlap->addArrow(localPoint, rank, remotePoint);
newRecvOverlap->addArrow(rank, localPoint, remotePoint);
} // for
More information about the CIG-COMMITS
mailing list