Commit 80b0524f authored by Matthieu Boileau's avatar Matthieu Boileau

Add tp6/exercice2

parent fb627c83
Pipeline #3997 passed with stages
in 2 minutes and 6 seconds
cmake_minimum_required(VERSION 3.3)
project(Maillage)
# Recherche des dépendances Python
message(STATUS "Find Python dependencies using CMake ${CMAKE_VERSION}")
if(CMAKE_VERSION VERSION_LESS 3.12) # Pour CMake < 3.12
find_package(PythonInterp)
find_package(PythonLibs)
set(Python_NumPy_INCLUDE_DIRS "/opt/conda/lib/python3.7/site-packages/numpy/core/include")
set(Python_LIBRARIES ${PYTHON_LIBRARIES})
set(Python_INCLUDE_DIRS ${PYTHON_INCLUDE_DIRS} ${Python_NumPy_INCLUDE_DIRS})
else()
if(APPLE) # Pour MacOS seulement
set(CMAKE_FIND_FRAMEWORK NEVER)
endif()
# Nécessite CMake >= 3.12
find_package(Python3 COMPONENTS Interpreter NumPy required)
set(Python_NumPy_INCLUDE_DIRS ${Python3_NumPy_INCLUDE_DIRS})
set(Python_LIBRARIES ${Python3_LIBRARIES})
set(Python_INCLUDE_DIRS ${Python3_INCLUDE_DIRS} ${Python_NumPy_INCLUDE_DIRS})
endif()
message(" Python_NumPy_INCLUDE_DIRS: ${Python_NumPy_INCLUDE_DIRS}")
message(" Python_INCLUDE_DIRS : ${Python_INCLUDE_DIRS}")
message(" Python_LIBRARIES : ${Python_LIBRARIES}")
# On ajoute les en-têtes des biliothèques internes
include_directories(ext ${Python_INCLUDE_DIRS})
# On ajoute les en-têtes du projet
include_directories(include)
# On construit la liste des fichiers sources à partir du contenu de src/
file(GLOB_RECURSE SOURCE_FILES src/*.cpp)
# On indique les options de compilations
set(CMAKE_CXX_FLAGS "-Wall -std=c++11")
# On crée la bibliothèque libgraphe.a
add_library(graphe ${SOURCE_FILES})
# On génère l'exécutable graphe.e à partir du fichier src/main.cpp
add_executable(graphe.e src/main.cpp)
# On lie l'exécutable avec la bibliothèque libgraphe.a
target_link_libraries(graphe.e graphe ${Python_LIBRARIES})
# On ajoute les test unitaires
enable_testing()
# On construit la liste des fichiers de test
file(GLOB_RECURSE TEST_FILES test/*.cpp)
foreach(TEST_FILE ${TEST_FILES})
# On récupère le nom sans l'extension
get_filename_component(TEST_NAME ${TEST_FILE} NAME_WE)
# On crée l'exécutable test_toto.e à partir du fichier test_toto.cpp
add_executable(${TEST_NAME}.e ${TEST_FILE})
target_link_libraries(${TEST_NAME}.e graphe ${Python_LIBRARIES})
# On ajoute le test "test_toto" qui utilise l'exécutable "test_toto.e"
add_test(${TEST_NAME} ${TEST_NAME}.e)
endforeach()
CC = g++
CC_FLAGS = -std=c++11 -Wall
EXEC_NAME = main
OBJ_FILES = main.o modele.o modele_sir.o modele_seir.o solveur_temps.o euler_ex.o rk2.o graphe.o
all: $(EXEC_NAME)
clean:
rm -f $(EXEC_NAME) $(OBJ_FILES)
$(EXEC_NAME) : $(OBJ_FILES)
$(CC) $(CC_FLAGS) -o $(EXEC_NAME) $(OBJ_FILES)
%.o: %.cpp %.hpp
$(CC) $(CC_FLAGS) -o $@ -c $<
main.o : main.cpp modele_graphe.hpp
This diff is collapsed.
This diff is collapsed.
#ifndef EULER_EX_HPP
#define EULER_EX_HPP
#include "solveur_temps.hpp"
class EulerEx: public SolveurTemps {
public:
EulerEx(): SolveurTemps() {}
EulerEx(double dt): SolveurTemps(dt) {}
// Calcule la solution pour un pas de temps supplémentaire
void applique();
};
#endif
#ifndef GRAPHE_HPP
#define GRAPHE_HPP
#include <cmath>
#include <iostream>
#include <list>
#include <vector>
using namespace std;
class voisin {
int number;
double poids;
public:
voisin(){
poids =0.0;
number = 0;
};
voisin(const voisin & n){
poids = n.poids;
number = n.number;
};
~voisin (){
};
int num(){
return number;
}
double p(){
return poids;
}
void set_num(int a){
number = a;
}
void set_p(double a){
poids = a;
}
};
template< class T > class noeud{
T val;
int number;
vector<voisin> voisins;
public:
noeud(){
val = T();
number = 0;
};
noeud(const noeud & n){
val = n.val;
number = n.number;
voisins = n.voisins;
};
~noeud (){
voisins.clear();
};
void ajout_val(T x){
val = x;
}
void ajout_num(int n){
number = n;
}
vector<voisin> get_voisins() {
return voisins;
}
T get_val() {
return val;
}
int num(){
return number;
}
void ajout_voisin(int nbnoeud, double p){
voisin v;
v.set_num(nbnoeud);
v.set_p(p);
voisins.push_back(v);
}
void supprime_voisin(int nbnoeud){
vector<voisin>::iterator pos;
for(pos = voisins.begin(); pos != voisins.end(); pos++){
if ( pos->num() == nbnoeud ) break;
}
if (pos != voisins.end() ) {
voisins.erase(pos);
}
}
noeud & operator = (const noeud & n){
if(this != &n){
voisins.clear();
val = n.val;
voisins = n.voisins;
number = n.number;
}
return *this;
}
bool operator == (const noeud & n){
bool res = false;
res = (val == n.val) && (n.number == number);
return res;
}
template< class U>
friend ostream & operator<<(ostream & out, noeud<U> & x);
};
template<class T>class graphe{
vector< noeud<T> * > noeuds;
int nbnodes;
public:
graphe(){
nbnodes = 0;
};
graphe(const graphe & a){
noeuds = a.noeuds;
nbnodes = a.nbnodes;
};
~graphe (){
noeuds.clear();
};
int get_nbnodes() {
return nbnodes;
}
noeud<T> * get_noeud(int i) {
return noeuds[i];
}
void ajout_noeud(noeud<T> * n){
n->ajout_num(nbnodes);
noeuds.push_back(n);
nbnodes++;
}
void supprime_noeud(noeud<T> & n){
typename vector<noeud<T>* >::iterator pos;
for(int i = 0; i< noeuds.size() ; i++){
noeuds[i]->supprime_voisin(n.num());
}
pos = find (noeuds.begin(), noeuds.end(), &n);
if (pos != noeuds.end() ) {
noeuds.erase(pos);
}
nbnodes--;
}
template< class U>
friend ostream & operator<<(ostream & out, graphe<U> & x);
};
template <class T> ostream & operator<<(ostream & out, noeud<T> & x){
vector<voisin>::iterator ii;
int nvoisins=0;
out<<" noeuf numero: "<<x.number<<" de valeur "<<x.val<<endl;
if(x.voisins.size() != 0){
out<<">>>>> liste des voisins:"<<endl;
for(ii = x.voisins.begin(); ii != x.voisins.end(); ii++){
out<<" voisins n"<<nvoisins<<" :"<<ii->num()<<endl;
nvoisins++;
}
}
return out;
}
template <class T> ostream & operator<<(ostream & out, graphe<T> & x){
int nnoeud=0;
int i=0;
if(x.noeuds.size() != 0){
cout<< ">>>> list des noeuds" <<endl;
for(i = 0; i< x.noeuds.size() ; i++){
out<<" node n"<<x.noeuds[i]->num()<<" :"<<*x.noeuds[i]<<endl;
nnoeud++;
}
} else { cout<< ">>>> graphe nul" <<endl; }
return out;
}
#endif
#ifndef MODELE_HPP
#define MODELE_HPP
#include <iostream>
#include <vector>
class Modele {
protected:
// Nombre de variables du modèle
int m_n;
// Conditions initiales du modèle
std::vector<double> m_initial_value;
public:
// Renvoie la valeur de m_n
int get_n();
// Renvoie a valeur de m_initial_value
std::vector<double> initial_value();
// Copie le paramètre v dans m_initial_value
void setInitialValue(std::vector<double> v);
// Calcule le flux associé au vecteur en paramètre
virtual std::vector<double> flux(std::vector<double>) = 0;
};
#endif
#ifndef MODELE_GRAPHE_HPP
#define MODELE_GRAPHE_HPP
#include <vector>
#include "modele.hpp"
#include "graphe.hpp"
// Le paramètre T correspond au type de modèle localement utilisé sur chaque noeud du graphe
template<class T>
class ModeleGraphe : public Modele {
private:
// Graphe sur lequel le modèle d'épidémie sera appliqué
graphe<T> *m_graphe;
public:
ModeleGraphe();
ModeleGraphe(graphe<T> *g);
~ModeleGraphe();
// Calcule les conditions initiales du modèle en fonction des conditions initiales
// du modèle associé à chaque noeud du graphe.
void setInitialValue();
// Calcule le flux associé au vecteur en paramètre
virtual std::vector<double> flux(std::vector<double>);
};
template<class T>
ModeleGraphe<T>::ModeleGraphe() {
m_graphe = NULL;
m_n = 0;
}
template<class T>
ModeleGraphe<T>::ModeleGraphe(graphe<T> *g) {
m_graphe = g;
T t;
if (g != NULL) {
m_n = g->get_nbnodes() * t.get_n();
} else {
m_n = 0;
}
}
template<class T>
ModeleGraphe<T>::~ModeleGraphe() {
}
template<class T>
void ModeleGraphe<T>::setInitialValue() {
// création du vecteur contenant les conditions initiale :
// On concatène les vecteurs de conditions initiales de chaque noeud ;
// Si ils sont invalides, on met la valeur initiale à 0 par défaut.
vector<double> res;
for (int i=0; i<m_graphe->get_nbnodes(); i++) {
T t;
if (m_graphe->get_noeud(i)->get_val().initial_value().size() == t.get_n()) {
for (int j=0; j<t.get_n(); j++) {
res.push_back(m_graphe->get_noeud(i)->get_val().initial_value()[j]);
}
} else {
for (int j=0; j<t.get_n(); j++) {
res.push_back(0);
}
}
}
m_initial_value = res;
}
template<class T>
vector<double> ModeleGraphe<T>::flux(vector<double> v) {
vector<double> res(m_n);
T t;
// Traitement des termes de flux pour chaque noeud
vector<double> flux_local;
// Repère placé sur le premier élément correspondant au noeud considéré
vector<double>::const_iterator debut_v_local;
for (int i=0; i<m_graphe->get_nbnodes(); i++) {
debut_v_local = v.begin() + i * t.get_n();
// Flux associé au noeud i
flux_local = m_graphe->get_noeud(i)->get_val().flux(vector<double>(debut_v_local, debut_v_local + t.get_n()));
// Ajout du flux local dans le vecteur de flux global
for (int j=0; j<t.get_n(); j++) {
res[j + t.get_n() * i] += flux_local[j];
}
// Parcours des voisins du noeud
// Numéro du noeud voisin considéré
int n_voisin;
// Poids associé à ce voisin
double poids_voisin;
// Ajout des termes de couplage
for (int ind_voisin=0; ind_voisin<m_graphe->get_noeud(i)->get_voisins().size(); ind_voisin++) {
n_voisin = m_graphe->get_noeud(i)->get_voisins()[ind_voisin].num();
poids_voisin = m_graphe->get_noeud(i)->get_voisins()[ind_voisin].p();
for (int j=0; j<t.get_n(); j++) {
res[j + t.get_n() * n_voisin] += m_graphe->get_noeud(i)->get_val().g() * poids_voisin * v[j + i * t.get_n()];
}
}
}
return res;
}
#endif
#ifndef MODELE_SEIR_HPP
#define MODELE_SEIR_HPP
#include <vector>
#include "modele.hpp"
class ModeleSeir : public Modele {
private:
double m_beta;
double m_gamma;
double m_alpha;
double m_nu;
double m_mu;
double m_g;
public:
ModeleSeir();
ModeleSeir(double beta, double gamma, double alpha, double nu, double mu, double g);
~ModeleSeir() {}
double g();
std::vector<double> flux(std::vector<double>);
};
#endif
#ifndef MODELE_SIR_HPP
#define MODELE_SIR_HPP
#include <vector>
#include "modele.hpp"
class ModeleSir : public Modele {
private:
double m_beta;
double m_gamma;
double m_nu;
double m_mu;
double m_g;
public:
ModeleSir();
ModeleSir(double beta, double gamma, double nu, double mu, double g);
~ModeleSir() {}
double beta();
double gamma();
double nu();
double mu();
double g();
virtual std::vector<double> flux(std::vector<double>);
};
#endif
#ifndef RK2_HPP
#define RK2_HPP
#include "solveur_temps.hpp"
class Rk2: public SolveurTemps {
public:
Rk2(): SolveurTemps() {}
Rk2(double dt): SolveurTemps(dt) {}
void applique();
};
#endif
#ifndef SOLVEUR_TEMPS_HPP
#define SOLVEUR_TEMPS_HPP
#include <vector>
#include "modele.hpp"
#include "matplotlibcpp.h"
class SolveurTemps {
protected:
// Pas de temps
double m_dt;
// Vecteur contenant les vecteurs solution à chaque pas de temps
std::vector< std::vector<double> > m_variablesTemps;
// Modèle d'épidémie utilisé
Modele *m_modele;
public:
SolveurTemps();
SolveurTemps(double dt);
~SolveurTemps();
// Calcule la solution pour un pas de temps supplémentaire
virtual void applique() = 0;
// Etablit le modèle utilisé
void setModele(Modele * m);
// Calcule la valeur initiale du vecteur solution
void setInitialValue();
// Renvoie le temps correspondant au dernier vecteur de m_variablesTemps
double getCurrentTime();
// Renvoie la dernière solution calculée
std::vector<double> getSolutionCurrentTime();
void plotSolution(const int inoeud, const std::string&);
};
#endif
This diff is collapsed.
import matplotlib.pyplot as plt
import csv
filename = "paris.txt"
S = []
I = []
R = []
with open(filename, 'r') as csvfile:
plots = csv.reader(csvfile, delimiter=',')
for row in plots:
S.append(float(row[0]))
I.append(float(row[1]))
R.append(float(row[2]))
x = [0.01 * i for i in range(len(S))]
plt.plot(x[:8000], S[:8000], label='S')
plt.plot(x[:8000], I[:8000], label='I')
plt.plot(x[:8000], R[:8000], label='R')
plt.title(filename)
plt.legend()
plt.show()
#include "euler_ex.hpp"
using namespace std;
void EulerEx::applique() {
// Solution au temps précédent
vector<double> sol_n = getSolutionCurrentTime();
// Solution calculée
vector<double> sol_np1;
// Flux évalué au temps t précédent
vector<double> flux = m_modele->flux(sol_n);
// Calcul de la solution à t + dt
for (int i=0; i<sol_n.size(); i++) {
sol_np1.push_back(sol_n[i] + m_dt * flux[i]);
}
// Ajout de la solution au vecteur des solutions
m_variablesTemps.push_back(sol_np1);
}
This diff is collapsed.
#include "graphe.hpp"
#include <cmath>
#include <list>
#include <vector>
using namespace std;
#include <iostream>
#include <fstream>
#include "graphe.hpp"
#include "euler_ex.hpp"
#include "rk2.hpp"
#include "modele_sir.hpp"
#include "modele_graphe.hpp"
using namespace std;
int main() {
// Création du graphe de villes
graphe<ModeleSir> france;
noeud<ModeleSir> Paris;
france.ajout_noeud(&Paris); // 0
noeud<ModeleSir> Lille;
france.ajout_noeud(&Lille); //1
noeud<ModeleSir> Bordeaux;
france.ajout_noeud(&Bordeaux); //2
noeud<ModeleSir> Strasbourg;
france.ajout_noeud(&Strasbourg); //3
noeud<ModeleSir> Nantes;
france.ajout_noeud(&Nantes); //4
noeud<ModeleSir> Lyon;
france.ajout_noeud(&Lyon); //5
noeud<ModeleSir> Marseille;
france.ajout_noeud(&Marseille); //6
// Ordre des paramètres dans le constructeur du modèle :
// beta, gamma, nu, mu, g
double g = 0.1;
double beta = 0.7; //0.5
double gamma = 0.1;
ModeleSir modeleParis(beta, gamma, 0, 0, g);
ModeleSir modeleLille(beta, gamma, 0, 0, g);
ModeleSir modeleBordeaux(beta, gamma, 0, 0, g);
ModeleSir modeleStrasbourg(beta, gamma, 0, 0, g);
ModeleSir modeleNantes(beta, gamma, 0, 0, g);
ModeleSir modeleLyon(beta, gamma, 0, 0, g);
ModeleSir modeleMarseille(beta, gamma, 0, 0, g);
// Création du vecteur contenant les valeurs initiales
vector<double> v;
v.push_back(0.34);
v.push_back(0.1);
v.push_back(0);
modeleParis.setInitialValue(v);
v[0] = 0.05;
v[1] = 0;
modeleLille.setInitialValue(v);
v[0] = 0.05;
modeleBordeaux.setInitialValue(v);