close all 
clc 
clear all 

% Paramètres physiques
c = 1.0;         % vitesse d'advection

% Paramètres numériques
L = 10.0;         % longueur du domaine
Nx = 500;        % nombre de points spatiaux
dx = L / (Nx-1); % pas spatial
x = linspace(0, L, Nx);

T = 10;         % temps final
CFL= 0.1;
dt = CFL*dx/c; % pas de temps
Nt = round(T/dt);

% Condition initiale : pic gaussien
u = (x<2).*(x>1);
u_new = zeros(size(u));

% Coefficients du schéma
lambda  = c * dt / dx;

% Boucle en temps
for n = 1:Nt
    for i = 2:Nx-1
        % Schéma upwind pour l'advection et centré pour la diffusion
        u_new(i) = u(i) ...
                   - lambda/2 * (u(i+1) - u(i-1)) ...
                   + lambda^2/2 * (u(i+1) - 2*u(i) + u(i-1));
    end

    % Conditions aux limites (Neumann homogènes ici)
    u_new(1) = u_new(2);
    u_new(end) = u_new(end-1);

    % Mise à jour
    u = u_new;

    % Affichage toutes les 1 itérations
    if mod(n,10)==0
        plot(x, u, 'b-', 'LineWidth', 2);
        %axis([0 1 0 1]);
        title(['t = ' num2str(n*dt)]);
        drawnow;
    end
end
