#include <iostream>
#include <vector>
#include <algorithm>
#include <cmath>
#include <time.h>
#include <ctime>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_randist.h>
#include <cstdio> // Pour popen
// Fonction pour simuler N_T (nombre de sauts)
unsigned int simulate_NT(gsl_rng* rng, double lambda_T) {
double sum_U = 0.0;
unsigned int k = 0;
while (sum_U <= lambda_T) {
double U = gsl_rng_uniform(rng); // Générer une variable uniforme [0, 1]
sum_U += U;
++k;
}
return k - 1; // Retourner N_T = k - 1
}
// Fonction pour simuler les instants d'arrivée et les sauts
void simulate_poisson_jumps(double T, double lambda, double mu_J, double sigma_J,
std::vector<double>& times, std::vector<double>& jumps)
{
gsl_rng *rng = gsl_rng_alloc(gsl_rng_mt19937); // Initialiser le générateur GSL
gsl_rng_set(rng, time(NULL)); // Graine pour les nombres aléatoires
// Simuler N_T
double lambda_T = lambda * T;
unsigned int N_T = simulate_NT(rng, lambda_T);
// Simuler les instants d'arrivée
std::vector<double> arrival_times(N_T);
for (unsigned int i = 0; i < N_T; ++i) {
arrival_times[i] = gsl_ran_flat(rng, 0, T);
}
std::sort(arrival_times.begin(), arrival_times.end());
// Simuler les sauts Y_i ~ N(mu_J, sigma_J^2)
for (unsigned int i = 0; i < N_T; ++i) {
double jump = gsl_ran_gaussian(rng, sigma_J) + mu_J;
times.push_back(arrival_times[i]);
jumps.push_back(jump);
}
gsl_rng_free(rng); // Libérer le générateur
}
// Fonction pour simuler la trajectoire d'un actif avec la dynamique Black-Scholes
avec sauts
void simulate_black_scholes_with_jumps(double T, double S0, double r, double sigma,
double lambda,
double mu_J, double sigma_J,
std::vector<double>& times,
std::vector<double>& values) {
gsl_rng *rng = gsl_rng_alloc(gsl_rng_mt19937); // Initialiser le générateur GSL
gsl_rng_set(rng, time(NULL)); // Graine pour les nombres aléatoires
// Simuler le mouvement Brownien
unsigned int steps = 10000; // Nombre de pas de temps
double dt = T / steps;
double sqrt_dt = sqrt(dt);
std::vector<double> W(steps + 1, 0.0);
for (unsigned int i = 1; i <= steps; ++i) {
W[i] = W[i - 1] + gsl_ran_gaussian(rng, sqrt_dt);
}
// Simuler les sauts
std::vector<double> jump_times, jump_magnitudes;
simulate_poisson_jumps(T, lambda, mu_J, sigma_J, jump_times, jump_magnitudes);
// Construire la trajectoire de S_t
times.push_back(0.0);
values.push_back(S0);
double current_S = S0;
double jump_sum = 0.0;
unsigned int jump_index = 0;
for (unsigned int i = 1; i <= steps; ++i) {
double t = i * dt;
// Ajouter les sauts jusqu'au temps t
while (jump_index < jump_times.size() && jump_times[jump_index] <= t) {
jump_sum += jump_magnitudes[jump_index]+1;
++jump_index;
}
// Calculer S_t
double drift = (r - 0.5 * sigma * sigma) * t;
double diffusion = sigma * W[i];
current_S = S0 * exp(drift + diffusion + jump_sum);
times.push_back(t);
values.push_back(current_S);
}
gsl_rng_free(rng); // Libérer le générateur
}
// Fonction pour afficher les données avec Gnuplot
void plot_with_gnuplot(const std::vector<double>& times, const std::vector<double>&
values, const char* title) {
FILE* gnuplot = popen("C:\\gnuplot\\bin\\[Link]", "w");
if (!gnuplot) {
std::cerr << "Erreur : Impossible de démarrer Gnuplot.\n";
return;
}
fprintf(gnuplot, "set title '%s'\n", title);
fprintf(gnuplot, "set xlabel 'Temps (t)'\n");
fprintf(gnuplot, "set ylabel 'S(t)'\n");
fprintf(gnuplot, "set grid\n");
fprintf(gnuplot, "plot '-' with lines title '%s'\n", title);
for (size_t i = 0; i < [Link](); ++i) {
fprintf(gnuplot, "%f %f\n", times[i], values[i]);
}
fprintf(gnuplot, "e\n");
fflush(gnuplot);
std::cout << "Appuyez sur Entrée pour fermer Gnuplot...\n";
std::[Link]();
pclose(gnuplot);
}
int main() {
double T = 1.0; // Intervalle de temps [0, T]
double S0 = 100.0; // Prix initial de l'actif
double r = 0.05; // Taux d'intérêt sans risque
double sigma = 0.2; // Volatilé
double lambda = 0.6; // Intensité du processus de Poisson
double mu_J = -0.1; // Moyenne des sauts
double sigma_J = 0.1; // Écart-type des sauts
std::vector<double> times, values;
// Simuler la trajectoire
simulate_black_scholes_with_jumps(T, S0, r, sigma, lambda, mu_J, sigma_J,
times, values);
// Tracer la trajectoire avec Gnuplot
plot_with_gnuplot(times, values, "Dynamique Black-Scholes avec sauts");
return 0;
}