%% serie de fourier 01

%        { 1       para -2 <= t <  1
% f(t) = {
%        { 0       para  1 <= t <= 2
       
clearvars;
close all;
clc;

%% definiciones
Tmin =  0;
Tmax =  20;
dt = 0.005;

%%tiempo
t = Tmin:dt:Tmax;    
Nt = length(t);


%%función

ft = zeros(1,Nt);

%% serie de fourier

%termino a0
a0 = 3/2;
ft(1:Nt) = a0/2;

figure(1)
plot(t,ft,'b','Linewidth',1);
title('serie de fourier');
xlabel('t (s)');
ylabel('f_1(t) (Pa)');
legend('f_1(t) n = 0');
axis([Tmin,Tmax,-1.1,1.1]);
grid on;
box on;

%sumatoria

N = 200;
nn = 1:N;
an = zeros(N,1);
bn = zeros(N,1);

for n = 1:N
    %coeficientes
    an(n) = (1/(n*pi))*(sin(n*pi/2)-sin(n*pi));
    bn(n) = (1/(n*pi))*(cos(n*pi)-cos(n*pi/2));
    ft = ft + an(n)*cos(n*pi*t/2) + bn(n)*sin(n*pi*t/2);
    
    figure(1)
    plot(t,ft,'b','Linewidth',1);
    title('serie de fourier');
    xlabel('t (s)');
    ylabel('f(t) (Pa)');
    legend(['f(t) n = ',num2str(n)]);
    axis([Tmin,Tmax,-0.5,1.5]);
    grid on;
    box on;
    
    pause(0.01);
    
end

%% figura final
figure(2)
plot(t,ft,'b','Linewidth',1);
title('serie de fourier');
xlabel('t (s)');
ylabel('f(t) (Pa)');
legend(['f(t) n = ',num2str(n)]);
axis([Tmin,Tmax,-0.5,1.5]);
grid on;
box on;

%% figura coeficientes
figure(3)

subplot(2,1,1)
stem(nn,an,'b','Linewidth',1);
title('coeficientes a_n');
xlabel('n');
ylabel('a_n');
legend('a_n');
grid on;
box on;

subplot(2,1,2)
stem(nn,bn,'b','Linewidth',1);
title('coeficientes b_n');
xlabel('n');
ylabel('b_n');
legend('b_n');
grid on;
box on;


%% solucion a ecuación diferencial

C = 4 - 3/pi;

for n = 1:N
    C = C -( -4/(n*pi*(4*pi*n^2))*sin(n*pi/2) +  8*n/(n*pi*(4*pi*n^2))*(cos(n*pi)-cos(n*pi/2)) );
end

%solucion complementaria
yc = C*exp(-pi*t/4);

figure(4)
plot(t,yc,'b','Linewidth',1);
title('solucion complementaria');
xlabel('t (s)');
ylabel('y_c(t) (Pa)');
legend('y_c(t)');
axis([Tmin,Tmax,-0.5,4.5]);
grid on;
box on;

%% solucion particular

yp = 3/pi;

for n = 1:N
    %coeficientes
    an(n) = 4/(n*pi)/(4*pi*n^2)*sin(n*pi/2);
    bn(n) = 4/(n*pi)/(4*pi*n^2)*(cos(n*pi) - cos(n*pi/2));
    yp = yp + an(n)*( 2*n*sin(n*pi*t/2) + cos(n*pi*t/2) ) + bn(n)*(-2*n*cos(n*pi*t/2) + sin(n*pi*t/2) ); 
    
    figure(5)
    plot(t,yp,'b','Linewidth',1);
    title('solucion particular');
    xlabel('t (s)');
    ylabel('y_p(t) (Pa)');
    legend(['y_p(t) n = ',num2str(n)]);
    axis([Tmin,Tmax,0,1.4]);
    grid on;
    box on;
    
    pause(0.01);
    
end

%% figura compuesta

figure(6)
plot(t,yp,'b','Linewidth',1);
hold on;
plot(t,ft,'k','Linewidth',1);
title('solucion particular');
xlabel('t (s)');
ylabel('y_p(t) (Pa)');
legend(['y_p(t) n = ',num2str(n)]);
axis([Tmin,Tmax,-0.5,1.5]);
grid on;
box on;

%% solucion completa

y = yp + yc;

figure(7)
plot(t,y,'b','Linewidth',1);
hold on;
plot(t,ft,'k','Linewidth',1);
title('solucion completa');
xlabel('t (s)');
ylabel('y(t) (Pa)');
legend(['y(t) n = ',num2str(n)]);
axis([Tmin,Tmax,-0.5,5]);
grid on;
box on;