%% programa solución de ecuacion de onda 01

clearvars;
close all;
clc;

%este programa grafica y anima la solucion de la ecuacion de onda
%mediante su solucion usando serie de fourier


%% datos inciales
L = 1;          %longitud de la cuerda
A = 1/4;        %amplitud
c = 4;          %veocidad de propagacion de ondas
x = 0:L/100:L;  %vector de distancia
Nx = length(x); %longitud de vector x

%% ciclo principal de animacion

N = 100;          %numero de terminos de la sumatoria

p(1:Nx) = 0;   %incializacion de presion

for t = 0:0.01:2
    
    p(1:Nx) = 0;   %incializacion de presion
    
    for n = 1:N
        p(1:Nx) = p(1:Nx) + (8*A)/(n^2*pi^2)*sin(n*pi/2)*cos(n*pi*c*t/L)*sin(n*pi*x/L);
    end
    
    figure(1)
    plot(x,p,'b','LineWidth',1.5)
    title('onda cuerda');
    xlabel('x (m)');
    ylabel('p(x,t) (m)')
    legend(['t = ',num2str(t),' (s)']);
    axis([0, L, -1.1*A, 1.1*A])
    grid on;
    box on;
    
    n= 1;
    f1 = (8*A)/(n^2*pi^2)*sin(n*pi/2)*sin(n*pi*x/L)*cos(n*pi*c*t/L);
    
    n= 2;
    f2 = (8*A)/(n^2*pi^2)*sin(n*pi/2)*sin(n*pi*x/L)*cos(n*pi*c*t/L);
    
    n= 3;
    f3 = (8*A)/(n^2*pi^2)*sin(n*pi/2)*sin(n*pi*x/L)*cos(n*pi*c*t/L);
    
    n= 4;
    f4 = (8*A)/(n^2*pi^2)*sin(n*pi/2)*sin(n*pi*x/L)*cos(n*pi*c*t/L);
    
    n= 5;
    f5 = (8*A)/(n^2*pi^2)*sin(n*pi/2)*sin(n*pi*x/L)*cos(n*pi*c*t/L);

    n= 6;
    f6 = (8*A)/(n^2*pi^2)*sin(n*pi/2)*sin(n*pi*x/L)*cos(n*pi*c*t/L);

    
    
    figure(2)
    plot(x,f1,x,f2,x,f3,x,f4,x,f5,x,f6)
    title('onda cuerda armónicos');
    xlabel('x (m)');
    ylabel('f(x,t) (m)')
    legend('f_1(t)','f_2(t)','f_3(t)','f_4(t)','f_5(t)','f_6(t)');
    axis([0, L, -1.1*A, 1.1*A])
    grid on;
    box on;
    
    pause(0.01);
    
end