Unverified Commit bb210cae authored by Guillaume Anciaux's avatar Guillaume Anciaux
Browse files

pushing exercise of this week

parent cbcfdbb8
step-*.csv
build
starting_point
\ No newline at end of file
cmake_minimum_required (VERSION 2.6)
project (Particles)
set(CMAKE_CXX_STANDARD 14)
add_library(part
vector.cc
compute_boundary.cc
compute_verlet_integration.cc
particle.cc
planet.cc
compute_gravity.cc
csv_reader.cc
particles_factory_interface.cc
planets_factory.cc
compute_contact.cc
compute_kinetic_energy.cc
csv_writer.cc
system.cc
compute_energy.cc
compute_potential_energy.cc
ping_pong_ball.cc
system_evolution.cc
ping_pong_balls_factory.cc
compute_interaction.cc)
add_subdirectory(googletest)
add_executable(particles main.cc)
target_link_libraries(particles gtest_main gtest pthread part)
add_executable(test_kepler test_kepler.cc)
target_link_libraries(test_kepler gtest_main gtest pthread part)
add_custom_target(test ./test_kepler DEPENDS part test_kepler)
This diff is collapsed.
0 0 0 0 0 0 1e-4 0 0 1 sun
1 0 0 0 1 0 -1e-4 0 0 1e-4 earth
#ifndef __COMPUTE__HH__
#define __COMPUTE__HH__
/* -------------------------------------------------------------------------- */
#include "system.hh"
/* -------------------------------------------------------------------------- */
//! Base class for all compute
class Compute {
public:
//! Virtual destructor needed because we have subclasses
virtual ~Compute() = default;
/* ------------------------------------------------------------------------ */
/* Methods */
/* ------------------------------------------------------------------------ */
public:
//! Compute is pure virtual
virtual void compute(System& system) = 0;
};
/* -------------------------------------------------------------------------- */
#endif //__COMPUTE__HH__
#include "compute_boundary.hh"
/* -------------------------------------------------------------------------- */
ComputeBoundary::ComputeBoundary(const Vector& xmin, const Vector& xmax) {}
/* -------------------------------------------------------------------------- */
void ComputeBoundary::compute(System& system) {}
/* -------------------------------------------------------------------------- */
#ifndef __COMPUTE_BOUNDARY__HH__
#define __COMPUTE_BOUNDARY__HH__
/* -------------------------------------------------------------------------- */
#include "compute.hh"
/* -------------------------------------------------------------------------- */
//! Compute interaction with simulation box
class ComputeBoundary : public Compute {
// Constructors/Destructors
public:
ComputeBoundary(const Vector& xmin, const Vector& xmax);
// Methods
public:
void compute(System& system) override;
// Members
protected:
Vector xmin, xmax;
};
/* -------------------------------------------------------------------------- */
#endif //__COMPUTE_BOUNDARY__HH__
#include "compute_contact.hh"
#include "ping_pong_ball.hh"
#include <cmath>
/* -------------------------------------------------------------------------- */
void ComputeContact::setPenalty(Real penalty) {}
/* -------------------------------------------------------------------------- */
void ComputeContact::computePairInteraction(Particle& par1, Particle& par2) {}
#ifndef __COMPUTE_CONTACT__HH__
#define __COMPUTE_CONTACT__HH__
/* -------------------------------------------------------------------------- */
#include "compute_interaction.hh"
//! Compute contact interaction between ping-pong balls
class ComputeContact : public ComputeInteraction {
// Virtual implementation
public:
//! Penalty contact implementation
void computePairInteraction(Particle& par1, Particle& par2) override;
// Accessors
public:
//! Set penalty
void setPenalty(Real penalty);
// Members
protected:
Real penalty;
};
/* -------------------------------------------------------------------------- */
#endif //__COMPUTE_CONTACT__HH__
#include "compute_energy.hh"
/* -------------------------------------------------------------------------- */
#ifndef __COMPUTE_ENERGY__HH__
#define __COMPUTE_ENERGY__HH__
/* -------------------------------------------------------------------------- */
#include "compute.hh"
//! Base class for energy computation
class ComputeEnergy : public Compute {
// Methods
public:
Real getValue() { return value; }
protected:
Real value;
};
/* -------------------------------------------------------------------------- */
#endif //__COMPUTE_ENERGY__HH__
#include "compute_gravity.hh"
#include <cmath>
/* -------------------------------------------------------------------------- */
void ComputeGravity::setG(Real G) { this->G = G; }
/* -------------------------------------------------------------------------- */
void ComputeGravity::computePairInteraction(Particle& par1, Particle& par2) {
// compute the distance vector and the square of distance
auto v_r = par2.getPosition() - par1.getPosition();
auto r2 = v_r.squaredNorm();
if (r2 == 0.)
return;
// compute the distance
auto r = sqrt(r2);
v_r *= 1 / r;
// compute the pair force
auto force = par1.getMass() * par2.getMass() * G / r2 * v_r;
// add up the force for both concerned particles
par2.getForce() -= force;
par1.getForce() += force;
}
#ifndef __COMPUTE_GRAVITY__HH__
#define __COMPUTE_GRAVITY__HH__
/* -------------------------------------------------------------------------- */
#include "compute_interaction.hh"
//! Compute Newton gravity interaction
class ComputeGravity : public ComputeInteraction {
// Virtual implementation
public:
//! Newton gravity implementation
void computePairInteraction(Particle& par1, Particle& par2) override;
// Accessors
public:
//! set the gravitational constant
void setG(Real G);
// Members
private:
//! newton constant
Real G = 1.;
};
/* -------------------------------------------------------------------------- */
#endif //__COMPUTE_GRAVITY__HH__
#include "compute_interaction.hh"
#include <cmath>
/* -------------------------------------------------------------------------- */
void ComputeInteraction::compute(System& system) {
UInt size = system.getNbParticles();
// first of all reset forces to zero
for (auto& particle : system)
particle.getForce() = 0;
for (UInt p1 = 0; p1 < size; ++p1) {
Particle& par1 = system.getParticle(p1);
for (UInt p2 = p1 + 1; p2 < size; ++p2) {
Particle& par2 = system.getParticle(p2);
computePairInteraction(par1, par2);
}
}
}
#ifndef __COMPUTE_INTERACTION__HH__
#define __COMPUTE_INTERACTION__HH__
/* -------------------------------------------------------------------------- */
#include "compute.hh"
//! Base class for interaction computation
class ComputeInteraction : public Compute {
public:
void compute(System& system) override;
//! Compute interaction for a pair of particles
virtual void computePairInteraction(Particle& par1, Particle& par2) = 0;
};
/* -------------------------------------------------------------------------- */
#endif //__COMPUTE_INTERACTION__HH__
#include "compute_kinetic_energy.hh"
/* -------------------------------------------------------------------------- */
void ComputeKineticEnergy::compute(System& system) {}
/* -------------------------------------------------------------------------- */
#ifndef __COMPUTE_KINETIC_ENERGY__HH__
#define __COMPUTE_KINETIC_ENERGY__HH__
/* -------------------------------------------------------------------------- */
#include "compute_energy.hh"
//! Compute kinetic energy of system
class ComputeKineticEnergy : public ComputeEnergy {
public:
void compute(System& system) override;
};
/* -------------------------------------------------------------------------- */
#endif //__COMPUTE_KINETIC_ENERGY__HH__
#include "compute_potential_energy.hh"
/* -------------------------------------------------------------------------- */
ComputePotentialEnergy::ComputePotentialEnergy(ComputeInteraction& cForces)
: cForces(cForces) {}
/* -------------------------------------------------------------------------- */
void ComputePotentialEnergy::compute(System& system) {}
#ifndef __COMPUTE_POTENTIAL_ENERGY__HH__
#define __COMPUTE_POTENTIAL_ENERGY__HH__
/* -------------------------------------------------------------------------- */
#include "compute_energy.hh"
#include "compute_interaction.hh"
/* -------------------------------------------------------------------------- */
//! Compute potential energy of system
class ComputePotentialEnergy : public ComputeEnergy {
/* ------------------------------------------------------------------------ */
/* Constructors/Destructors */
/* ------------------------------------------------------------------------ */
public:
ComputePotentialEnergy(ComputeInteraction& cForces);
/* ------------------------------------------------------------------------ */
/* Methods */
/* ------------------------------------------------------------------------ */
public:
void compute(System& system) override;
/* ------------------------------------------------------------------------ */
/* Members */
/* ------------------------------------------------------------------------ */
protected:
ComputeInteraction& cForces;
};
/* -------------------------------------------------------------------------- */
#endif //__COMPUTE_POTENTIAL_ENERGY__HH__
#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);
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment