MercuryDPM  Trunk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
SubcriticalMaserBoundaryTEST.cc
Go to the documentation of this file.
1 //Copyright (c) 2013-2020, The MercuryDPM Developers Team. All rights reserved.
2 //For the list of developers, see <http://www.MercuryDPM.org/Team>.
3 //
4 //Redistribution and use in source and binary forms, with or without
5 //modification, are permitted provided that the following conditions are met:
6 // * Redistributions of source code must retain the above copyright
7 // notice, this list of conditions and the following disclaimer.
8 // * Redistributions in binary form must reproduce the above copyright
9 // notice, this list of conditions and the following disclaimer in the
10 // documentation and/or other materials provided with the distribution.
11 // * Neither the name MercuryDPM nor the
12 // names of its contributors may be used to endorse or promote products
13 // derived from this software without specific prior written permission.
14 //
15 //THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
16 //ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
17 //WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
18 //DISCLAIMED. IN NO EVENT SHALL THE MERCURYDPM DEVELOPERS TEAM BE LIABLE FOR ANY
19 //DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
20 //(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
21 //LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
22 //ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
23 //(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24 //SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25 
26 
27 #include "DPMBase.h"
29 
34 {
35  distanceLeft_ = std::numeric_limits<Mdouble>::quiet_NaN();
36  distanceRight_ = std::numeric_limits<Mdouble>::quiet_NaN();
37  maserIsActivated_ = false;
38  activationTime_ = std::numeric_limits<Mdouble>::quiet_NaN();
39  copyFlowParticles_ = false;
40 }
41 
46 {
47  logger(DEBUG, "SubcriticalMaserBoundaryTEST::~SubcriticalBoundary() finished");
48 }
49 
56 {
57  return new SubcriticalMaserBoundaryTEST(*this);
58 }
59 
64 void SubcriticalMaserBoundaryTEST::read(std::istream& is)
65 {
67  std::string dummy;
68  is >> dummy >> maserIsActivated_
69  >> dummy >> activationTime_;
70  logger(INFO, "maser is activated: %", maserIsActivated_);
71 }
72 
77 void SubcriticalMaserBoundaryTEST::write(std::ostream& os) const
78 {
80  os << " maserIsActivated " << maserIsActivated_
81  << " activationTime " << activationTime_;
82 }
83 
89 {
90  return "SubcriticalMaserBoundaryTEST";
91 }
92 
94 {
95  copyFlowParticles_ = copyFlowParticles;
96 }
97 
102 {
103  if (maserIsActivated_)
104  {
105  for (BaseParticle* p : getHandler()->getDPMBase()->particleHandler)
106  {
107  if (getDistance(p->getPosition()) > 0)
108  {
109  p->setMaserParticle(true);
110  }
111  }
112  }
113 }
114 
128 {
130  const Mdouble proximityDistance = p->getMaxInteractionRadius() + 2.0* pH.getLargestParticle()->getMaxInteractionRadius();
131  if (maserIsActivated_)
132  {
133  if (p->isMaserParticle())
134  {
135  //Only if the particle is close to the right wall a ghost can be generated
136  if (getDistanceFromRight(p->getPosition()) < proximityDistance)
137  {
139  }
140  }
141  }
142  else
143  {
144  //If a particle is at any wall, the ghost can be generated
145  if (getDistance(p->getPosition()) < proximityDistance)
146  {
148  }
149  }
150 }
151 
162 {
163  // check if particle passed either of the boundary walls
164  if (maserIsActivated_)
165  {
166  if (p->isMaserParticle() && (getDistanceFromRight(p->getPosition()) < 0))
167  {
168  //Do the maser stuff
169  BaseParticle* pCopy = p->copy();
170  pCopy->setMaserParticle(false);
171  pH.addObject(pCopy);
172  /*
173  logger(INFO, "copying particle");
174  p->write(std::cout);
175  */
176  logger(VERBOSE, "copying particle %", p);
177  //Shift position of the particle
178  shiftPosition(p);
179  }
180  }
181  else
182  {
183  //note, if we already made a rough bottom on beforehand, we do not want the bottom particles to be eaten up,
184  //that is, all shifted inside the domain.
185  if (getDistance(p->getPosition()) < 0)
186  {
187  shiftPosition(p);
188  }
189  }
190 
191  return false;
192 }
193 
200 {
201 #ifdef MERCURY_USE_MPI
202  if (NUMBER_OF_PROCESSORS == 1)
203  {
204 #endif
205  //Activate the maser at the given time
206  if (getHandler()->getDPMBase()->getTime() > activationTime_)
207  {
208  if (!maserIsActivated_)
209  {
210  logger(INFO, "Activating the maser at time: %", activationTime_);
211  activateMaser();
212  }
213  }
214 
215  //Update particles
216  for (BaseParticle* p : pH)
217  {
219  }
220 #ifdef MERCURY_USE_MPI
221  }
222 
223 #endif
224 }
225 
226 
232 {
233  logger(INFO, "SubcriticalMaserBoundaryTEST::activateMaser ");
234 #ifdef MERCURY_USE_MPI
235  //Clear the periodic boundaryHandler, only want real particles
237 #endif
238  //Flag that the maser is active!
239  maserIsActivated_ = true;
240 
241  extendBottom();
242  if (copyFlowParticles_)
243  {
245  }
246 
247  //Loop over all particles to see if they are within the maser boundary
249  for (BaseParticle* particle : pH)
250  {
251  particle->setMaserParticle((getDistance(particle->getPosition()) > 0));
252  }
253 
254 #ifdef MERCURY_USE_MPI
255  //Generate ghost particles
257 #endif
258 
259 }
260 
262 {
263  if (maserIsActivated_)
264  {
265  for (BaseParticle* particle : getHandler()->getDPMBase()->particleHandler)
266  {
267  if (getDistance(particle->getPosition()) > 0)
268  {
269  particle->setMaserParticle(false);
270  }
271  }
272 
273  maserIsActivated_ = false;
274 
275 
276  // TODO JMFT: @Marnix In activateMaser(), what do extendBottom() and
277  // addNewParticles() do, and how do I undo their effects?
278  }
279  else
280  logger(WARN, "[SubcriticalMaserBoundaryTEST::deactivateMaser()] Maser is not activated, so can't deactivate");
281 }
282 
284 {
285  return maserIsActivated_;
286 }
287 
294 {
295  activationTime_ = time;
296 }
297 
307 {
308  if (maserIsActivated_)
309  {
310  return getDistanceFromRight(position);
311  }
312  else
313  {
314  return PeriodicBoundary::getDistance(position);
315  }
316 }
317 
324 {
325  Mdouble distanceFromPlaneThroughOrigin = Vec3D::dot(position, normal_);
326  return distanceRight_ - distanceFromPlaneThroughOrigin;
327 }
328 
340 void SubcriticalMaserBoundaryTEST::modifyPeriodicComplexity(std::vector<int>& complexity, int& totalPeriodicComplexity,
341  BaseParticle* particle, int i) const
342 {
343  if (maserIsActivated_)
344  {
345  if (particle->isMaserParticle())
346  {
347  //Check if this particle is flagged periodic in this boundary
348  if (complexity[i] < 0)
349  {
350  //Make sure that only a ghost close to the right boundary is flagged as real in this boundary
351  if (getDistanceFromRight(particle->getPosition()) < 0)
352  {
353  complexity[i] = 3;
354  }
355 
356  }
357  if (complexity[i] == -3)
358  {
359  logger(INFO, "Something went wrong in SubcriticalMaserBoundaryTEST::modifyPeriodicComplexity, "
360  "complexity[%] = -3", i);
361  }
362  }
363  else
364  {
365  if (complexity[i] == 1)
366  {
367  totalPeriodicComplexity--;
368  }
369  complexity[i] = 3;
370  }
371  }
372 }
373 
375 {
376  if (maserIsActivated_)
377  {
378  if (!particle->isMaserParticle())
379  {
380  std::vector<int> periodicComplexity = particle->getPeriodicComplexity();
381  periodicComplexity[i] = 3;
382  particle->setPeriodicComplexity(periodicComplexity);
383  }
384  }
385 }
386 
387 
394 {
395  //Activate the maser at the given time
396  if (getHandler()->getDPMBase()->getTime() > activationTime_)
397  {
398  if (!maserIsActivated_)
399  {
400  logger(INFO, "Activating the maser at time: %", activationTime_);
401  activateMaser();
402  }
403  }
404 }
405 
411 {
412  logger(INFO, "extending bottom");
413 #ifdef MERCURY_USE_MPI
414  MPIContainer& communicator = MPIContainer::Instance();
415  std::vector<unsigned int> numberOfParticlesPerCore(NUMBER_OF_PROCESSORS);
416  std::vector<std::vector<BaseParticle*>> particlesToCores(NUMBER_OF_PROCESSORS);
417  if (PROCESSOR_ID == 0)
418  {
419 #endif
420  std::vector<BaseParticle*> fixedParticles;
421  std::vector<BaseParticle*> newParticles;
422  for (BaseParticle* p : getHandler()->getDPMBase()->particleHandler)
423  {
424  if (p->isFixed() && !(p->isPeriodicGhostParticle()) && !(p->isMPIParticle()))
425  {
426  fixedParticles.push_back(p);
427  }
428  }
429 
430  for (BaseParticle* p : fixedParticles)
431  {
432  Vec3D newPosition = p->getPosition() + shift_;
433  Vec3D maxDomain = getHandler()->getDPMBase()->getMax();
434  Vec3D minDomain = getHandler()->getDPMBase()->getMin();
435  while (newPosition.X < maxDomain.X && newPosition.Y < maxDomain.Y && newPosition.Z < maxDomain.Z
436  && newPosition.X > minDomain.X && newPosition.Y > minDomain.Y && newPosition.Z > minDomain.Z)
437  {
438  p = p->copy();
439  p->setPosition(newPosition);
440  newParticles.push_back(p);
441  newPosition += shift_;
442  }
443  }
444 #ifndef MERCURY_USE_MPI
445  for (BaseParticle* p : newParticles)
446  {
448  }
449 #else
450  if (NUMBER_OF_PROCESSORS == 1)
451  {
452  for (BaseParticle* p : newParticles)
453  {
455  }
456  }
457  else
458  {
459  //count how many particles go to each core
460  for (BaseParticle* p : newParticles)
461  {
462  int targetGlobalIndex = getHandler()->getDPMBase()->domainHandler.getParticleDomainGlobalIndex(p);
463  int targetProcessor = getHandler()->getDPMBase()->domainHandler.getParticleProcessor(
464  targetGlobalIndex);
465  particlesToCores[targetProcessor].push_back(p);
466  }
467  for (int i = 0; i < NUMBER_OF_PROCESSORS; ++i)
468  {
469  numberOfParticlesPerCore[i] = (particlesToCores[i].size());
470  }
471  }
472 }
473 if (NUMBER_OF_PROCESSORS > 1)
474 {
475  //broadcast numberOfParticlesPerCore to other processors
476  communicator.broadcast(numberOfParticlesPerCore.data(), NUMBER_OF_PROCESSORS, 0);
477  for (unsigned i = 0; i < numberOfParticlesPerCore.size(); ++i)
478  {
479  for (unsigned p = 0; p < numberOfParticlesPerCore[i]; ++p)
480  {
481  BaseParticle* particle;
482  if (PROCESSOR_ID == 0)
483  {
484  particle = particlesToCores[i][p];
485  }
486  else
487  {
488  particle = new SphericalParticle();
489  }
490 
491  getHandler()->getDPMBase()->particleHandler.addObject(0, particle);
492  }
493  }
494 }
495 #endif
496 }
497 
503 {
504  logger(INFO, "copying flow particles");
505 #ifdef MERCURY_USE_MPI
506  MPIContainer& communicator = MPIContainer::Instance();
507  std::vector<unsigned int> numberOfParticlesPerCore(NUMBER_OF_PROCESSORS);
508  std::vector<std::vector<BaseParticle*>> particlesToCores(NUMBER_OF_PROCESSORS);
509  if (PROCESSOR_ID == 0)
510  {
511 #endif
512  std::vector<BaseParticle*> flowParticles;
513  std::vector<BaseParticle*> newParticles;
514  for (BaseParticle* p : getHandler()->getDPMBase()->particleHandler)
515  {
516  if (!p->isFixed() && !(p->isPeriodicGhostParticle()) && !(p->isMPIParticle()))
517  {
518  flowParticles.push_back(p);
519  }
520  }
521 
522  for (BaseParticle* p : flowParticles)
523  {
524  Vec3D newPosition = p->getPosition() + shift_;
525  for (unsigned i = 0; i < 4; ++i)
526  {
527  p = p->copy();
528  p->setPosition(newPosition);
529  newParticles.push_back(p);
530  newPosition += shift_;
531  }
532  }
533 #ifndef MERCURY_USE_MPI
534  for (BaseParticle* p : newParticles)
535  {
537  }
538 #else
539  if (NUMBER_OF_PROCESSORS == 1)
540  {
541  for (BaseParticle* p : newParticles)
542  {
544  }
545  }
546  else
547  {
548  //count how many particles go to each core
549  for (BaseParticle* p : newParticles)
550  {
551  int targetGlobalIndex = getHandler()->getDPMBase()->domainHandler.getParticleDomainGlobalIndex(p);
552  int targetProcessor = getHandler()->getDPMBase()->domainHandler.getParticleProcessor(
553  targetGlobalIndex);
554  particlesToCores[targetProcessor].push_back(p);
555  }
556  for (int i = 0; i < NUMBER_OF_PROCESSORS; ++i)
557  {
558  numberOfParticlesPerCore[i] = (particlesToCores[i].size());
559  }
560  }
561 }
562 if (NUMBER_OF_PROCESSORS > 1)
563 {
564  //broadcast numberOfParticlesPerCore to other processors
565  communicator.broadcast(numberOfParticlesPerCore.data(), NUMBER_OF_PROCESSORS, 0);
566  for (unsigned i = 0; i < numberOfParticlesPerCore.size(); ++i)
567  {
568  for (unsigned p = 0; p < numberOfParticlesPerCore[i]; ++p)
569  {
570  BaseParticle* particle;
571  if (PROCESSOR_ID == 0)
572  {
573  particle = particlesToCores[i][p];
574  }
575  else
576  {
577  particle = new SphericalParticle();
578  }
579 
580  getHandler()->getDPMBase()->particleHandler.addObject(0, particle);
581  //particle->setPeriodicComplexity(getPeriodicHandler()->computePeriodicComplexity(particle->getPosition()));
582  }
583  }
584 }
585 #endif
586 
587  logger(INFO, "completed copying particles");
588 }
This class contains all information and functions required for communication between processors...
Definition: MpiContainer.h:125
void shiftPosition(BaseParticle *p) const override
shifts the particle
BaseParticle * getLargestParticle() const
Returns the pointer of the largest particle in the particle handler. When mercury is running in paral...
SubcriticalMaserBoundaryTEST * copy() const override
Creates a copy of this maser on the heap.
const Vec3D & getPosition() const
Returns the position of this BaseInteractable.
static MPIContainer & Instance()
fetch the instance to be used for communication
Definition: MpiContainer.h:130
Mdouble X
the vector components
Definition: Vector.h:65
A basic particle.
void addObject(BaseParticle *P) override
Adds a BaseParticle to the ParticleHandler.
bool isMaserParticle() const
Indicates if this particle belongs to the maser boundary.
bool isActivated() const
Returns whether the maser is activated or not.
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
double Mdouble
Definition: GeneralDefine.h:34
void write(std::ostream &os) const override
writes boundary properties to ostream
Mdouble getMaxInteractionRadius() const
Returns the particle's interaction radius, which might be different from radius_ (e.g., when dealing with wet particles)
Definition: BaseParticle.h:359
Mdouble getDistanceFromRight(const Vec3D &position) const
returns the distance to the right wall
Vec3D getMin() const
Definition: DPMBase.h:623
void activateMaser()
Activates the maser functionaly of this periodic boundary.
const std::complex< Mdouble > i
Definition: ExtendedMath.h:50
~SubcriticalMaserBoundaryTEST() override
destructor
PeriodicBoundaryHandler * getPeriodicHandler() const
Returns the periodic boundary handler.
DomainHandler domainHandler
An object of the class DomainHandler which deals with parallel code.
Definition: DPMBase.h:1354
static Mdouble dot(const Vec3D &a, const Vec3D &b)
Calculates the dot product of two Vec3D: .
Definition: Vector.cc:76
virtual BaseParticle * copy() const =0
Particle copy method. It calls to copy constructor of this Particle, useful for polymorfism.
Mdouble distanceRight_
position of right edge, s.t. normal*x = distanceRight_
bool maserIsActivated_
Flag whether or not the gap is created and particles transformed already.
Mdouble getDistance(const BaseParticle &p) const override
Returns the distance of the edge to the particle.
Vec3D normal_
outward unit normal vector for right wall
void modifyPeriodicComplexity(std::vector< int > &complexity, int &totalPeriodicComplexity, BaseParticle *particle, int i) const override
modifies the periodic complexity to support a maser boundary
void deactivateMaser()
Stops copying particles, and act merely as a periodic domain.
void checkBoundaryAfterParticlesMove(ParticleHandler &pH) override
Evaluates what the particles have to do after they have changed position.
void modifyGhostAfterCreation(BaseParticle *particle, int i) override
std::string getName() const override
Returns the name of the object.
void read(std::istream &is) override
reads boundary properties from istream
void addNewParticles()
Adds new particles to the periodic particle lists.
ParticleHandler particleHandler
An object of the class ParticleHandler, contains the pointers to all the particles created...
Definition: DPMBase.h:1329
void createGhostParticle(BaseParticle *pReal)
Creates and adds a ghost particle from a give real particle.
void performActionsBeforeAddingParticles() override
Checks before adding particles if the maser needs to be activated.
#define PROCESSOR_ID
Definition: GeneralDefine.h:63
BoundaryHandler * getHandler() const
Returns the boundary's BoundaryHandler.
Mdouble distanceLeft_
position of left edge, s.t. normal*x = distanceLeft_
void setActivationTime(Mdouble time)
sets the activate time of the maser
bool checkBoundaryAfterParticleMoved(BaseParticle *p, ParticleHandler &pH) const
Shifts the particle to its 'periodic' position if it is a maser particle and has crossed either of th...
SubcriticalMaserBoundaryTEST()
MaserBoundary constructor.
Vec3D getMax() const
Definition: DPMBase.h:629
Container to store all BaseParticle.
#define NUMBER_OF_PROCESSORS
For the MPI communication routines this quantity is often required. defining this macro makes the cod...
Definition: GeneralDefine.h:62
Mdouble Y
Definition: Vector.h:65
Vec3D shift_
shift from left to right boundary
std::enable_if< std::is_scalar< T >::value, void >::type broadcast(T &t, int fromProcessor=0)
Broadcasts a scalar from the root to all other processors.
Definition: MpiContainer.h:429
int getParticleDomainGlobalIndex(BaseParticle *particle)
void setMaserParticle(bool flag)
Flags the status of the particle if it belongs to the maser boundary or not.
const std::vector< int > & getPeriodicComplexity()
Obtains the periodic communication complexity of the particle.
int getParticleProcessor(int globalIndex)
void setPeriodicComplexity(std::vector< int > complexity)
Set the periodic communication complexity of the particle.
void createPeriodicParticle(BaseParticle *p, ParticleHandler &pH) override
Creates periodic particles when the particle is a maser particle and is sufficiently close to one of ...
Definition: Vector.h:49
void write(std::ostream &os) const override
writes boundary properties to ostream
Mdouble activationTime_
Time at which the maser opens.
DPMBase * getDPMBase()
Gets the problem that is solved using this handler.
Definition: BaseHandler.h:725
Mdouble getDistance(const Vec3D &position) const override
gets the distance to the closest wall if maser is inactive, otherwise distance to right wall ...
Mdouble Z
Definition: Vector.h:65
void read(std::istream &is) override
reads boundary properties from istream
void clearCommunicationLists()
Removes all ghost particles and bookkeeping for a fresh start.
bool copyFlowParticles_
Flag for whether or not we copy a few blocks of flow particles in the front when activating the maser...
void setCopyFlowParticles(bool copyFlowParticles)