Commit de963aa6 authored by Manon Eugénie Voisin--Leprince's avatar Manon Eugénie Voisin--Leprince
Browse files

exercises week11

parent 0e67bed0
---
Language: Cpp
BasedOnStyle: LLVM
Standard: Cpp11
AlwaysBreakTemplateDeclarations: true
PointerAlignment: Left
SpacesBeforeTrailingComments: 2
#FixNamespaceComments: true
...
step-*.csv
build
starting_point
\ No newline at end of file
[submodule "googletest"]
path = googletest
url = https://github.com/google/googletest.git
cmake_minimum_required (VERSION 3.1)
project (Particles)
set(CMAKE_CXX_STANDARD 14)
################################################################
# libpart
################################################################
add_library(part
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
material_point.cc
system_evolution.cc
ping_pong_balls_factory.cc
compute_interaction.cc
compute_temperature.cc
material_points_factory.cc
)
add_executable(particles main.cc)
target_link_libraries(particles gtest_main gtest pthread part)
if(EXISTS "${PROJECT_SOURCE_DIR}/googletest/CMakeLists.txt")
add_subdirectory(googletest)
add_executable(test_kepler test_kepler.cc)
add_executable(test_fft test_fft.cc)
target_link_libraries(test_kepler part gtest_main gtest pthread )
target_link_libraries(test_fft part gtest_main gtest pthread ${FFTW_LIBRARIES})
add_custom_target(test ./test_kepler && ./test_fft DEPENDS part test_kepler test_fft)
else()
message("
*********************************************
google tests is missing.
Did you forget `git submodule update --init`
*********************************************
")
message(FATAL_ERROR "exit")
endif()
################################################################
# Doxygen
################################################################
find_package(Doxygen)
if (DOXYGEN_FOUND)
# to set other options, read: https://cmake.org/cmake/help/v3.9/module/FindDoxygen.html
set(DOXYGEN_EXCLUDE_PATTERNS */googletest/*)
doxygen_add_docs(
doxygen
${PROJECT_SOURCE_DIR}
COMMENT "Generate html pages"
)
add_custom_target(doc DEPENDS doxygen)
endif(DOXYGEN_FOUND)
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)
: xmin(xmin), xmax(xmax) {
Vector d = xmax - xmin;
for (UInt i = 0; i < Vector::dim; ++i)
if (d[i] < 0) {
std::cout << "XMax and XMin do not form a domain range" << std::endl;
std::exit(1);
}
}
/* -------------------------------------------------------------------------- */
void ComputeBoundary::compute(System& system) {
for (auto& par : system) {
auto& x = par.getPosition();
auto& v = par.getVelocity();
Vector d1 = xmax - x, d2 = x - xmin;
for (UInt i = 0; i < Vector::dim; ++i) {
if (d1[i] < 0 or d2[i] < 0)
v[i] *= -1;
}
}
}
/* -------------------------------------------------------------------------- */
#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) { this->penalty = penalty; }
/* -------------------------------------------------------------------------- */
void ComputeContact::compute(System& system) {
auto penalty_contact = [this](Particle& par1, Particle& par2) {
auto& ppb1 = static_cast<PingPongBall&>(par1);
auto& ppb2 = static_cast<PingPongBall&>(par2);
auto v_r = ppb2.getPosition() - ppb1.getPosition();
auto r = std::sqrt(v_r.squaredNorm());
v_r /= r;
auto penetration = ppb1.getRadius() + ppb2.getRadius() - r;
Vector force;
if (penetration > 0) {
force = v_r;
force *= penalty * penetration;
ppb1.getForce() -= force;
ppb2.getForce() += force;
}
};
applyOnPairs(penalty_contact, system);
}
#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 compute(System& system) 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::compute(System& system) {
// define a lambda function for newtonian gravity
auto newton_gravity = [this](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;
};
// apply newtonian gravity on all pairs of particles
applyOnPairs(newton_gravity, system);
}
#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 compute(System& system) 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>
/* -------------------------------------------------------------------------- */
#ifndef __COMPUTE_INTERACTION__HH__
#define __COMPUTE_INTERACTION__HH__
/* -------------------------------------------------------------------------- */
#include "compute.hh"
//! Base class for interaction computation
class ComputeInteraction : public Compute {
public:
//! Apply a functor on all particle pairs
template <typename Functor>
void applyOnPairs(Functor func, System& system);
};
/* -------------------------------------------------------------------------- */
/* Template implementation */
/* -------------------------------------------------------------------------- */
template <typename Functor>
void ComputeInteraction::applyOnPairs(Functor func, System &system) {
UInt size = system.getNbParticles();
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);
func(par1, par2);
}
}
}
/* -------------------------------------------------------------------------- */
#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) {}
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