C++ ParticleList类代码示例

void ParticleFilterLocalizer::resample() {
    ParticleList resampledParticles;
	int particleCount = (int)particles.size();
    int index = Math::randomInt(0, particleCount - 1);
    float beta = 0.0f;
    float maxProbability = 1.0f;

    for (int i = 0; i < particleCount; i++) {
        beta += Math::randomFloat() * 2.0f * maxProbability;

        while (beta > particles[index]->probability) {
            beta -= particles[index]->probability;
            index = (index + 1) % particleCount;
        }

        resampledParticles.push_back(new Particle(
            particles[index]->x,
            particles[index]->y,
            particles[index]->orientation,
            particles[index]->probability
        ));
    }

    for (ParticleList::const_iterator it = particles.begin(); it != particles.end(); it++) {
        delete *it;
    }

    particles.clear();

    particles = resampledParticles;
}

void MDcalculation::calculateMeanDisplacement(){
	//Resetting
	for(int t = 0; t < 3; t++){
		rMeanSquared[t] = 0;
	}
	ParticleList* list;
	for(int i = 1; i < Nx + 1; i++){
	for(int j = 1; j < Ny + 1; j++){
	for(int k = 1; k < Nz + 1; k++){
		list = cells[i][j][k];
		for(Particle* p = list->getFirst(); list->hasNext(); p = list->getNext()){
			for(int d = 0; d < 3; d++){
				rMeanSquared[p->type] += p->rReal[d]*p->rReal[d];
			}
		}
	}
	}
	}
	MPI_Reduce(rMeanSquared, rMeanSquaredSum, 3, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);
	if(rank == 0){
		for(int t = 0; t < 3; t++){
			if(NparticlesType[t] == 0){
				rMeanSquaredSum[t] = 0;
			}else{
				rMeanSquaredSum[t] = rMeanSquaredSum[t]/NparticlesType[t];
			}
		}
	}
}

void AdvanceTime(ParticleList& list, float dt, float timeStep)
{
    Time += dt;

    for (size_t i = 0; i < list.size(); ++i) {
        Point3 p(list[i].Px, list[i].Py, list[i].Pz);
        Vector3 v = ComputeCurl(p);
        Point3 midx = p + 0.5f * timeStep * v;
        p += timeStep * ComputeCurl(midx);
        list[i].Px = p.getX();
        list[i].Py = p.getY();
        list[i].Pz = p.getZ();
        list[i].Vx = v.getX();
        list[i].Vy = v.getY();
        list[i].Vz = v.getZ();
    }

    for (ParticleList::iterator i = list.begin(); i != list.end();) {
        if (i->Py > PlumeCeiling) {
            i = list.erase(i);
        } else {
            ++i;
        }
    }

   noise.set_time(0.5f*NoiseGain[0]/NoiseLengthScale[0]*Time);

   if (!ShowStreamlines || Time < 0.1f)
        SeedParticles(list, dt);
}

void add_particles_mpi(ParticleList &particles, CellMatrix &cells)
{
    ParticleList::iterator iter = particles.begin();
    while(iter != particles.end())
    {
        particle_t tmp = (*iter);
        Point p = get_cell_index(tmp);
        cells[p.y][p.x].push_back(&(*iter));
        ++iter;
    }
}

ParticleList * ReadParticleListFromFile(char * file)
{
  ParticleList * ReadParticleList = 0; 
  PLNode       * CurrentNode      = 0;

  char   character;
  int    pid, itotal; 
  double E, px, py, pz, weight;

  itotal = 0;
  ReadParticleList = new ParticleList;
  CurrentNode      = ReadParticleList->GetHeadNode();

  cout << "Reading particle list from file: " 
       << file << endl;
  cout << endl;

  ifstream fin(file);
  while (fin.get(character)) 
  {
    fin >> pid;
    fin >> E;
    fin >> px;
    fin >> py;
    fin >> pz;
    fin >> weight;
    Particle * pParticle = new Particle;
    pParticle->Set4mom(E,px,py,pz);
    pParticle->SetID(pid);
    pParticle->SetWeight(weight);
    pParticle->SetDecaysum(0.0);
    pParticle->SetGeneration(1);
    if ( E>0.0 ) 
    { 
      itotal++;
      if ( fmod((double)itotal,10000.0)==0.0 )
	cout << itotal << " particles read from file" << endl;
      ReadParticleList->InsertAfter(CurrentNode, pParticle);
      CurrentNode = CurrentNode->GetNextNode();
     }
    E = 0.0;
  }

  fin.close();

  cout << endl;

  return ReadParticleList;
}

void beforeResort(ParticleList &pl, OutBuffer& out) {
  if (pl.size() == 1) {
    beforeResortCalled = true;
    int res = 42;
    out.write(res);
  }
}

void ParticleRender::drawParticles(const ParticleList &particles) {
	ci::gl::ScopedTextureBind	stb(mTexture);
	// Prevent writing to the depth buffer, which will block out
	// pixels that are supposed to be transparent.
	ci::gl::ScopedDepthWrite	sdw(false);
	ci::gl::ScopedBlendAlpha	sba;
	ci::gl::color(1.0f, 1.0f, 1.0f, 1.0f);

	size_t						k = 0, size = particles.size();
	while (size >= BUFFER_SIZE) {
		drawParticles(k, k + BUFFER_SIZE, particles);
		k += BUFFER_SIZE;
		size -= BUFFER_SIZE;
	}
	drawParticles(k, particles.size(), particles);
}

void afterResort(ParticleList &pl, InBuffer& inbuf) {
  if (pl.size() == 1) {
    afterResortCalled = true;
    int res;
    inbuf.read(res);
    BOOST_CHECK_EQUAL(res, 42);
  }
}

void MDcalculation::thermostat(double dt, double T0, double tau, double T){
	double gamma = sqrt(1 + dt/tau*(T0/T - 1));
	//Scaling velocities
	ParticleList* list;
	for(int i = 1; i <= Nx; i++){
	for(int j = 1; j <= Ny; j++){
	for(int k = 1; k <= Nz; k++){
		list = cells[i][j][k];
		for(Particle* p = list->getFirst(); list->hasNext(); p = list->getNext()){
			for(int d = 0; d < 3; d++){
				p->v[d] = p->v[d]*gamma;
			}
		}
	}
	}
	}
} 

void MDcalculation::calculateForces(bool monitoring){
	V = 0;
	ParticleList* list;
	//Calculating weights, finding 2BF and 3BF particle lists
	findNeighbours();
	// return;
	//Calculating 2BF and 3BF
	for(int i = 0; i < Nx + 2; i++){
	for(int j = 0; j < Ny + 2; j++){
	for(int k = 0; k < Nz + 2; k++){
		list = cells[i][j][k];
		if(i != 0 && i != Nx + 1 && j != 0 && j != Ny + 1 && k != 0 && k != Nz + 1){
			//Internal 3BF calculation
			for(Particle* p = list->getFirst(); list->hasNext(); p = list->getNext()){
				V += potential->force2(p);
				V += potential->force3(p);
			}
		}else{
			//3BF for boundary with internal
			for(Particle* p = list->getFirst(); list->hasNext(); p = list->getNext()){
				V += potential->boundaryForce3(p);
			}
		}
	}
	}
	}
}

  void FixedSingleList::afterRecvParticles(ParticleList &pl, InBuffer& buf) {
	longint size = pl.size();
    std::vector< longint > received(size);
    buf.read(received);
    for (std::vector< longint >::iterator it=received.begin(); it!=received.end(); it++) {
    	globalSingles.insert(*it);
    }
    LOG4ESPP_INFO(theLogger, "received fixed single list after receive particles");
  }

void Cell::removeLostParticles(ParticleList &lost) {
    for (ParticleList::iterator it = particles.begin(); it != particles.end();) {
        const Particle &particle = *it;
        if (!isInside(particle)) { // particle moved outside of cell
            lost.push_back(particle); // add to lost
            it = particles.erase(it); // erase particle
        } else
            it++;
    }
}

  void FixedSingleList::beforeSendParticles(ParticleList& pl, OutBuffer& buf) {
	std::vector< longint > toSend(pl.size());
    // loop over the particle list
    for (ParticleList::Iterator pit(pl); pit.isValid(); ++pit) {
      longint pid = pit->id();
      toSend.push_back(pid);
      globalSingles.erase(pid);
      LOG4ESPP_DEBUG(theLogger, "erase and send particle with pid from FixedSingleList" << pid);
    }
    // send the list
    buf.write(toSend);
    LOG4ESPP_INFO(theLogger, "prepared fixed single list before send particles");
  }

static void SeedParticles(ParticleList& list, float dt)
{
    static float time = 0;
    time += dt;
    unsigned int num_new = (unsigned int) (time * ParticlesPerSecond);
    if (num_new == 0)
        return;

    time = 0;
    list.reserve(list.size() + num_new);
    for (unsigned int i = 0; i < num_new; ++i) {
        float theta = randhashf(Seed++, 0, M_PI * 2);
        float r = randhashf(Seed++, 0, SeedRadius);
        float y = randhashf(Seed++, 0, InitialBand);
        Particle p;
        p.Px = r*std::cos(theta);
        p.Py = PlumeBase + y;
        p.Pz = r*std::sin(theta) + 0.125f; // Nudge the emitter towards the viewer
        p.ToB = Time;
        list.push_back(p);
    }
}

// apply force from actors to reactors
void Bin::apply_forces(ParticleList& reactors, ParticleList& actors) {
	for( ParticleIterator reactor = reactors.begin(); reactor != reactors.end(); reactor++ ) {
	if ((*reactor)->tag < 0)
	    continue;

        for (ParticleIterator actor = actors.begin(); actor != actors.end(); actor++ ){
            if ( *reactor == *actor)
                continue;
            double dx = (*actor)->x - (*reactor)->x;
            double dy = (*actor)->y - (*reactor)->y;
            double r2 = dx * dx + dy * dy;
            if( r2 > cutoff * cutoff )
                continue;
            r2 = fmax( r2, min_r*min_r );
            double r = sqrt( r2 );

            //  very simple short-range repulsive force
            double coef = ( 1 - cutoff / r ) / r2 / mass;
            (*reactor)->ax += coef * dx;
            (*reactor)->ay += coef * dy;
        }
    }
}