43 : m_trace(trace), m_interface(interfaceShPtr)
55 : m_graph(meshGraph), m_movement(meshGraph->GetMovement()), m_trace(trace)
57 auto comm =
m_trace->GetComm()->GetSpaceComm();
58 auto interfaceCollection =
m_movement->GetInterfaces();
64 std::map<int, int> myIndxLRMap;
69 std::map<int, std::pair<InterfaceTraceSharedPtr, InterfaceTraceSharedPtr>>
79 for (
const auto &interface : interfaceCollection)
81 indxToInterfaceID[cnt] = interface.first.first;
82 myIndxLRMap[interface.first.first] = 0;
84 if (!interface.second->GetLeftInterface()->IsEmpty())
86 myIndxLRMap[interface.first.first] += 1;
88 localInterfaces[interface.first.first].first =
90 trace, interface.second->GetLeftInterface());
92 localInterfaces[interface.first.first].first);
95 if (!interface.second->GetRightInterface()->IsEmpty())
97 myIndxLRMap[interface.first.first] += 2;
99 localInterfaces[interface.first.first].second =
101 trace, interface.second->GetRightInterface());
103 localInterfaces[interface.first.first].second);
110 int nRanks = comm->GetSize();
115 for (
auto pres : myIndxLRMap)
117 interfaceEdges[cnt++] = pres.second;
120 comm->AllGather(interfaceEdges, rankLocalInterfaceIds);
123 std::map<int, std::vector<InterfaceTraceSharedPtr>> oppRankSharedInterface;
127 size_t myRank = comm->GetRank();
128 size_t numInterfaces = interfaceCollection.size();
129 for (
int i = 0; i < nRanks; ++i)
131 for (
size_t j = 0; j < numInterfaces; ++j)
133 int otherId = indxToInterfaceID[j];
137 int otherCode = rankLocalInterfaceIds[i * numInterfaces + j];
141 int myCode = myIndxLRMap[otherId];
149 localInterfaces[otherId].first->SetCheckLocal(
true);
150 localInterfaces[otherId].second->SetCheckLocal(
true);
157 if ((myCode == 1 && otherCode == 2) ||
158 (myCode == 1 && otherCode == 3) ||
159 (myCode == 3 && otherCode == 2))
161 oppRankSharedInterface[i].emplace_back(
162 localInterfaces[otherId].first);
165 else if ((myCode == 2 && otherCode == 1) ||
166 (myCode == 2 && otherCode == 3) ||
167 (myCode == 3 && otherCode == 1))
169 oppRankSharedInterface[i].emplace_back(
170 localInterfaces[otherId].second);
173 else if (myCode == 3 && otherCode == 3)
175 oppRankSharedInterface[i].emplace_back(
176 localInterfaces[otherId].first);
177 oppRankSharedInterface[i].emplace_back(
178 localInterfaces[otherId].second);
185 for (
auto &rank : oppRankSharedInterface)
198 auto comm =
m_trace->GetComm()->GetSpaceComm();
203 interfaceTrace->CalcLocalMissing();
209 auto requestSend = comm->CreateRequest(
m_exchange.size());
210 auto requestRecv = comm->CreateRequest(
m_exchange.size());
214 m_exchange[i]->RankFillSizes(requestSend, requestRecv, i);
216 comm->WaitAll(requestSend);
217 comm->WaitAll(requestRecv);
221 m_exchange[i]->SendMissing(requestSend, requestRecv, i);
223 comm->WaitAll(requestSend);
224 comm->WaitAll(requestRecv);
228 i->CalcRankDistances();
250 for (
auto &childId : childEdge)
252 auto childElmt =
m_trace->GetExpFromGeomId(childId.first);
253 size_t nq = childElmt->GetTotPoints();
255 childElmt->GetCoords(xc, yc, zc);
259 for (
int i = 0; i < nq; ++i, ++cnt)
283 for (
auto &childId : childEdge)
285 auto childElmt =
m_trace->GetExpFromGeomId(childId.first);
286 size_t nq = childElmt->GetTotPoints();
288 childElmt->GetCoords(xc, yc, zc);
292 for (
int i = 0; i < nq; ++i)
301 auto parentEdge =
m_interface->GetOppInterface()->GetEdge();
302 for (
auto &edge : parentEdge)
305 if (!edge.second->MinMaxCheck(xs))
311 edge.second->FindDistance(xs, foundLocCoord);
316 edge.second->GetGlobalID(), foundLocCoord);
331 if (
m_trace->GetComm()->IsSerial())
335 " coordinates on interface ID " +
337 " linked to interface ID " +
338 std::to_string(
m_interface->GetOppInterface()->GetId()) +
339 ". Check both sides of the interface line up.");
358 m_zones[localInterface->GetInterface()->GetId()]->GetMoved();
367 if (
m_zones[interface->GetInterface()->GetId()]->GetMoved())
373 interface->GetMissingCoords().size() * 3;
379 requestSend, requestNum);
381 requestRecv, requestNum);
401 if (
m_zones[interface->GetInterface()->GetId()]->GetMoved())
403 auto missing = interface->GetMissingCoords();
404 for (
auto coord : missing)
408 for (
int k = 0; k < 3; ++k, ++cnt)
433 auto comm =
m_trace->GetComm()->GetSpaceComm();
441 m_localInterface->FillLocalBwdTrace(Fwd, Bwd);
446 auto requestSend = comm->CreateRequest(
m_exchange.size());
447 auto requestRecv = comm->CreateRequest(
m_exchange.size());
450 m_exchange[i]->SendFwdTrace(requestSend, requestRecv, i, Fwd);
456 m_localInterface->FillLocalBwdTrace(Fwd, Bwd);
459 comm->WaitAll(requestSend);
460 comm->WaitAll(requestRecv);
465 i->FillRankBwdTraceExchange(Bwd);
478 int traceId =
m_trace->GetElmtToExpId(foundLocCoord.second.first);
482 Fwd +
m_trace->GetPhys_Offset(traceId);
484 Bwd[foundLocCoord.first] =
485 m_trace->GetExp(traceId)->StdPhysEvaluate(locCoord, edgePhys);
510 if (!std::isnan(trace[i]))
529 int traceId =
m_trace->GetElmtToExpId(i.second.first);
533 Fwd +
m_trace->GetPhys_Offset(traceId);
536 m_trace->GetExp(traceId)->StdPhysEvaluate(locCoord, edgePhys);
570 for (
int j = disp[i]; j < disp[i + 1]; j += 3)
578 for (
auto &edge : localEdge)
581 if (!edge.second->MinMaxCheck(xs))
586 NekDouble dist = edge.second->FindDistance(xs, foundLocCoord);
591 edge.second->GetGlobalID(), foundLocCoord);
#define ASSERTL0(condition, msg)
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::map< int, Array< OneD, int > > m_sendSize
Map of rank to array of size of send buffer for each interface.
void RankFillSizes(LibUtilities::CommRequestSharedPtr &requestSend, LibUtilities::CommRequestSharedPtr &requestRecv, int requestNum)
std::map< int, SpatialDomains::ZoneBaseShPtr > m_zones
Map of zone IDs to zone bases.
const ExpListSharedPtr m_trace
Trace expansion list.
void FillRankBwdTraceExchange(Array< OneD, NekDouble > &Bwd)
const std::vector< InterfaceTraceSharedPtr > m_interfaceTraces
Vector of interface traces i.e. every interface side present on m_rank.
const LibUtilities::CommSharedPtr m_comm
Communicator.
void SendFwdTrace(LibUtilities::CommRequestSharedPtr &requestSend, LibUtilities::CommRequestSharedPtr &requestRecv, int requestNum, Array< OneD, NekDouble > &Fwd)
Array< OneD, NekDouble > m_send
Send buffer for coord exchange.
std::map< int, std::map< int, std::pair< int, Array< OneD, NekDouble > > > > m_foundRankCoords
std::map< int, int > m_totSendSize
Map of rank to total size of send buffer for all interfaces.
void CalcRankDistances()
Populates m_foundRankCoords using the FindDistance function.
std::map< int, int > m_totRecvSize
Map of rank to total size of receive buffer for all interfaces.
std::map< int, Array< OneD, int > > m_recvSize
Map of rank to array of size of receive buffer for each interface.
Array< OneD, NekDouble > m_recv
Receive buffer for coord exchange.
void SendMissing(LibUtilities::CommRequestSharedPtr &requestSend, LibUtilities::CommRequestSharedPtr &requestRecv, int requestNum)
Array< OneD, NekDouble > m_recvTrace
Receive buffer for trace exchange.
Array< OneD, NekDouble > m_sendTrace
Send buffer for trace exchange.
const ExpListSharedPtr m_trace
Trace expansion list.
void ExchangeTrace(Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
Perform the trace exchange between processors, given the forwards and backwards spaces.
std::vector< InterfaceExchangeSharedPtr > m_exchange
Vector of interface exchanges, i.e. every rank-to-rank comm needed.
void ExchangeCoords()
Perform the coordinate exchange between processors. This is where the missing coordinates on the inte...
SpatialDomains::MovementSharedPtr m_movement
Movement object associated with the non-conformal interfaces.
InterfaceMapDG(const SpatialDomains::MeshGraphSharedPtr &graph, const ExpListSharedPtr &trace)
std::vector< InterfaceTraceSharedPtr > m_localInterfaces
Interface sides present on current process.
ExpListSharedPtr m_trace
Trace expansion list.
std::map< int, std::pair< int, Array< OneD, NekDouble > > > m_foundLocalCoords
Map of found coordinates present locally.
std::vector< int > m_mapMissingCoordToTrace
Vector of indices corresponding to m_missingCoord locations in trace.
SpatialDomains::InterfaceShPtr m_interface
Local interface object.
std::vector< Array< OneD, NekDouble > > m_missingCoords
Vector of coordinates on interface missing from the other side locally.
InterfaceTrace(const ExpListSharedPtr &trace, const SpatialDomains::InterfaceShPtr &interfaceShPtr)
Constructor.
void FillRankBwdTrace(Array< OneD, NekDouble > &trace, Array< OneD, NekDouble > &Bwd)
Fills the Bwd trace from partitioned trace.
void CalcLocalMissing()
Calculates what coordinates on the interface are missing locally.
void FillLocalBwdTrace(Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
Fills the Bwd trace by interpolating from the Fwd for local interfaces.
bool m_checkLocal
Flag whether the opposite side of the interface is present locally.
std::shared_ptr< CommRequest > CommRequestSharedPtr
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
static const NekDouble kFindDistanceMin
std::shared_ptr< Interface > InterfaceShPtr
std::shared_ptr< MeshGraph > MeshGraphSharedPtr