adding the solution of week8

exercises/week8/particles/solution/compute_verlet_integration.cc

+#include "compute_verlet_integration.hh"
+
+ComputeVerletIntegration::ComputeVerletIntegration(Real dt) : dt(dt) {}
+
+/* -------------------------------------------------------------------------- */
+
+void ComputeVerletIntegration::setDeltaT(Real dt) { this->dt = dt; }
+
+/* -------------------------------------------------------------------------- */
+
+void ComputeVerletIntegration::compute(System& system) {
+  UInt size = system.getNbParticles();
+
+  for (auto& par : system) {
+    par.getVelocity() += .5 * dt * par.getForce() / par.getMass();
+    par.getPosition() += dt * par.getVelocity();
+  }
+
+  auto& sun = system.getParticle(0);
+  sun.getPosition() = 0;
+
+  // first of all reset forces to zero
+  for (auto& particle : system)
+    particle.getForce() = 0;
+
+  for (auto& interaction : interactions)
+    interaction->compute(system);
+
+  for (auto& par : system) {
+    par.getVelocity() += .5 * dt * par.getForce() / par.getMass();
+  }
+}
+
+/* -------------------------------------------------------------------------- */
+
+void ComputeVerletIntegration::addInteraction(
+    std::shared_ptr<ComputeInteraction> interaction) {
+  interactions.push_back(interaction);
+}

exercises/week8/particles/solution/compute_verlet_integration.hh

+#ifndef __COMPUTE_VERLET_INTEGRATION__HH__
+#define __COMPUTE_VERLET_INTEGRATION__HH__
+
+/* -------------------------------------------------------------------------- */
+#include "compute.hh"
+#include "compute_interaction.hh"
+/* -------------------------------------------------------------------------- */
+#include <cmath>
+
+//! Integrate equation of motion
+class ComputeVerletIntegration : public Compute {
+  using InteractionList = std::vector<std::shared_ptr<ComputeInteraction>>;
+
+  // Constructors/Destructors
+public:
+  ComputeVerletIntegration(Real timestep);
+
+  // Methods
+public:
+  //! Set time step
+  void setDeltaT(Real dt);
+  //! Evolve positions and velocities
+  void compute(System& system) override;
+  //! Add an interaction to the computation of forces
+  void addInteraction(std::shared_ptr<ComputeInteraction> interaction);
+
+private:
+  Real dt;
+  InteractionList interactions;
+};
+
+/* -------------------------------------------------------------------------- */
+#endif  //__COMPUTE_VERLET_INTEGRATION__HH__

exercises/week8/particles/solution/csv_reader.cc

+#include "csv_reader.hh"
+#include "particles_factory_interface.hh"
+#include <fstream>
+#include <sstream>
+/* -------------------------------------------------------------------------- */
+CsvReader::CsvReader(const std::string& filename) : filename(filename) {}
+/* -------------------------------------------------------------------------- */
+void CsvReader::read(System& system) { this->compute(system); }
+/* -------------------------------------------------------------------------- */
+
+void CsvReader::compute(System& system) {
+  std::ifstream is(filename.c_str());
+  std::string line;
+
+  if (is.is_open() == false) {
+    std::cerr << "cannot open file " << filename << std::endl;
+    throw;
+  }
+
+  while (is.good()) {
+    getline(is, line);
+
+    if (line[0] == '#' || line.size() == 0)
+      continue;
+
+    auto p = ParticlesFactoryInterface::getInstance().createParticle();
+    std::stringstream sstr(line);
+    sstr >> *p;
+    system.addParticle(std::move(p));
+  }
+
+  is.close();
+}
+
+/* -------------------------------------------------------------------------- */

exercises/week8/particles/solution/csv_reader.hh

+#ifndef __CSV_READER__HH__
+#define __CSV_READER__HH__
+
+/* -------------------------------------------------------------------------- */
+#include "compute.hh"
+#include <string>
+
+//! Read system from csv input file
+class CsvReader : public Compute {
+
+  /* ------------------------------------------------------------------------ */
+  /* Constructors/Destructors                                                 */
+  /* ------------------------------------------------------------------------ */
+public:
+  //! Construct from filename
+  CsvReader(const std::string& filename);
+
+  /* ------------------------------------------------------------------------ */
+  /* Methods                                                                  */
+  /* ------------------------------------------------------------------------ */
+public:
+  //! Left here for convenience (calls compute)
+  void read(System& system);
+  //! Read input file into system
+  void compute(System& system);
+
+protected:
+  std::string filename;
+};
+
+/* -------------------------------------------------------------------------- */
+#endif  //__CSV_READER__HH__

exercises/week8/particles/solution/csv_writer.cc

+#include "csv_writer.hh"
+#include <cstdlib>
+#include <fstream>
+/* -------------------------------------------------------------------------- */
+CsvWriter::CsvWriter(const std::string& filename) : filename(filename) {}
+
+/* -------------------------------------------------------------------------- */
+void CsvWriter::write(System& system) { this->compute(system); }
+
+/* -------------------------------------------------------------------------- */
+
+void CsvWriter::compute(System& system) {
+  std::ofstream os(filename.c_str());
+
+  if (os.is_open() == false) {
+    std::cerr << "cannot open file " << filename << std::endl
+              << "check that the dumps folder exists" << std::endl;
+    std::exit(1);
+  }
+
+  UInt nb_particles = system.getNbParticles();
+
+  for (UInt p = 0; p < nb_particles; ++p) {
+    os << system.getParticle(p) << std::endl;
+  }
+
+  os.close();
+}

exercises/week8/particles/solution/csv_writer.hh

+#ifndef __CSV_WRITER__HH__
+#define __CSV_WRITER__HH__
+
+/* -------------------------------------------------------------------------- */
+#include "compute.hh"
+#include <string>
+
+//! Write system state to csv file
+class CsvWriter : public Compute {
+  // Constructors/Destructors
+public:
+  //! Construct from filename
+  CsvWriter(const std::string& filename);
+
+  // Methods
+public:
+  //! Write file (calls compute)
+  void write(System& system);
+  //! Write file
+  void compute(System& system);
+
+protected:
+  std::string filename;
+};
+
+/* -------------------------------------------------------------------------- */
+#endif  //__CSV_WRITER__HH__

exercises/week8/particles/solution/main.cc

+#include "compute_gravity.hh"
+#include "compute_verlet_integration.hh"
+#include "csv_reader.hh"
+#include "csv_writer.hh"
+#include "my_types.hh"
+#include "ping_pong_balls_factory.hh"
+#include "planets_factory.hh"
+#include "system.hh"
+/* -------------------------------------------------------------------------- */
+#include <cstdlib>
+#include <iostream>
+#include <sstream>
+/* -------------------------------------------------------------------------- */
+
+int main(int argc, char** argv) {
+  if (argc != 6) {
+    std::cout << "Usage: " << argv[0]
+              << " nsteps dump_freq input.csv particle_type timestep"
+              << std::endl;
+    std::cout << "\tparticle type can be: planet, ping_pong" << std::endl;
+    std::exit(EXIT_FAILURE);
+  }
+
+  // the number of steps to perform
+  Real nsteps;
+  std::stringstream(argv[1]) >> nsteps;
+  // freq to dump
+  int freq;
+  std::stringstream(argv[2]) >> freq;
+  // init file
+  std::string filename = argv[3];
+  // particle type
+  std::string type = argv[4];
+  // timestep
+  Real timestep;
+  std::stringstream(argv[5]) >> timestep;
+
+  if (type == "planet")
+    PlanetsFactory::getInstance();
+  else if (type == "ping_pong")
+    PingPongBallsFactory::getInstance();
+  else {
+    std::cout << "Unknown particle type: " << type << std::endl;
+    std::exit(EXIT_FAILURE);
+  }
+
+  ParticlesFactoryInterface& factory = ParticlesFactoryInterface::getInstance();
+
+  SystemEvolution& evol = factory.createSimulation(filename, timestep);
+
+  evol.setNSteps(nsteps);
+  evol.setDumpFreq(freq);
+
+  evol.evolve();
+
+  return EXIT_SUCCESS;
+}

exercises/week8/particles/solution/my_types.hh

+#ifndef __MY_TYPES_HH__
+#define __MY_TYPES_HH__
+
+/* -------------------------------------------------------------------------- */
+#include <iostream>
+
+typedef unsigned int UInt;
+typedef double Real;
+
+/* -------------------------------------------------------------------------- */
+#define TO_IMPLEMENT                                                           \
+  {                                                                            \
+    std::cerr << std::endl                                                     \
+              << std::endl                                                     \
+              << "*********************************************************\n" \
+              << "To be filled @ " << __FILE__ << ":" << __LINE__ << std::endl \
+              << "*********************************************************\n" \
+              << std::endl                                                     \
+              << std::endl;                                                    \
+    throw;                                                                     \
+  }
+/* -------------------------------------------------------------------------- */
+
+#define EXERCISE_BEGIN_CORRECTION
+#define EXERCISE_END_CORRECTION
+
+#endif /* __MY_TYPES_HH__ */

exercises/week8/particles/solution/particle.cc

+#include "particle.hh"
+
+void Particle::printself(std::ostream& stream) const {
+  stream << " " << position;
+  stream << " " << velocity;
+  stream << " " << force;
+  stream << " " << mass;
+}
+
+/* -------------------------------------------------------------------------- */
+
+void Particle::initself(std::istream& sstr) {
+  sstr >> position;
+  sstr >> velocity;
+  sstr >> force;
+  sstr >> mass;
+}

exercises/week8/particles/solution/particle.hh

+#ifndef __PARTICLE__HH__
+#define __PARTICLE__HH__
+
+/* -------------------------------------------------------------------------- */
+#include "vector.hh"
+/* -------------------------------------------------------------------------- */
+
+//! Particle base class
+class Particle {
+  /* ------------------------------------------------------------------------ */
+  /* Constructors/Destructors                                                 */
+  /* ------------------------------------------------------------------------ */
+public:
+  //! Default destructor (virtual because of virtual functions in this class)
+  virtual ~Particle() = default;
+
+  // Accessors
+public:
+  //! Get mass
+  Real& getMass() { return mass; }
+  //! Get position
+  Vector& getPosition() { return position; }
+  //! Get force
+  Vector& getForce() { return force; }
+  //! Get velocity
+  Vector& getVelocity() { return velocity; }
+
+  // I/O functions
+public:
+  /// function to print the content of the class
+  virtual void printself(std::ostream& stream) const;
+
+  /// function to initilize values from input files
+  virtual void initself(std::istream& sstr);
+
+  // Members
+protected:
+  Real mass;
+  Vector position, force, velocity;
+};
+
+/* -------------------------------------------------------------------------- */
+/* Inline functions */
+/* -------------------------------------------------------------------------- */
+
+inline std::istream& operator>>(std::istream& sstr, Particle& _this) {
+  _this.initself(sstr);
+  return sstr;
+}
+
+/* -------------------------------------------------------------------------- */
+inline std::ostream& operator<<(std::ostream& sstr, Particle& _this) {
+  _this.printself(sstr);
+  return sstr;
+}
+
+
+/* -------------------------------------------------------------------------- */
+#endif  //__PARTICLE__HH__

exercises/week8/particles/solution/particles_factory_interface.cc

+#include "particles_factory_interface.hh"
+#include "planets_factory.hh"
+/* -------------------------------------------------------------------------- */
+ParticlesFactoryInterface& ParticlesFactoryInterface::getInstance() {
+  return *factory;
+}
+
+/* -------------------------------------------------------------------------- */
+ParticlesFactoryInterface* ParticlesFactoryInterface::factory = nullptr;

exercises/week8/particles/solution/particles_factory_interface.hh

+#ifndef __PARTICLES_FACTORY_INTERFACE__HH__
+#define __PARTICLES_FACTORY_INTERFACE__HH__
+
+/* -------------------------------------------------------------------------- */
+#include "system_evolution.hh"
+/* -------------------------------------------------------------------------- */
+
+//! Abstract factory defining interface
+class ParticlesFactoryInterface {
+  // Constructors/Destructors
+protected:
+  //! Instance constructor (protected)
+  ParticlesFactoryInterface() = default;
+
+public:
+  virtual ~ParticlesFactoryInterface() = default;
+
+  // Methods
+public:
+  //! Create a whole simulation from file
+  virtual SystemEvolution& createSimulation(const std::string& fname,
+                                            Real timestep) = 0;
+  //! Create a new particle
+  virtual std::unique_ptr<Particle> createParticle() = 0;
+
+  //! Get singleton instance
+  static ParticlesFactoryInterface& getInstance();
+
+  // Members
+protected:
+  std::vector<Particle*> list_particles;
+  std::unique_ptr<SystemEvolution> system_evolution = nullptr;
+
+  // Standard pointer because constructor is protected (cannot use make_unique)
+  static ParticlesFactoryInterface* factory;
+};
+
+/* -------------------------------------------------------------------------- */
+#endif  //__PARTICLES_FACTORY_INTERFACE__HH__

exercises/week8/particles/solution/ping_pong_ball.cc

+#include "ping_pong_ball.hh"
+
+/* -------------------------------------------------------------------------- */
+void PingPongBall::printself(std::ostream& stream) const {
+  Particle::printself(stream);
+
+  stream << " " << radius;
+}
+
+/* -------------------------------------------------------------------------- */
+
+void PingPongBall::initself(std::istream& sstr) {
+  Particle::initself(sstr);
+
+  sstr >> radius;
+}

exercises/week8/particles/solution/ping_pong_ball.hh

+#ifndef __PING_PONG_BALL__HH__
+#define __PING_PONG_BALL__HH__
+
+/* -------------------------------------------------------------------------- */
+#include "particle.hh"
+
+//! Class for ping-pong ball
+class PingPongBall : public Particle {
+  /* ------------------------------------------------------------------------ */
+  /* Methods                                                                  */
+  /* ------------------------------------------------------------------------ */
+
+public:
+  //! Get contact dissipation
+  Real& getContactDissipation() { return contact_dissipation; }
+  //! Get ball radius
+  Real& getRadius() { return radius; }
+
+  void printself(std::ostream& stream) const override;
+  void initself(std::istream& sstr) override;
+
+private:
+  Real radius, contact_dissipation;
+};
+
+/* -------------------------------------------------------------------------- */
+#endif  //__PING_PONG_BALL__HH__

exercises/week8/particles/solution/ping_pong_balls_factory.cc

+#include "ping_pong_balls_factory.hh"
+#include "compute_contact.hh"
+#include "compute_verlet_integration.hh"
+#include "csv_reader.hh"
+#include "csv_writer.hh"
+#include "ping_pong_ball.hh"
+#include <cmath>
+#include <iostream>
+/* -------------------------------------------------------------------------- */
+
+std::unique_ptr<Particle> PingPongBallsFactory::createParticle() {
+  return std::make_unique<PingPongBall>();
+}
+
+/* -------------------------------------------------------------------------- */
+
+SystemEvolution&
+PingPongBallsFactory::createSimulation(const std::string& fname,
+                                       Real timestep) {
+  this->system_evolution =
+      std::make_unique<SystemEvolution>(std::make_unique<System>());
+
+  CsvReader reader(fname);
+  reader.read(this->system_evolution->getSystem());
+
+  auto contact = std::make_shared<ComputeContact>();
+  auto verlet = std::make_shared<ComputeVerletIntegration>(timestep);
+
+  contact->setPenalty(1.);
+  verlet->addInteraction(contact);
+  this->system_evolution->addCompute(verlet);
+
+  return *system_evolution;
+}
+
+/* -------------------------------------------------------------------------- */
+
+ParticlesFactoryInterface& PingPongBallsFactory::getInstance() {
+  if (not ParticlesFactoryInterface::factory)
+    ParticlesFactoryInterface::factory = new PingPongBallsFactory;
+
+  return *factory;
+}
+
+/* -------------------------------------------------------------------------- */

exercises/week8/particles/solution/ping_pong_balls_factory.hh

+#ifndef __PING_PING_BALLS_FACTORY__HH__
+#define __PING_PING_BALLS_FACTORY__HH__
+
+/* -------------------------------------------------------------------------- */
+#include "particles_factory_interface.hh"
+#include "ping_pong_ball.hh"
+/* -------------------------------------------------------------------------- */
+
+//! Factory for ping-pong balls
+class PingPongBallsFactory : public ParticlesFactoryInterface {
+  /* ------------------------------------------------------------------------ */
+  /* Constructors/Destructors                                                 */
+  /* ------------------------------------------------------------------------ */
+private:
+  PingPongBallsFactory() = default;
+
+  /* ------------------------------------------------------------------------ */
+  /* Methods                                                                  */
+  /* ------------------------------------------------------------------------ */
+public:
+  SystemEvolution& createSimulation(const std::string& fname,
+				    Real timestep) override;
+
+  std::unique_ptr<Particle> createParticle() override;
+
+  static ParticlesFactoryInterface& getInstance();
+};
+
+/* -------------------------------------------------------------------------- */
+#endif  //__PING_PING_BALLS_FACTORY__HH__

exercises/week8/particles/solution/planet.cc

+#include "planet.hh"
+
+void Planet::initself(std::istream &stream) {
+  Particle::initself(stream);
+  stream >> name;
+}
+
+/* -------------------------------------------------------------------------- */
+
+void Planet::printself(std::ostream &stream) const {
+  Particle::printself(stream);
+  stream << " " << name;
+}

exercises/week8/particles/solution/planet.hh

+#ifndef __PLANET__HH__
+#define __PLANET__HH__
+
+/* -------------------------------------------------------------------------- */
+#include "particle.hh"
+/* -------------------------------------------------------------------------- */
+#include <iostream>
+/* -------------------------------------------------------------------------- */
+
+//! Class for planet
+class Planet : public Particle {
+  // Methods
+public:
+  //! Get name
+  std::string& getName() { return name; }
+
+  // Inherited methods
+public:
+  void initself(std::istream& stream) override;
+  void printself(std::ostream& stream) const override;
+
+private:
+  std::string name;
+};
+
+/* -------------------------------------------------------------------------- */
+#endif  //__PLANET__HH__

exercises/week8/particles/solution/planets_factory.cc

+#include "planets_factory.hh"
+#include "compute_gravity.hh"
+#include "compute_verlet_integration.hh"
+#include "csv_reader.hh"
+#include "csv_writer.hh"
+#include "planet.hh"
+#include <cmath>
+/* -------------------------------------------------------------------------- */