% ============================================================
% Simulazione di onde viaggianti tumore-sistema immunitario
% Modello di reazione-diffusione senza infezione
% ============================================================

% Parametri del modello
Du = 0.01; Dz = 0.05;
p = 1; K = 1; zeta = 1; qz = 1; Sz = 0.5;

% Spazio e tempo
L = 100; N = 400; dx = 2*L/N; x = linspace(-L, L, N);
dt = 0.1; T = 1000; Nt = round(T/dt);

% Inizializzazione
u = zeros(N,1);
z = ones(N,1) * Sz/qz; % inizialmente tutto immunitario
u(x < 0) = 1;          % tumore inizialmente a sinistra

% Matrice di diffusione (Laplaciano discreto 1D)
e = ones(N,1);
Lapl = spdiags([e -2*e e], -1:1, N, N) / dx^2;
Lapl(1,:) = 0; Lapl(end,:) = 0; % condizioni di Neumann (flusso nullo)

% Tempi ai quali salvare i profili
tempi_salvati = [200, 400, 600];
u_salvati = zeros(N, length(tempi_salvati));
z_salvati = zeros(N, length(tempi_salvati));

% Evoluzione temporale
for n = 1:Nt
    % Aggiornamento delle equazioni di reazione-diffusione
    u = u + dt * (Du * Lapl * u + p*u.*(1 - u/K) - (zeta/K)*u.*z);
    z = z + dt * (Dz * Lapl * z - qz*z + Sz);

    % Salva i profili ai tempi desiderati
    t_attuale = round(n*dt);
    idx = find(tempi_salvati == t_attuale, 1);
    if ~isempty(idx)
        u_salvati(:,idx) = u;
        z_salvati(:,idx) = z;
    end

    % (Opzionale) aggiornamento grafico in tempo reale
    if mod(n,100)==0
        plot(x,u,'r',x,z,'b');
        legend('Tumore','Immunitario');
        title(['Tempo t = ' num2str(n*dt)]);
        xlabel('x'); ylabel('densità');
        drawnow;
    end
end

% ============================================================
% Grafico finale: famiglia di curve del tumore a t = 200,400,600
% ============================================================
figure;
hold on;
colori = lines(length(tempi_salvati));
for k = 1:length(tempi_salvati)
    plot(x, u_salvati(:,k), 'Color', colori(k,:), 'LineWidth', 1.5);
end
legend(arrayfun(@(t)['t = ' num2str(t)], tempi_salvati, 'UniformOutput', false));
xlabel('x');
ylabel('u(x,t)');
title('Evoluzione temporale del tumore');
grid on;
hold off;

% Salvataggio dell'immagine finale
saveas(gcf, 'famiglia_curve_tumore.png');

% ============================================================
% FINE DEL PROGRAMMA
% ============================================================
