41 template <
class T,
unsigned Dim,
class Mesh >
51 template <
class T,
unsigned Dim,
class Mesh >
61 template <
class T,
unsigned Dim,
class Mesh >
71 template <
class T,
unsigned Dim,
class Mesh >
80 template <
class T,
unsigned Dim,
class Mesh >
85 InterNodeList =
new bool[N];
86 SentToNodeList =
new bool[N];
92 MaxGlobalInterRadius = 0;
98 template <
class T,
unsigned Dim,
class Mesh >
101 delete [] InterNodeList;
102 delete [] SentToNodeList;
104 for (
int i=(PairList.size() - 1); i >= 0; --i)
105 delete (PairList[i]);
111 template <
class T,
unsigned Dim,
class Mesh >
114 if (InterRadiusArray != 0) {
115 if (InterRadiusArray->size() > 0)
116 return *(max_element(InterRadiusArray->begin(),
117 InterRadiusArray->end()));
132 template <
class T,
unsigned Dim,
class Mesh >
141 bool periodicBC[2*
Dim];
142 unsigned numPeriodicBC = 0;
144 for (
unsigned d=0; d<2*
Dim; ++d) {
145 periodicBC[d] = (pBConds[d] == periodicBCond);
146 if (periodicBC[d]) ++numPeriodicBC;
148 if (numPeriodicBC>0) {
151 swap_ghost_particles(PData.getLocalNum(), PData, periodicBC);
154 swap_ghost_particles(PData.getLocalNum(), PData);
160 swap_ghost_particles(PData.getLocalNum(), PData);
164 bpi = PairList[
n]->begin();
165 epi = PairList[
n]->end();
175 template <
class T,
unsigned Dim,
class Mesh >
180 DEBUGMSG(
"ParticleInteractLayout: rebuilding interaction node data."<<
endl);
183 InteractionNodes = 0;
184 T interRad = 2.0 * getMaxInteractionRadius();
188 for (j=0; j < N; ++j)
189 InterNodeList[j] =
false;
198 for (localVN = this->RLayout.
begin_iv(); localVN != endLocalVN; ++localVN) {
202 for (d=0; d <
Dim; ++d) {
203 chkDom[d] =
PRegion<T>(chkDom[d].first() - interRad,
204 chkDom[d].last() + interRad);
214 for ( ; tVN != touchingVN.second; ++tVN) {
217 DEBUGMSG(
" vnode: node = " << ((*tVN).second)->getNode());
218 DEBUGMSG(
", domain = " << ((*tVN).second)->getDomain() <<
endl);
219 unsigned vn = ((*tVN).second)->getNode();
220 if ( ! InterNodeList[vn] ) {
221 InterNodeList[vn] =
true;
227 DEBUGMSG(
" There are " << InteractionNodes <<
" inter. nodes." <<
endl);
232 NeedGhostSwap =
true;
242 template <
class T,
unsigned Dim,
class Mesh >
244 const bool periodicBC[2*
Dim])
249 DEBUGMSG(
"ParticleInteractLayout: rebuilding interaction node data."<<
endl);
252 InteractionNodes = 0;
253 T interRad = 2.0 * getMaxInteractionRadius();
257 for (j=0; j < N; ++j)
258 InterNodeList[j] =
false;
270 bool flipBit, activeBit[
Dim], refBit[
Dim];
280 for (localVN = this->RLayout.
begin_iv(); localVN != endLocalVN; ++localVN) {
283 chkDom = (*localVN).second->getDomain();
284 for (d=0; d<
Dim; ++d) {
285 chkDom[d] =
PRegion<T>(chkDom[d].first() - interRad,
286 chkDom[d].last() + interRad);
294 for (tVN = touchingVN.first; tVN != touchingVN.second; ++tVN) {
297 DEBUGMSG(
" vnode: node = " << ((*tVN).second)->getNode());
298 DEBUGMSG(
", domain = " << ((*tVN).second)->getDomain() <<
endl);
299 unsigned vn = ((*tVN).second)->getNode();
300 if ( ! InterNodeList[vn] ) {
301 InterNodeList[vn] =
true;
308 for (d=0; d<
Dim; ++d) {
309 if (periodicBC[2*d] && chkDom[d].first()<globalDom[d].first()) {
311 offset[d] = globalDom[d].length();
313 numRef = 2 * numRef + 1;
315 else if (periodicBC[2*d+1] && chkDom[d].last()>globalDom[d].last()) {
317 offset[d] = -globalDom[d].length();
319 numRef = 2 * numRef + 1;
323 activeBit[d] =
false;
329 for (j=0; j<numRef; ++j) {
335 while (d<Dim && !flipBit) {
354 for (d=0; d<
Dim; ++d) {
355 if (refBit[d]) refDom[d] = refDom[d] + offset[d];
363 for (tVN = touchingVN.first; tVN != touchingVN.second; ++tVN) {
366 DEBUGMSG(
" vnode: node = " << ((*tVN).second)->getNode());
367 DEBUGMSG(
", domain = " << ((*tVN).second)->getDomain() <<
endl);
368 unsigned vn = ((*tVN).second)->getNode();
369 if ( ! InterNodeList[vn] ) {
370 InterNodeList[vn] =
true;
375 if (!InterNodeList[pe]) {
377 bool interact =
false;
378 localVN2 = this->RLayout.
begin_iv();
379 while (localVN2 != endLocalVN && !interact) {
380 interact = refDom.
touches((*localVN2).second->getDomain());
384 InterNodeList[pe] =
true;
392 DEBUGMSG(
" There are " << InteractionNodes <<
" inter. nodes." <<
endl);
397 NeedGhostSwap =
true;
407 template <
class T,
unsigned Dim,
class Mesh >
414 unsigned LocalNum = PData.getLocalNum();
415 unsigned DestroyNum = PData.getDestroyNum();
417 T maxrad = getMaxLocalInteractionRadius();
421 PData.performDestroy();
422 LocalNum -= DestroyNum;
426 if ( ! this->RLayout.initialized())
427 rebuild_layout(LocalNum,PData);
431 apply_bconds(LocalNum, PData.R, this->getBConds(), this->RLayout.getDomain());
436 if (N > 1 && getUpdateFlag(this->
SWAP)) {
438 LocalNum = swap_particles(LocalNum, PData);
440 LocalNum = swap_particles(LocalNum, PData, *canSwap);
444 NeedGhostSwap =
true;
447 TotalNum = this->NodeCount[myN] = LocalNum;
472 msg->
get(this->NodeCount);
477 int notrecvd = N - 1;
479 while (notrecvd > 0) {
483 int remNodeCount = 0;
485 msg->
get(remNodeCount);
491 TotalNum += remNodeCount;
492 this->NodeCount[node] = remNodeCount;
493 if (remMaxRad > maxrad)
499 msg->
put(this->NodeCount, this->NodeCount + N);
507 PData.setTotalNum(TotalNum);
508 PData.setLocalNum(LocalNum);
511 if (maxrad != getMaxInteractionRadius()) {
512 setMaxInteractionRadius(maxrad);
517 bool periodicBC[2*
Dim];
518 unsigned numPeriodicBC = 0;
520 for (
unsigned d=0; d<2*
Dim; ++d) {
521 periodicBC[d] = (pBConds[d] == periodicBCond);
522 if (periodicBC[d]) ++numPeriodicBC;
524 if (numPeriodicBC>0) {
527 rebuild_interaction_data(periodicBC);
530 rebuild_interaction_data();
534 rebuild_interaction_data();
547 template <
class T,
unsigned Dim,
class Mesh >
555 if ( ! NeedGhostSwap )
return;
558 NeedGhostSwap =
false;
563 PData.ghostDestroy(PData.getGhostNum(), 0);
567 unsigned sendnum = 0;
568 for (i=0; i < N; i++)
569 if (InterNodeList[i] && this->NodeCount[i] > 0)
574 if (sendnum == 0 || LocalNum == 0) {
575 find_pairs(LocalNum, 0, LocalNum,
true, PData);
581 T interRad = 2.0 * getMaxInteractionRadius();
591 for (i=0; i < N; i++)
592 if (InterNodeList[i] && this->NodeCount[i] > 0)
593 this->SwapMsgList[i] =
new Message;
653 find_pairs(LocalNum, 0, LocalNum,
true, PData);
682 template <
class T,
unsigned Dim,
class Mesh >
686 const bool periodicBC[2*
Dim])
692 if ( ! NeedGhostSwap )
return;
695 NeedGhostSwap =
false;
700 PData.ghostDestroy(PData.getGhostNum(), 0);
705 unsigned sendnum = 0;
706 for (i=0; i < N; i++)
707 if (InterNodeList[i] && this->NodeCount[i] > 0)
712 if (sendnum == 0 || LocalNum == 0) {
713 find_pairs(LocalNum, 0, LocalNum,
true, PData);
719 T interRad = 2.0 * getMaxInteractionRadius();
727 bool flipBit, activeBit[
Dim], refBit[
Dim];
741 for (i=0; i < N; i++)
742 if (InterNodeList[i] && this->NodeCount[i] > 0)
743 this->SwapMsgList[i] =
new Message;
749 for (i=0; i < LocalNum; ++i) {
752 memset((
void *)SentToNodeList, 0, N *
sizeof(
bool));
756 for (j=0; j < (
unsigned int) Dim; ++j)
757 pLoc[j] =
PRegion<T>(PData.R[i][j] - interRad,
758 PData.R[i][j] + interRad);
762 for (tVNit = touchingVN.first; tVNit != touchingVN.second; ++tVNit) {
764 unsigned node = tVN->
getNode();
767 if (this->NodeCount[node] > 0 && ! SentToNodeList[node]) {
768 if (! InterNodeList[node]) {
769 ERRORMSG(
"ParticleInteractLayout: Cannot send ghost " << i);
770 ERRORMSG(
" to node " << node <<
" ... skipping." <<
endl);
773 PData.ghostPutMessage(*(this->SwapMsgList[node]), 1, i);
774 SentToNodeList[node] =
true;
781 for (d=0; d<
Dim; ++d) {
782 if (periodicBC[2*d] && pLoc[d].first()<globalDom[d].first()) {
784 offset[d] = globalDom[d].length();
786 numRef = 2 * numRef + 1;
788 else if (periodicBC[2*d+1] && pLoc[d].last()>globalDom[d].last()) {
790 offset[d] = -globalDom[d].length();
792 numRef = 2 * numRef + 1;
796 activeBit[d] =
false;
801 if (numRef>0) savePos = PData.R[i];
804 for (j=0; j<numRef; ++j) {
807 PData.R[i] = savePos;
811 while (d<Dim && !flipBit) {
830 for (d=0; d<
Dim; ++d) {
832 refLoc[d] = refLoc[d] + offset[d];
833 PData.R[i][d] = PData.R[i][d] + offset[d];
838 memset((
void *)SentToNodeList, 0, N *
sizeof(
bool));
842 for (tVNit = touchingVN.first; tVNit != touchingVN.second; ++tVNit) {
844 unsigned node = tVN->
getNode();
847 if (this->NodeCount[node] > 0 && ! SentToNodeList[node]) {
848 if (! InterNodeList[node]) {
849 ERRORMSG(
"ParticleInteractLayout: Cannot send ghost " << i);
850 ERRORMSG(
" to node " << node <<
" ... skipping." <<
endl);
853 PData.ghostPutMessage(*(this->SwapMsgList[node]), 1, i);
854 SentToNodeList[node] =
true;
859 if (InterNodeList[pe]) {
861 bool interact =
false;
863 while (localVN != endLocalVN && !interact) {
864 interact = refLoc.
touches((*localVN).second->getDomain());
868 PData.ghostPutMessage(*(this->SwapMsgList[pe]), 1, i);
872 if (numRef>0) PData.R[i] = savePos;
878 for (i=0; i < N; ++i) {
879 if (InterNodeList[i] && this->NodeCount[i] > 0) {
881 PData.ghostPutMessage(*(this->SwapMsgList[i]), (
unsigned)0, (
unsigned)0);
889 find_pairs(LocalNum, 0, LocalNum,
true, PData);
892 while (sendnum-- > 0) {
894 unsigned oldGN = PData.getGhostNum();
897 while (PData.ghostGetMessage(*recmsg, node) > 0);
901 find_pairs(LocalNum, LocalNum + oldGN, LocalNum + PData.getGhostNum(),
910 template <
class T,
unsigned Dim,
class Mesh >
912 const unsigned a1,
const unsigned a2,
const bool initLists,
919 unsigned vlen = PairList.size();
922 for (i=0; i < vlen; ++i)
923 PairList[i]->erase(PairList[i]->begin(), PairList[i]->end());
926 if (PairList.size() < LocalNum) {
927 int newamt = LocalNum - PairList.size();
928 PairList.reserve(newamt);
929 for (
int k=0; k < newamt; ++k)
930 PairList.push_back(
new std::vector<pair_t>);
935 if (a2 <= a1)
return;
938 for (i=0; i < LocalNum; ++i) {
940 T intrad1 = getInteractionRadius(i);
943 j = (a1 > i ? a1 : i + 1);
947 for (; j < LocalNum; ++j) {
949 T intrad2 = intrad1 + getInteractionRadius(j);
955 T sep2 =
dot(rsep, rsep);
959 if (sep2 < intrad2) {
960 PairList[i]->push_back(
pair_t(j, sep2));
961 PairList[j]->push_back(
pair_t(i, sep2));
967 for (; j < a2; ++j) {
969 T intrad2 = intrad1 + getInteractionRadius(j);
975 T sep2 =
dot(rsep, rsep);
980 if (sep2 < intrad2) {
981 PairList[i]->push_back(
pair_t(j, sep2));
990 template <
class T,
unsigned Dim,
class Mesh >
991 std::ostream& operator<<(std::ostream& out, const ParticleInteractLayout<T,Dim,Mesh>& L) {
993 out <<
"ParticleInteractLayout, with particle distribution:\n ";
995 out << L.getNodeCount(i) <<
" ";
996 out <<
"\nInteractLayout decomposition = " << L.getLayout();
1005 template <
class T,
unsigned Dim,
class Mesh >
1012 this->rebuild_neighbor_data();
1017 InteractionNodes = 0;
1018 setMaxInteractionRadius(0);
1019 NeedGhostSwap =
true;
void swap_ghost_particles(unsigned, IpplParticleBase< ParticleInteractLayout< T, Dim, Mesh > > &)
T ParticlePeriodicBCond(const T t, const T minval, const T maxval)
ID_t getUserListID() const
#define P_SPATIAL_LAYOUT_TAG
virtual void Repartition(UserList *)
double dot(const Vector3D &lhs, const Vector3D &rhs)
Vector dot product.
int next_tag(int t, int s=1000)
void update(IpplParticleBase< ParticleInteractLayout< T, Dim, Mesh > > &p, const ParticleAttrib< char > *canSwap=0)
#define P_SPATIAL_RETURN_TAG
bool touches(const NDRegion< T, Dim > &nr) const
std::vector< pair_t >::iterator pair_iterator
#define P_SPATIAL_GHOST_TAG
void getPairlist(unsigned, pair_iterator &, pair_iterator &, IpplParticleBase< ParticleInteractLayout< T, Dim, Mesh > > &)
virtual int broadcast_others(Message *, int, bool delmsg=true)
Message & get(const T &cval)
Message & put(const T &val)
touch_range_dv touch_range_rdv(const NDRegion< T, Dim > &domain)
~ParticleInteractLayout()
Message * receive_block(int &node, int &tag)
static Communicate * Comm
bool send(Message *, int node, int tag, bool delmsg=true)
void find_pairs(const unsigned LocalNum, const unsigned a1, const unsigned a2, const bool initLists, IpplParticleBase< ParticleInteractLayout< T, Dim, Mesh > > &PData)
void rebuild_interaction_data()
Inform & endl(Inform &inf)
T getMaxLocalInteractionRadius()