Migration des particules vers les faisceau.

Il reste une subtilite a gerer: pour des raisons de compatibilite, un accelerateur contient toujours une methode getParticles(). Celle-ci renvoi une concatenation de tous les particules de tous les faisceaux de l'accelerateur.
Or il y a un probleme d'allocation de la memoire avec cette methode (temporairement repare dans cette revision mais degeu!). Je vois deux solutions:
1) changer tous les test qui appellent la methode getParticles()
2) essayer de reparer le probleme d'allocation (difficile et je ne sais pas si c'est possible)

14 files changed, 328 insertions, 167 deletions

diff --git a/Makefile b/Makefile
index da10e0f..72ab508 100644
--- a/Makefile
+++ b/Makefile
@@ -25,7 +25,7 @@ export BINDIR = $(BASEDIR)/bin
 #
 # CXXFLAGS += -Wall
 # CXXFLAGS += -ansi -pedantic -Wall   # pour les purs et durs
-# CXXFLAGS += -g                      # pour debugger
+CXXFLAGS += -g                      # pour debugger
 # CXXFLAGS += -pg                     # pour profiler
 # LDFLAGS  += -pg                     # pour profiler
  CXXFLAGS += -O2                     # pour optimiser la vitesse
diff --git a/src/gui/Main.cc b/src/gui/Main.cc
index 28ba09f..94f68ee 100644
--- a/src/gui/Main.cc
+++ b/src/gui/Main.cc
@@ -21,6 +21,7 @@
 #include "Vector3D.h"
 #include "constants.h"
 #include <vector>
+#include "CircularBeam.h"
 
 using namespace std;
 using namespace vhc;
@@ -129,7 +130,6 @@ Une particule :
 	//acc->add(e7);
 	acc->add(e8);
 
-	acc->close();
 
 	//proton
 	Particle p1 = Particle(Vector3D(3.01, 0, 0), constants::PROTON_MASS, constants::E, 2 * constants::GeV, -Vector3D::j);
@@ -143,7 +143,9 @@ Une particule :
 	acc->add(ap1);
 	acc->add(ap2);
 
+	acc->add(CircularBeam(p1, 10, 1));
 
+	acc->close();
 
 	/*std::vector< Particle > ps = createParticles(e1.getEntryPosition(), 1000);
 
diff --git a/src/gui/Stage.cc b/src/gui/Stage.cc
index 81c405d..020f31a 100644
--- a/src/gui/Stage.cc
+++ b/src/gui/Stage.cc
@@ -118,7 +118,7 @@ void Stage::paintGL() {
 		glBlendFunc(GL_SRC_ALPHA, GL_ONE);
 		glEnable(GL_BLEND);
 		glDisable(GL_DEPTH_TEST);
-		for (list<Element*>::const_iterator i = accelerator->getElements().begin(); i != accelerator->getElements().end(); ++i) {
+		for (Accelerator::ElementCollection::const_iterator i = accelerator->getElements().begin(); i != accelerator->getElements().end(); ++i) {
 			elementRenderer.render(**i);
 		}
 
@@ -127,7 +127,8 @@ void Stage::paintGL() {
 
 		glColor3d(0, 0, 1);
 
-		for (list<Particle*>::const_iterator i = accelerator->getParticles().begin(); i != accelerator->getParticles().end(); ++i) {
+		Accelerator::ParticleCollection particles = accelerator->getParticles();
+		for (Accelerator::ParticleCollection::const_iterator i = particles.begin(); i != particles.end(); ++i) {
 			particleRenderer.render(**i);
 		}
 	}
@@ -144,7 +145,7 @@ void Stage::paintGL() {
 		camera.move(mv);
 
 
-	if (!paused) accelerator->step(h);
+	if (!paused) for (int i = 0; i < 10; ++i) accelerator->step(h);
 
 	glColor3d(1,1,0);
 	util::crosshair();
diff --git a/src/gui/gui.pro b/src/gui/gui.pro
index 74c6a72..03b954b 100644
--- a/src/gui/gui.pro
+++ b/src/gui/gui.pro
@@ -10,8 +10,8 @@ LIBS += -L$(BINDIR)/main -lvhc
 QT += opengl
 
 #CONFIG += DEBUG
-#QMAKE_CXXFLAGS_DEBUG += -pg
-#QMAKE_LFLAGS_DEBUG += -pg
+#QMAKE_CXXFLAGS_DEBUG += -g
+#QMAKE_LFLAGS_DEBUG += -g
 
 
 # Input
diff --git a/src/main/Accelerator.cc b/src/main/Accelerator.cc
index ca5de8a..40b9ff7 100644
--- a/src/main/Accelerator.cc
+++ b/src/main/Accelerator.cc
@@ -8,14 +8,15 @@
 #include <iostream>
 #include "Accelerator.h"
 #include "exceptions.h"
+#include "SingleBeam.h"
 
 using namespace std;
 
 namespace vhc {
 
 Accelerator::Accelerator():
-		elementCollec(0),
-		particleCollec(0),
+		elements(0),
+		beams(0),
 		allowLinear(false),
 		closed(false)
 	{};
@@ -24,51 +25,95 @@ Accelerator::~Accelerator() {
 	clear();
 }
 
-void Accelerator::initializeParticles() {
-	//rajouter les particules dans leurs elements respectifs
-	for (ParticleIterator i = particleCollec.begin(); i != particleCollec.end(); ++i) {
+void Accelerator::initializeBeams() {
+	//pour chauqe faisceau (vide a priori)
+	for (BeamCollection::iterator i = beams.begin(); i != beams.end(); ++i) {
 
-		for (ElementIterator j = elementCollec.begin(); j != elementCollec.end(); ++j) {
-			if ((**j).contains(**i)) {
-				(**i).setElement(*j);
+		(**i).getReferenceParticle().setElement(NULL);
+
+		//rajouter les particules de reference dans leurs elements respectifs
+		for (ElementIterator j = elements.begin(); j != elements.end(); ++j) {
+			if ((**j).contains((**i).getReferenceParticle())) {
+				(**i).getReferenceParticle().setElement(*j);
 				break;
 			}
 		}
 
-		//si une particule n'est pas contenue dans un element elle est supprimee
-		if ((**i).getElement() == NULL) {
+		//si une particule de reference n'est pas contenue dans un element le faisceau est supprimee
+		//et une exception est lancee
+		if ((**i).getReferenceParticle().getElement() == NULL) {
 			delete *i;
-			i = particleCollec.erase(i);
+			i = beams.erase(i);
 			--i;
-			//std::cout << "Particle hit wall. Removed from simulation" << std::endl;
+		} else {
+			//sinon, le faisceau est initialise
+			(**i).clear();
+			(**i).initializeParticles();
+
+			//mettre chaque particule dans l'element respectif
+			for (ParticleCollection::iterator j = (**i).getParticles().begin(); j != (**i).getParticles().end(); ++j) {
+				(**j).setElement(NULL);
+				for (ElementIterator k = elements.begin(); k != elements.end(); ++k) {
+					if ((**k).contains(**j)) {
+						(**j).setElement(*k);
+						break;
+					}
+				}
+
+				//supprimer la particule si elle n'est pas contenue dans un element
+				if ((**j).getElement() == NULL) {
+					delete *j;
+					j = (**i).getParticles().erase(j);
+					--j;
+				}
+			}
 		}
+
 	}
 }
 
+
 Element& Accelerator::add(const Element& element) {
 	Element* e = element.clone();
-	elementCollec.push_back(e);
+	elements.push_back(e);
 	closed = false;
 	return *e;
 }
 
-Particle& Accelerator::add(const Particle& particle) {
-	Particle* p = particle.clone();
-	particleCollec.push_back(p);
+Beam& Accelerator::add(const Beam& beam) {
+	Beam* b = beam.clone();
+	beams.push_back(b);
 	closed = false;
-	return *p;
+	return *b;
 }
 
-const Accelerator::ElementCollection& Accelerator::getElements() const { return elementCollec;}
+Beam& Accelerator::add(const Particle& particle) {
+	return add(SingleBeam(particle));
+}
 
-const Accelerator::ParticleCollection& Accelerator::getParticles() const { return particleCollec;}
+const Accelerator::ElementCollection& Accelerator::getElements() const {
+	return elements;
+}
+
+const Accelerator::BeamCollection& Accelerator::getBeams() const {
+	return beams;
+}
+
+//TODO the segmentation fault came from here! find a solution
+const Accelerator::ParticleCollection& Accelerator::getParticles() const {
+	ParticleCollection* particles = new ParticleCollection(0);
+	for (BeamCollection::const_iterator i = beams.begin(); i != beams.end(); ++i) {
+		particles->insert(particles->end(), (**i).getParticles().begin(), (**i).getParticles().end());
+	}
+	return *particles;
+}
 
 void Accelerator::close() {
-	for (ElementIterator current = elementCollec.begin(); current != elementCollec.end(); ++current) {
+	for (ElementIterator current = elements.begin(); current != elements.end(); ++current) {
 
 		ElementIterator next = current;
 		++next;
-		if (next == elementCollec.end()) next = elementCollec.begin();
+		if (next == elements.end()) next = elements.begin();
 
 
 		// est-ce que les elements se suivent (sont connectes)?
@@ -77,63 +122,31 @@ void Accelerator::close() {
 			(**next).setPrevious(*current);
 
 		//sinon est-ce qu'il s'agit du dernier element?
-		} else if (next == elementCollec.begin()) {
+		} else if (next == elements.begin()) {
 			if (!allowLinear) throw UnsupportedOperationException("Cannot close accelerator. Linear Accelerators are not allowed.");
 
 		//sinon
 		} else throw UnsupportedOperationException("Cannot close accelerator. Two succeeding elements are not physically connected. (not close enough)");
 	}
 
-	initializeParticles();
+	initializeBeams();
 
 	closed = true;
 
 }
 
-void Accelerator::updateParticles() {
-	for (ParticleIterator i = particleCollec.begin(); i != particleCollec.end(); ++i) {
-		Particle& particle = **i;
-		if (particle.getElement()->isAfter(particle)) {
-			if (particle.getElement()->getNext() == NULL)
-				if (allowLinear) {
-					delete *i;
-					i = particleCollec.erase(i);
-					--i;
-					//cout << "Particle reached end of accelerator. Removed from simulation" << std::endl;
-				}
-				else throw Exception("Element in accelerator not connected to next.");
-			else particle.setElement(particle.getElement()->getNext());
-		} else if (particle.getElement()->isBefore(particle)) {
-			if (particle.getElement()->getPrevious() == NULL)
-				if (allowLinear) {
-					delete *i;
-					i = particleCollec.erase(i);
-					--i;
-					//cout << "Particle reached beginning of accelerator. Removed from simulation" << std::endl;
-				}
-				else throw Exception("Element in accelerator not connected to previous.");
-			else particle.setElement(particle.getElement()->getPrevious());
-		} else if (particle.getElement()->isBeside(particle)) {
-			//std::cout << "Particle hit wall. Removed from simulation" << std::endl;
-			delete *i;
-			i = particleCollec.erase(i);
-			--i;
-		}
-	}
-}
-
 void Accelerator::clear() {
-	for (ParticleIterator i = particleCollec.begin(); i != particleCollec.end(); ++i) {
+	for (BeamCollection::iterator i = beams.begin(); i != beams.end(); ++i) {
 		delete *i;
 		*i = NULL;
 	}
-	particleCollec.clear();
+	beams.clear();
 
-	for (ElementIterator i = elementCollec.begin(); i != elementCollec.end(); ++i) {
+	for (ElementCollection::iterator i = elements.begin(); i != elements.end(); ++i) {
 		delete *i;
 		*i = NULL;
 	}
-	elementCollec.clear();
+	elements.clear();
 
 	closed = false;
 }
@@ -141,22 +154,9 @@ void Accelerator::clear() {
 void Accelerator::step(double dt) {
 	if (!closed) close();
 
-	for (ParticleIterator i = particleCollec.begin(); i != particleCollec.end(); ++i) {
-		Particle& particle = **i;
-
-		particle.setForce(Vector3D::Null);
-
-		particle.applyMagneticForce(particle.getElement()->magneticFieldAt(particle.getPosition()), dt);
-
-		Vector3D a = particle.getForce() / (particle.getGamma() * particle.getMass());
-		particle.setVelocity(particle.getVelocity() + a * dt);
-
-		particle.translate(particle.getVelocity() * dt);
-
+	for (BeamCollection::iterator i = beams.begin(); i != beams.end(); ++i) {
+		(**i).step(dt);
 	}
-
-
-	updateParticles();
 }
 
 void Accelerator::enableLinear(bool value) {
@@ -167,26 +167,28 @@ void Accelerator::enableLinear(bool value) {
 std::string Accelerator::toString() const {
 	std::stringstream s;
 
-	if (elementCollec.size() == 0) {
+	if (elements.size() == 0) {
 		s << "This accelerator doesn't contain anything."<<"\n";
-	} else if (elementCollec.size()==1) {
+	} else if (elements.size()==1) {
 		s << "This accelerator is made of the following element :" << "\n";
-		s << elementCollec.front()->toString()<<"\n";
+		s << elements.front()->toString()<<"\n";
 	} else {
-		s << "This accelerator is made of the " << elementCollec.size() << " following elements :" << "\n";
-		for (ElementCollection::const_iterator i = elementCollec.begin(); i != elementCollec.end(); ++i) {
+		s << "This accelerator is made of the " << elements.size() << " following elements :" << "\n";
+		for (ElementCollection::const_iterator i = elements.begin(); i != elements.end(); ++i) {
 			s << (*i)->toString() << "\n";
 		}
 	}
 
-	if (particleCollec.size() == 0) {
+	const ParticleCollection& particles = getParticles();
+
+	if (particles.size() == 0) {
 		s << "This accelerator doesn't contain any particle." << "\n";
-	} else if (particleCollec.size() == 1) {
+	} else if (particles.size() == 1) {
 		s << "This accelerator contains the following particle :" << "\n";
-		s << particleCollec.front()->toString() << "\n";
+		s << particles.front()->toString() << "\n";
 	}else{
-		s << "This accelerator contains the " << particleCollec.size() << " following particles :"<<"\n";
-		for (list<Particle*>::const_iterator i = particleCollec.begin(); i != particleCollec.end(); ++i) {
+		s << "This accelerator contains the " << particles.size() << " following particles :"<<"\n";
+		for (ParticleCollection::const_iterator i = particles.begin(); i != particles.end(); ++i) {
 					s << (*i)->toString() << "\n";
 		}
 	}
diff --git a/src/main/Accelerator.h b/src/main/Accelerator.h
index 175529b..9f990af 100644
--- a/src/main/Accelerator.h
+++ b/src/main/Accelerator.h
@@ -10,6 +10,7 @@
 #include <list>
 #include "Vector3D.h"
 #include "Particle.h"
+#include "Beam.h"
 #include "Element.h"
 
 namespace vhc {
@@ -22,6 +23,7 @@ public:
 
 	typedef std::list<Particle*> ParticleCollection;
 	typedef std::list<Element*>  ElementCollection;
+	typedef std::list<Beam*> BeamCollection;
 
 	typedef ParticleCollection::iterator ParticleIterator;
 	typedef ElementCollection::iterator ElementIterator;
@@ -35,8 +37,9 @@ public:
 	 *  L'accelerateur est ouvert en ajoutant un element. */
 	Element& add(const Element& element);
 
-	/** Copie une particule dans l'accélérateur. */
-	Particle& add(const Particle& particle);
+	Beam& add(const Beam& beam);
+
+	Beam& add(const Particle& particle);
 
 	/** Retourne la liste d'elements contenus dans cet accelerateur.
 	 *  <b>ATTENTION:</b> les elements peuvent etre supprimes sans preavis par l'accelerateur! */
@@ -46,6 +49,8 @@ public:
 	 *  <b>ATTENTION:</b> les particules peuvent etre supprimes sans preavis par l'accelerateur! */
 	const ParticleCollection & getParticles() const;
 
+	const BeamCollection& getBeams() const;
+
 	/** Ferme l'accelerateur.
 	 *  En invoquant cette methode, la continuite des elements est verifiee et les particules sont attribues leurs elements respectifs.
 	 *  @throws UnsupportedOperationException si les elements sauf le dernier ne sont pas continus
@@ -78,10 +83,9 @@ private:
 protected:
 
 	/** Collection d'elements contenus dans cet accelerateur. */
-	ElementCollection elementCollec;
+	ElementCollection elements;
 
-	/** Collection de particules contenus danc cet accelerateur. */
-	ParticleCollection particleCollec;
+	BeamCollection beams;
 
 	/** Autorise les accelerateurs lineaires.
 	 *  @see enableLinear */
@@ -91,19 +95,11 @@ protected:
 	 *  @see close */
 	bool closed;
 
-	/** Initialize les particules en leur attribuant l'element dans lequel ils sont contenus.
-	 *  Les particules non-contenus sont supprimes de l'accelerateur. */
-	void initializeParticles();
-
-	/** Met a jour les particules en leur attribuant l'element dans lequel ils sont contenus.
-	 *  Contrairement a <code>initializeParticles()</code>, les elements consideres sont:
-	 *  - l'element actuel de la particule
-	 *  - l'element precedent
-	 *  - l'element suivant
-	 *  Si la particule se situe a cote de son element, elle est supprimee de l'accelerateur.
-	 *  Attention: si la particule saute un element, elle est tout de meme consideree comme etant dans l'element suivant (ou precedent)! Ceci
-	 *  peut survenir si un element est trop petit ou si la simulation est faite avec des pas de temps trop grands. */
-	void updateParticles();
+	/** Initialize les faisceaux.
+	 *  Les faisceaux dont la particule de refernce n'est pas contenue dans un element
+	 *  sont supprimes de l'accelerateur. */
+	void initializeBeams();
+
 };
 
 }
diff --git a/src/main/Beam.cc b/src/main/Beam.cc
index f01368e..97c2688 100644
--- a/src/main/Beam.cc
+++ b/src/main/Beam.cc
@@ -13,18 +13,27 @@ using namespace std;
 
 namespace vhc {
 
-Beam::Beam(const Particle& referenceParticle, int quantity, int lambda): referenceParticle(referenceParticle), quantity(quantity), lambda(lambda) {
-	init(referenceParticle, quantity, lambda);
-}
+Beam::Beam(const Particle& referenceParticle, int quantity, int lambda):
+		referenceParticle(referenceParticle),
+		quantity(quantity),
+		lambda(lambda) {}
 
 Beam::~Beam() {
-	for (list<Particle*>::iterator i = particles.begin(); i != particles.end(); ++i) {
+	clear();
+}
+
+void Beam::clear() {
+	for (ParticleCollection::iterator i = particles.begin(); i != particles.end(); ++i) {
 		delete *i;
 		*i = NULL;
 	}
 	particles.clear();
 }
 
+int Beam::getQuantity() const {
+	return quantity;
+}
+
 int Beam::getSize() const {
 	return particles.size();
 }
@@ -33,10 +42,55 @@ int Beam::getLambda() const {
 	return lambda;
 }
 
-const Particle& Beam::getReferenceParticle() const {
+Particle& Beam::getReferenceParticle() {
 	return referenceParticle;
 }
 
+Beam::ParticleCollection& Beam::getParticles() {return particles;}
+
+void Beam::updateParticles() {
+	for (ParticleCollection::iterator i = particles.begin(); i != particles.end(); ++i) {
+		Particle& particle = **i;
+		if (particle.getElement()->isAfter(particle)) {
+			if (particle.getElement()->getNext() == NULL) {
+				delete *i;
+				i = particles.erase(i);
+				--i;
+				//cout << "Particle reached end of accelerator. Removed from simulation" << std::endl;
+			} else particle.setElement(particle.getElement()->getNext());
+		} else if (particle.getElement()->isBefore(particle)) {
+			if (particle.getElement()->getPrevious() == NULL) {
+				delete *i;
+				i = particles.erase(i);
+				--i;
+				//cout << "Particle reached beginning of accelerator. Removed from simulation" << std::endl;
+			} else particle.setElement(particle.getElement()->getPrevious());
+		} else if (particle.getElement()->isBeside(particle)) {
+			//std::cout << "Particle hit wall. Removed from simulation" << std::endl;
+			delete *i;
+			i = particles.erase(i);
+			--i;
+		}
+	}
+}
+
+void Beam::step(double dt) {
+	for (ParticleCollection::iterator i = particles.begin(); i != particles.end(); ++i) {
+		Particle& particle = **i;
+
+		particle.applyMagneticForce(particle.getElement()->magneticFieldAt(particle.getPosition()), dt);
+
+		Vector3D a = particle.getForce() / (particle.getGamma() * particle.getMass());
+		particle.setVelocity(particle.getVelocity() + a * dt);
+
+		particle.translate(particle.getVelocity() * dt);
+
+		particle.setForce(Vector3D::Null);
+	}
+
+	updateParticles();
+}
+
 double Beam::getVR2() const {
 	double z2 = 0;
 	for (list<Particle*>::const_iterator i = particles.begin(); i != particles.end(); ++i) {
diff --git a/src/main/Beam.h b/src/main/Beam.h
index f03af97..70f9f86 100644
--- a/src/main/Beam.h
+++ b/src/main/Beam.h
@@ -14,52 +14,12 @@
 
 namespace vhc {
 
-class Beam {
-
-protected:
-
-	/** Particule de reference. */
-	Particle referenceParticle;
-
-	/** Particules contenus dans ce faisceau. */
-	std::list<Particle*> particles;
-
-	/** Coefficient des macrosparticules. */
-	int lambda;
-
-	int quantity;
-
-	/** Retourne la moyenne de la distribution horizontale de la position des particules.
-	 *  (<r^2> horizontal du complement) */
-	double getHR2() const;
-
-	/** Retourne la moyenne de la distribution horizontale de la vitesse des particules.
-	 *  (<v^2> horizontal du complement) */
-	double getHV2() const;
-
-	/** Retourne la moyenne de la distribution horizontale du produit de la position et de la vitesse des particules.
-	 *  (<r*v>^2 horizontal du complement) */
-	double getHRV2() const;
-
-	/** Retourne la moyenne de la distribution verticale de la position des particules.
-	 *  (<r^2> vertical du complement) */
-	double getVR2() const;
-
-	/** Retourne la moyenne de la distribution verticale de la vitesse des particules.
-	 *  (<v^2> vertical du complement) */
-	double getVV2() const;
-
-	/** Retourne la moyenne de la distribution verticale du produit de la position et de la vitesse des particules.
-	 *  (<r*v>^2 vertical du complement) */
-	double getVRV2() const;
-
-	/** Initialise ce faisceau a partir de la particule de reference. Comme il s'agit d'une sorte de factory method (au sens vraiment large),
-	 *  la responsabilite de gestion des particules appartient a ce faisceau (et non a la sous-classe). */
-	virtual void init(const Particle& referenceParticle, int quantity, int lambda) = 0;
-
+class Beam: public Cloneable {
 
 public:
 
+	typedef std::list<Particle*> ParticleCollection;
+
 	/** Cree un nouveau faisceaux.
 	 *  <b>ATTENTION:</b> un `Beam' est abstrait et n'initialise pas les particules a partir de la reference.
 	 *  C'est au client de faire cela! */
@@ -67,14 +27,34 @@ public:
 
 	virtual ~Beam();
 
+
+	ParticleCollection& getParticles();
+
+	/** Met a jour les particules en leur attribuant l'element dans lequel ils sont contenus.
+	 *  Contrairement a <code>Accelerator::initializeBeams()</code>, les elements consideres sont:
+	 *  - l'element actuel de la particule
+	 *  - l'element precedent
+	 *  - l'element suivant
+	 *  Si la particule se situe a cote de son element, elle est supprimee de l'accelerateur.
+	 *  Attention: si la particule saute un element, elle est tout de meme consideree comme etant dans l'element suivant (ou precedent)! Ceci
+	 *  peut survenir si un element est trop petit ou si la simulation est faite avec des pas de temps trop grands. */
+	void updateParticles();
+
+	virtual void initializeParticles() = 0;
+
+	void step(double dt);
+
 	/** Retourne la quantite de particules contenus dans ce faisceau. */
+	int getQuantity() const;
+
+	/** Retourne la quantite de macroparticules contenus dans ce faisceau. */
 	int getSize() const;
 
 	/** Retourne le coefficient des macroparticules. */
 	int getLambda() const;
 
 	/** Retourne la particule de reference. */
-	const Particle& getReferenceParticle() const;
+	Particle& getReferenceParticle();
 
 	/** Retourne l'energie moyenne des particules dans ce faisceau. */
 	double getAverageEnergy() const;
@@ -85,6 +65,15 @@ public:
 	/** Retourne l'emmitance horizontale de ce faisceau. */
 	double getHorizontalEmittance() const;
 
+	void clear();
+
+	virtual Beam* clone() const = 0;
+
+
+
+
+	//-------------------------------------------------------------------
+
 	/** Retourne coefficient des ellipses de phases vertical. */
 	double getVerticalA11() const;
 
@@ -103,6 +92,45 @@ public:
 	/** Retourne coefficient des ellipses de phases horizontal. */
 	double getHorizontalA22() const;
 
+protected:
+
+	/** Particule de reference. */
+	Particle referenceParticle;
+
+	/** (Macro-)Particules contenus dans ce faisceau. */
+	ParticleCollection particles;
+
+	/** Coefficient des macrosparticules. */
+	int lambda;
+
+	int quantity;
+
+	/** Retourne la moyenne de la distribution horizontale de la position des particules.
+	 *  (<r^2> horizontal du complement) */
+	double getHR2() const;
+
+	/** Retourne la moyenne de la distribution horizontale de la vitesse des particules.
+	 *  (<v^2> horizontal du complement) */
+	double getHV2() const;
+
+	/** Retourne la moyenne de la distribution horizontale du produit de la position et de la vitesse des particules.
+	 *  (<r*v>^2 horizontal du complement) */
+	double getHRV2() const;
+
+	/** Retourne la moyenne de la distribution verticale de la position des particules.
+	 *  (<r^2> vertical du complement) */
+	double getVR2() const;
+
+	/** Retourne la moyenne de la distribution verticale de la vitesse des particules.
+	 *  (<v^2> vertical du complement) */
+	double getVV2() const;
+
+	/** Retourne la moyenne de la distribution verticale du produit de la position et de la vitesse des particules.
+	 *  (<r*v>^2 vertical du complement) */
+	double getVRV2() const;
+
+
+
 };
 
 }
diff --git a/src/main/CircularBeam.cc b/src/main/CircularBeam.cc
index ee2fa9d..fefc3f9 100644
--- a/src/main/CircularBeam.cc
+++ b/src/main/CircularBeam.cc
@@ -6,6 +6,8 @@
  */
 
 #include "CircularBeam.h"
+#include "Particle.h"
+#include "Element.h"
 
 namespace vhc {
 
@@ -15,7 +17,26 @@ CircularBeam::CircularBeam(const Particle& referenceParticle, int quantity, int
 
 CircularBeam::~CircularBeam() {}
 
-CircularBeam::init(const Particle& referenceParticle, int quantity, int lambda) {
+void CircularBeam::initializeParticles() {
+	Element* element = referenceParticle.getElement();
+
+	for (int i = 0; i < quantity / lambda; ++i) {
+
+		particles.push_back(new Particle(
+				element->getEntryPosition(),
+				referenceParticle.getMass() * lambda,
+				referenceParticle.getCharge() * lambda,
+				referenceParticle.getEnergy(),
+				element->getDiagonal()
+				));
+
+		element = element->getNext();
+	}
 
 }
+
+CircularBeam* CircularBeam::clone() const {
+	return new CircularBeam(referenceParticle, quantity, lambda);
+}
+
 }
diff --git a/src/main/CircularBeam.h b/src/main/CircularBeam.h
index 086abfb..853774d 100644
--- a/src/main/CircularBeam.h
+++ b/src/main/CircularBeam.h
@@ -14,14 +14,13 @@ namespace vhc {
 
 class CircularBeam: public Beam {
 
-protected:
-
-	virtual void init(const Particle& referenceParticle, int quantity, int lambda) = 0;
 
 public:
 	CircularBeam(const Particle& referenceParticle, int quantity, int lambda);
 	virtual ~CircularBeam();
 
+	virtual void initializeParticles();
+
 	virtual CircularBeam* clone() const;
 };
 
diff --git a/src/main/Makefile b/src/main/Makefile
index 83385fd..96f4d98 100644
--- a/src/main/Makefile
+++ b/src/main/Makefile
@@ -14,7 +14,8 @@ LOCALDIR = main
 # Si un objet nécessite une compilation non-standard (i.e. pas de règle du style Foo.o : Foo.cc Foo.h), rajouter
 # cette règle.
 LOCALOBJS = Vector3D.o Particle.o Printable.o Element.o CurvedElement.o StraightElement.o \
-			CompositeElement.o Dipole.o Quadrupole.o FODO.o ElementVisitor.o Cloneable.o Accelerator.o exceptions.o Stepper.o Beam.o
+			CompositeElement.o Dipole.o Quadrupole.o FODO.o ElementVisitor.o Cloneable.o \
+			Accelerator.o exceptions.o Stepper.o Beam.o SingleBeam.o CircularBeam.o
 
 OBJS=$(addprefix $(BINDIR)/$(LOCALDIR)/,$(LOCALOBJS))
 .PHONY = all checkdirs lib
diff --git a/src/main/SingleBeam.cc b/src/main/SingleBeam.cc
new file mode 100644
index 0000000..fc5f63e
--- /dev/null
+++ b/src/main/SingleBeam.cc
@@ -0,0 +1,27 @@
+/*
+ * SingleBeam.cc
+ *
+ *  Created on: May 18, 2011
+ *      Author: jakob
+ */
+
+#include "SingleBeam.h"
+#include "iostream"
+
+namespace vhc {
+
+SingleBeam::SingleBeam(const Particle& particle): Beam(particle, 1, 1) {}
+
+SingleBeam::~SingleBeam() {
+
+}
+
+void SingleBeam::initializeParticles() {
+	particles.push_back(new Particle(referenceParticle));
+}
+
+SingleBeam* SingleBeam::clone() const {
+	return new SingleBeam(referenceParticle);
+}
+
+}
diff --git a/src/main/SingleBeam.h b/src/main/SingleBeam.h
new file mode 100644
index 0000000..95451de
--- /dev/null
+++ b/src/main/SingleBeam.h
@@ -0,0 +1,28 @@
+/*
+ * SingleBeam.h
+ *
+ *  Created on: May 18, 2011
+ *      Author: jakob
+ */
+
+#ifndef SINGLEBEAM_H_
+#define SINGLEBEAM_H_
+
+#include "Beam.h"
+#include "Particle.h"
+
+namespace vhc {
+
+class SingleBeam: public vhc::Beam {
+public:
+	SingleBeam(const Particle& particle);
+	virtual ~SingleBeam();
+
+	virtual void initializeParticles();
+
+	virtual SingleBeam* clone() const;
+};
+
+}
+
+#endif /* SINGLEBEAM_H_ */
diff --git a/src/test/AcceleratorBenchmarkTest.cc b/src/test/AcceleratorBenchmarkTest.cc
index 671cf03..520d3fe 100644
--- a/src/test/AcceleratorBenchmarkTest.cc
+++ b/src/test/AcceleratorBenchmarkTest.cc
@@ -17,6 +17,7 @@
 #include "FODO.h"
 #include "Vector3D.h"
 #include "constants.h"
+#include "SingleBeam.h"
 
 using namespace vhc;
 using namespace std;
@@ -57,8 +58,6 @@ Accelerator* standard() {
 	acc->add(e7);
 	acc->add(e8);
 
-	acc->close();
-
 	Particle p1 = Particle(Vector3D(3.00, 0, 0), constants::PROTON_MASS, constants::E, 2 * constants::GeV, -Vector3D::j);
 	Particle p2 = Particle(Vector3D(2.99, 0, 0), constants::PROTON_MASS, constants::E, 2 * constants::GeV, -Vector3D::j);
 	acc->add(p1);
@@ -70,7 +69,7 @@ Accelerator* standard() {
 		acc->add(ps[i]);
 	}
 
-
+	acc->close();
 	return acc;
 }
 
@@ -82,7 +81,10 @@ int main() {
 	int steps = 1000;
 	double dt = 1E-11;
 
-	cout << "Simulating " << steps << " steps with " << accelerator->getParticles().size() << " particles...";
+	cout << "Simulating " << steps << " steps with " << accelerator->getParticles().size() << " particles in " << accelerator->getBeams().size() << " beams...";
+	for (Accelerator::BeamCollection::const_iterator i = accelerator->getBeams().begin(); i != accelerator->getBeams().end(); ++i) {
+		cout << (**i).getParticles().size();
+		}
 	cout.flush();
 	int t0 = clock();
 	for (int i = 0; i < steps; ++i) {