moab
|
#include <Intx2Mesh.hpp>
Definition at line 41 of file Intx2Mesh.hpp.
moab::Intx2Mesh::Intx2Mesh | ( | Interface * | mbimpl | ) |
Definition at line 22 of file Intx2Mesh.cpp.
: mb(mbimpl) #ifdef USE_MPI , parcomm(NULL), remote_cells(NULL) #endif { dbg_1=0; box_error=0; my_rank=0; BlueFlagTag=0; RedFlagTag=0; redParentTag =0; blueParentTag = 0; countTag = 0; }
moab::Intx2Mesh::~Intx2Mesh | ( | ) | [virtual] |
Definition at line 37 of file Intx2Mesh.cpp.
{ // TODO Auto-generated destructor stub #ifdef USE_MPI if (remote_cells) { delete remote_cells; remote_cells=NULL; } #endif }
ErrorCode moab::Intx2Mesh::build_processor_euler_boxes | ( | EntityHandle | euler_set, |
Range & | local_verts | ||
) |
void moab::Intx2Mesh::clean | ( | ) |
Definition at line 425 of file Intx2Mesh.cpp.
{ // int indx = 0; for (Range::iterator eit = RedEdges.begin(); eit != RedEdges.end(); eit++, indx++) { //EntityHandle edge = *eit; //delete extraNodesMap[edge]; delete extraNodesVec[indx]; } //extraNodesMap.clear(); extraNodesVec.clear(); // also, delete some bit tags, used to mark processed reds and blues mb->tag_delete(RedFlagTag);// to mark blue quads already considered mb->tag_delete(BlueFlagTag); /*mb->tag_delete(redParentTag); mb->tag_delete(blueParentTag); mb->tag_delete(countTag); RedEdges.clear(); localEnts.clear();*/ }
virtual int moab::Intx2Mesh::computeIntersectionBetweenRedAndBlue | ( | EntityHandle | red, |
EntityHandle | blue, | ||
double * | P, | ||
int & | nP, | ||
double & | area, | ||
int | markb[MAXEDGES], | ||
int | markr[MAXEDGES], | ||
int & | nsidesBlue, | ||
int & | nsidesRed, | ||
bool | check_boxes_first = false |
||
) | [pure virtual] |
Implemented in moab::Intx2MeshOnSphere, and moab::Intx2MeshInPlane.
void moab::Intx2Mesh::correct_polygon | ( | EntityHandle * | foundIds, |
int & | nP | ||
) |
Definition at line 452 of file Intx2Mesh.cpp.
{ int i = 0; while(i<nP) { int nextIndex = (i+1)%nP; if (nodes[i]==nodes[nextIndex]) { // we need to reduce nP, and collapse nodes if (dbg_1) { std::cout<<" nodes duplicated in list: " ; for (int j=0; j<nP; j++) std::cout<<nodes[j] << " " ; std::cout<<"\n"; std::cout<<" node " << nodes[i] << " at index " << i << " is duplicated" << "\n"; } // this will work even if we start from 1 2 3 1; when i is 3, we find nextIndex is 0, then next thing does nothing // (nP-1 is 3, so k is already >= nP-1); it will result in nodes -> 1, 2, 3 for (int k=i; k<nP-1; k++) nodes[k] = nodes[k+1]; nP--; // decrease the number of nodes; also, decrease i, just if we may need to check again i--; } i++; } return; }
ErrorCode moab::Intx2Mesh::create_departure_mesh_2nd_alg | ( | EntityHandle & | euler_set, |
EntityHandle & | covering_lagr_set | ||
) |
ErrorCode moab::Intx2Mesh::create_departure_mesh_3rd_alg | ( | EntityHandle & | lagr_set, |
EntityHandle & | covering_set | ||
) |
void moab::Intx2Mesh::createTags | ( | ) | [virtual] |
Definition at line 48 of file Intx2Mesh.cpp.
{ if (redParentTag) mb->tag_delete(redParentTag); if(blueParentTag) mb->tag_delete(blueParentTag); if (countTag) mb->tag_delete(countTag); /*RedEdges.clear(); localEnts.clear()*/ unsigned char def_data_bit = 0; // unused by default ErrorCode rval = mb->tag_get_handle("blueFlag", 1, MB_TYPE_BIT, BlueFlagTag, MB_TAG_CREAT, &def_data_bit); if (MB_SUCCESS != rval) return; // maybe the red tag is better to be deleted every time, and recreated; // or is it easy to set all values to something again? like 0? rval = mb->tag_get_handle("redFlag", 1, MB_TYPE_BIT, RedFlagTag, MB_TAG_CREAT, &def_data_bit); ERRORV(rval, "can't get red flag tag"); // assume that the edges are on the red triangles Range redElements; //Range redEdges; rval = mb->get_entities_by_dimension(mbs2, 2, redElements, false); // so all tri, quad, poly ERRORV(rval, "can't get ents by dimension"); // create red edges if they do not exist yet; so when they are looked upon, they are found // this is the only call that is potentially NlogN, in the whole method rval = mb->get_adjacencies(redElements, 1, true, RedEdges, Interface::UNION); ERRORV(rval, "can't get adjacent red edges"); // now, create a map from each edge to a list of potential new nodes on a red edge // this memory has to be cleaned up // change it to a vector, and use the index in range of red edges int indx = 0; extraNodesVec.reserve(RedEdges.size()); for (Range::iterator eit = RedEdges.begin(); eit != RedEdges.end(); eit++, indx++) { //EntityHandle edge = *eit; //extraNodesMap[edge] = new std::vector<EntityHandle>; std::vector<EntityHandle> * nv = new std::vector<EntityHandle>; extraNodesVec.push_back(nv); } int defaultInt = 0; rval = mb->tag_get_handle("RedParent", 1, MB_TYPE_INTEGER, redParentTag, MB_TAG_SPARSE | MB_TAG_CREAT, &defaultInt); ERRORV(rval, "can't create positive tag"); rval = mb->tag_get_handle("BlueParent", 1, MB_TYPE_INTEGER, blueParentTag, MB_TAG_SPARSE | MB_TAG_CREAT, &defaultInt); ERRORV(rval, "can't create negative tag"); rval = mb->tag_get_handle("Counting", 1, MB_TYPE_INTEGER, countTag, MB_TAG_SPARSE | MB_TAG_CREAT, &defaultInt); ERRORV(rval, "can't create Counting tag"); return; }
void moab::Intx2Mesh::disable_debug | ( | ) | [inline] |
Definition at line 102 of file Intx2Mesh.hpp.
{dbg_1 = 0;}
void moab::Intx2Mesh::enable_debug | ( | ) | [inline] |
Definition at line 101 of file Intx2Mesh.hpp.
{dbg_1 = 1;}
virtual int moab::Intx2Mesh::findNodes | ( | EntityHandle | red, |
int | nsRed, | ||
EntityHandle | blue, | ||
int | nsBlue, | ||
double * | iP, | ||
int | nP | ||
) | [pure virtual] |
Implemented in moab::Intx2MeshOnSphere, and moab::Intx2MeshInPlane.
ErrorCode moab::Intx2Mesh::GetOrderedNeighbors | ( | EntityHandle | set, |
EntityHandle | quad, | ||
EntityHandle | neighbors[MAXEDGES] | ||
) |
Definition at line 116 of file Intx2Mesh.cpp.
{ int nnodes = 3; // will get the nnodes ordered neighbors; // first cell is for nodes 0, 1, second to 1, 2, third to 2, 3, last to nnodes-1, const EntityHandle * conn4; ErrorCode rval = mb->get_connectivity(cell, conn4, nnodes); int nsides = nnodes; // account for possible padded polygons while (conn4[nsides-2]==conn4[nsides-1] && nsides>3) nsides--; ERRORR(rval, "can't get connectivity on an element"); for (int i = 0; i < nsides; i++) { EntityHandle v[2]; v[0] = conn4[i]; v[1] = conn4[(i + 1) % nsides]; // get quads adjacent to vertices std::vector<EntityHandle> cells; std::vector<EntityHandle> cellsInSet; rval = mb->get_adjacencies(v, 2, 2, false, cells, Interface::INTERSECT); ERRORR(rval, "can't get adjacencies on 2 nodes"); size_t siz = cells.size(); for (size_t j = 0; j < siz; j++) if (mb->contains_entities(set, &(cells[j]), 1)) cellsInSet.push_back(cells[j]); siz = cellsInSet.size(); if (siz > 2) { std::cout << "non manifold mesh, error" << mb->list_entities(&(cellsInSet[0]), cellsInSet.size()) << "\n"; return MB_FAILURE; // non-manifold } if (siz == 1) { // it must be the border, neighbors[i] = 0; // we are guaranteed that ids are !=0; this is marking a border // borders do not appear for a sphere in serial, but they do appear for // parallel processing anyway continue; } // here siz ==2, it is either the first or second if (cell == cellsInSet[0]) neighbors[i] = cellsInSet[1]; else neighbors[i] = cellsInSet[0]; } return MB_SUCCESS; }
ErrorCode moab::Intx2Mesh::intersect_meshes | ( | EntityHandle | mbs1, |
EntityHandle | mbs2, | ||
EntityHandle & | outputSet | ||
) |
Definition at line 169 of file Intx2Mesh.cpp.
{ ErrorCode rval; mbs1 = mbset1; // set 1 is departure, and it is completely covering the euler set on proc mbs2 = mbset2; outSet = outputSet; // really, should be something from t1 and t2; blue is 1 (lagrange), red is 2 (euler) createTags(); // EntityHandle startBlue, startRed; mb->get_entities_by_dimension(mbs1, 2, rs1); mb->get_entities_by_dimension(mbs2, 2, rs2); Range rs22=rs2; // a copy of the initial range; we will remove from it elements as we // advance ; rs2 is needed for marking the polygon to the red parent while (!rs22.empty()) { for (Range::iterator it = rs1.begin(); it != rs1.end(); it++) { startBlue = *it; int found = 0; for (Range::iterator it2 = rs22.begin(); it2 != rs22.end() && !found; it2++) { startRed = *it2; double area = 0; // if area is > 0 , we have intersections double P[10*MAXEDGES]; // max 8 intx points + 8 more in the polygon // int nP = 0; int nb[MAXEDGES], nr[MAXEDGES]; // sides 3 or 4? also, check boxes first int nsRed, nsBlue; computeIntersectionBetweenRedAndBlue(startRed, startBlue, P, nP, area, nb, nr, nsBlue, nsRed, true); if (area > 0) { found = 1; break; // found 2 elements that intersect; these will be the seeds } } if (found) break; } std::queue<EntityHandle> blueQueue; // these are corresponding to Ta, blueQueue.push(startBlue); std::queue<EntityHandle> redQueue; redQueue.push(startRed); Range toResetBlues; // will be used to reset blue flags for every red quad // processed /*if (my_rank==0) dbg_1 = 1;*/ unsigned char used = 1; unsigned char unused = 0; // for red flags // mark the start blue quad as used, so it will not come back again mb->tag_set_data(RedFlagTag, &startRed, 1, &used); while (!redQueue.empty()) { // flags for the side : 0 means a blue quad not found on side // a paired blue not found yet for the neighbors of red EntityHandle n[MAXEDGES] = { EntityHandle(0) }; int nsidesRed; // will be initialized later, when we compute intersection EntityHandle currentRed = redQueue.front(); redQueue.pop(); // for (k=0; k<m_numPos; k++) // redFlag[k] = 0; // redFlag[m_numPos] = 1; // to guard for the boundary // all reds that were tagged, are now cleared for (Range::iterator itr = toResetBlues.begin(); itr != toResetBlues.end(); itr++) { EntityHandle ttt = *itr; rval = mb->tag_set_data(BlueFlagTag, &ttt, 1, &unused); ERRORR(rval, "can't set blue unused tag"); } //rval = mb2->tag_set_data(RedFlagTag, toResetReds, &unused); /*if (dbg_1) { std::cout << "reset blues: "; mb->list_entities(toResetBlues); }*/ //rval = mb2->tag_set_data(RedFlagTag, toResetReds, &unused); if (dbg_1) { std::cout << "reset blues: "; for (Range::iterator itr = toResetBlues.begin(); itr != toResetBlues.end(); itr++) std::cout << mb->id_from_handle(*itr) << " "; std::cout << std::endl; } EntityHandle currentBlue = blueQueue.front(); // where do we check for redQueue???? // red and blue queues are parallel blueQueue.pop(); // mark the current red //redFlag[currentRed] = 1; // toResetBlues.clear(); // empty the range of used blues, will have to be set unused again, // at the end of red element processing toResetBlues.insert(currentBlue); rval = mb->tag_set_data(BlueFlagTag, ¤tBlue, 1, &used); ERRORR(rval, "can't set blue tag"); //mb2->set_tag_data std::queue<EntityHandle> localBlue; localBlue.push(currentBlue); while (!localBlue.empty()) { // EntityHandle blueT = localBlue.front(); localBlue.pop(); double P[10*MAXEDGES], area; // int nP = 0; int nb[MAXEDGES] = {0}; int nr[MAXEDGES] = {0}; int nsidesBlue; // area is in 2d, points are in 3d (on a sphere), back-projected, or in a plane // intersection points could include the vertices of initial elements // nb [j] = 0 means no intersection on the side k for element blue (markers) // nb [j] = 1 means that the side j (from j to j+1) of blue poly intersects the // red poly. A potential next poly is the red poly that is adjacent to this side computeIntersectionBetweenRedAndBlue(/* red */currentRed, blueT, P, nP, area, nb, nr, nsidesBlue, nsidesRed); if (nP > 0) { if (dbg_1) { for (int k=0; k<3; k++) { std::cout << " nb, nr: " << k << " " << nb[k] << " " << nr[k] << "\n"; } } // intersection found: output P and original triangles if nP > 2 EntityHandle neighbors[MAXEDGES]; rval = GetOrderedNeighbors(mbs1, blueT, neighbors); if (rval != MB_SUCCESS) { std::cout << " can't get the neighbors for blue element " << mb->id_from_handle(blueT); return MB_FAILURE; } // add neighbors to the localBlue queue, if they are not marked for (int nn = 0; nn < nsidesBlue; nn++) { EntityHandle neighbor = neighbors[nn]; if (neighbor > 0 && nb[nn]>0) // advance across blue boundary n { //n[nn] = redT; // start from 0!! unsigned char status = 0; mb->tag_get_data(BlueFlagTag, &neighbor, 1, &status); if (status == 0) { localBlue.push(neighbor); /*if (dbg_1) { std::cout << " local blue elem " << mb->list_entities(&neighbor, 1) << "\n for red:" << mb->list_entities(¤tRed, 1) << "\n"; }*/ if (dbg_1) { std::cout << " local blue elem " << mb->id_from_handle(neighbor) << " for red:" << mb->id_from_handle(currentRed) << "\n"; mb->list_entities(&neighbor, 1); } rval = mb->tag_set_data(BlueFlagTag, &neighbor, 1, &used); //redFlag[neighbor] = 1; // flag it to not be added anymore toResetBlues.insert(neighbor); // this is used to reset the red flag } } // n(find(nc>0))=ac; % ac is starting candidate for neighbor if (nr[nn] > 0) n[nn] = blueT; } if (nP > 1) // this will also construct triangles/polygons in the new mesh, if needed findNodes(currentRed, nsidesRed, blueT, nsidesBlue, P, nP); } else if (dbg_1) /*{ std::cout << " red, blue, do not intersect: " << mb->list_entities(¤tRed, 1) << " " << mb->list_entities(&blueT, 1) << "\n"; }*/ { std::cout << " red, blue, do not intersect: " << mb->id_from_handle(currentRed) << " " << mb->id_from_handle(blueT) << "\n"; } } // end while (!localBlue.empty()) // here, we are finished with redCurrent, take it out of the rs22 range (red, arrival mesh) rs22.erase(currentRed); // also, look at its neighbors, and add to the seeds a next one EntityHandle redNeighbors[MAXEDGES]; rval = GetOrderedNeighbors(mbs2, currentRed, redNeighbors); ERRORR(rval, "can't get neighbors"); if (dbg_1) { std::cout << "Next: neighbors for current red "; for (int kk = 0; kk < nsidesRed; kk++) { if (redNeighbors[kk] > 0) std::cout << mb->id_from_handle(redNeighbors[kk]) << " "; else std::cout << 0 << " "; } std::cout << std::endl; } for (int j = 0; j < nsidesRed; j++) { EntityHandle redNeigh = redNeighbors[j]; unsigned char status = 1; if (redNeigh == 0) continue; mb->tag_get_data(RedFlagTag, &redNeigh, 1, &status); // status 0 is unused if (status == 0 && n[j] > 0) // not treated yet and marked as a neighbor { // we identified red quad n[j] as intersecting with neighbor j of the blue quad redQueue.push(redNeigh); blueQueue.push(n[j]); if (dbg_1) std::cout << "new polys pushed: blue, red:" << mb->id_from_handle(redNeigh) << " " << mb->id_from_handle(n[j]) << std::endl; mb->tag_set_data(RedFlagTag, &redNeigh, 1, &used); } } } // end while (!redQueue.empty()) } if (dbg_1) { for (int k = 0; k < 6; k++) mout_1[k].close(); } // before cleaning up , we need to settle the position of the intersection points // on the boundary edges // this needs to be collective, so we should maybe wait something #ifdef USE_MPI rval = correct_intersection_points_positions(); if (rval!=MB_SUCCESS) { std::cout << "can't correct position, Intx2Mesh.cpp \n"; } #endif clean(); return MB_SUCCESS; }
virtual bool moab::Intx2Mesh::is_inside_element | ( | double | xyz[3], |
EntityHandle | eh | ||
) | [pure virtual] |
Implemented in moab::Intx2MeshOnSphere, and moab::Intx2MeshInPlane.
void moab::Intx2Mesh::set_box_error | ( | double | berror | ) | [inline] |
Definition at line 81 of file Intx2Mesh.hpp.
{box_error = berror;}
void moab::Intx2Mesh::SetErrorTolerance | ( | double | eps | ) | [inline] |
Definition at line 72 of file Intx2Mesh.hpp.
{ epsilon_1=eps; epsilon_area = eps*eps/2;}
std::vector<double> moab::Intx2Mesh::allBoxes [protected] |
Definition at line 146 of file Intx2Mesh.hpp.
const EntityHandle* moab::Intx2Mesh::blueConn [protected] |
Definition at line 129 of file Intx2Mesh.hpp.
CartVect moab::Intx2Mesh::blueCoords[MAXEDGES] [protected] |
Definition at line 131 of file Intx2Mesh.hpp.
double moab::Intx2Mesh::blueCoords2D[MAXEDGES2] [protected] |
Definition at line 133 of file Intx2Mesh.hpp.
Tag moab::Intx2Mesh::BlueFlagTag [protected] |
Definition at line 114 of file Intx2Mesh.hpp.
Tag moab::Intx2Mesh::blueParentTag [protected] |
Definition at line 123 of file Intx2Mesh.hpp.
double moab::Intx2Mesh::box_error [protected] |
Definition at line 147 of file Intx2Mesh.hpp.
Tag moab::Intx2Mesh::countTag [protected] |
Definition at line 124 of file Intx2Mesh.hpp.
int moab::Intx2Mesh::dbg_1 [protected] |
Definition at line 136 of file Intx2Mesh.hpp.
double moab::Intx2Mesh::epsilon_1 [protected] |
Definition at line 143 of file Intx2Mesh.hpp.
double moab::Intx2Mesh::epsilon_area [protected] |
Definition at line 144 of file Intx2Mesh.hpp.
std::vector<std::vector<EntityHandle> *> moab::Intx2Mesh::extraNodesVec [protected] |
Definition at line 141 of file Intx2Mesh.hpp.
Range moab::Intx2Mesh::localEnts [protected] |
Definition at line 150 of file Intx2Mesh.hpp.
EntityHandle moab::Intx2Mesh::localRoot [protected] |
Definition at line 149 of file Intx2Mesh.hpp.
int moab::Intx2Mesh::max_edges [protected] |
Definition at line 158 of file Intx2Mesh.hpp.
Interface* moab::Intx2Mesh::mb [protected] |
Definition at line 104 of file Intx2Mesh.hpp.
EntityHandle moab::Intx2Mesh::mbs1 [protected] |
Definition at line 106 of file Intx2Mesh.hpp.
EntityHandle moab::Intx2Mesh::mbs2 [protected] |
Definition at line 107 of file Intx2Mesh.hpp.
std::ofstream moab::Intx2Mesh::mout_1[6] [protected] |
Definition at line 135 of file Intx2Mesh.hpp.
unsigned int moab::Intx2Mesh::my_rank [protected] |
Definition at line 151 of file Intx2Mesh.hpp.
EntityHandle moab::Intx2Mesh::outSet [protected] |
Definition at line 111 of file Intx2Mesh.hpp.
const EntityHandle* moab::Intx2Mesh::redConn [protected] |
Definition at line 128 of file Intx2Mesh.hpp.
CartVect moab::Intx2Mesh::redCoords[MAXEDGES] [protected] |
Definition at line 130 of file Intx2Mesh.hpp.
double moab::Intx2Mesh::redCoords2D[MAXEDGES2] [protected] |
Definition at line 132 of file Intx2Mesh.hpp.
Range moab::Intx2Mesh::RedEdges [protected] |
Definition at line 118 of file Intx2Mesh.hpp.
Tag moab::Intx2Mesh::RedFlagTag [protected] |
Definition at line 116 of file Intx2Mesh.hpp.
Tag moab::Intx2Mesh::redParentTag [protected] |
Definition at line 122 of file Intx2Mesh.hpp.
Range moab::Intx2Mesh::rs1 [protected] |
Definition at line 108 of file Intx2Mesh.hpp.
Range moab::Intx2Mesh::rs2 [protected] |
Definition at line 109 of file Intx2Mesh.hpp.