Resposta de um sistema discreto (PDS)

Esta é uma daquelas postagem somente para os iniciados, mas se você tiver coragem suficiente siga em frente!

***************

Seja um sistema discreto definido por:
$$y[n]  + 0.1y[n-1]  - 0.41y[n-2]  - 0.105y[n-3] = x[n] - 0.5x[n-1]$$
com $y[-1]=y[-2]=y[-3]=0$ e $x[n] = \cos(\pi k /6) u[n]$. Calcule a saída $y[n]$.

Solução. Sendo o sistema estável, sabemos que para $n$ grande a saída será dada por
$$y[n] = A \cos(\pi n /6 + \phi)$$
Sabemos também que durante transitório as raízes da equação característica são importantes. A equação característica é $$\lambda^3 + 0.1\lambda^2 - 0.41\lambda - 0.105 = 0$$ cujas raízes (podemos usar o comando roots do Scilab) são $\lambda_1 = 0.7$, $\lambda_2 = -0.5$ e $\lambda_3 = -0.3$. Logo: $$y[n] = C_1 (0.7)^n + C_2 (-0.5)^n + C_3(-0.3)^n + A \cos(\pi n /6 + \phi)$$
Determinação das constantes $C_1, C_2, C_3, A$ e $\phi$. Podemos calcular recursivamente os primeiros valores de $y[n]$, são eles: $[1.0,  1.2660254,  1.2164102,  0.7524294, 0.0564179]$. Logo, podemos montar o seguinte sistema de equações:
$$C_1   +  C_2  +  C_3  + A\cos(\phi)  = 1.0$$
$$0.7C_1   - 0.5C_2   - 0.3C_3  + A\cos(\pi/6 + \phi) = 1.2660254$$
$$0.49C_1  +  0.25C_2 +   0.09C_3 +   A\cos(\pi/3 + \phi) = 1.2164102$$
$$0.343C_1 - 0.125C_2 - 0.027C_3  + A\cos(\pi/2 + \phi)  = 0.7524294$$
$$0.2401C_1 + 0.0625C_2 + 0.0081C_3 + A\cos(4\pi/6 + \phi)  = 0.0564179$$

Expandindo o $\cos(a+b)$:
$$C_1   +  C_2  +  C_3  + A\cos(\phi)  = 1.0$$
$$0.7C_1   - 0.5C_2   - 0.3C_3  + A\sqrt{3}\cos(\phi)/2 - A\sin(\phi)/2 = 1.2660254$$
$$0.49C_1  +  0.25C_2 +   0.09C_3 +   A\cos(\phi)/2 - A\sqrt{3}\sin(\phi)/2 = 1.2164102$$
$$0.343C_1 - 0.125C_2 - 0.027C_3  - A\sin(\phi)  = 0.7524294$$
$$0.2401C_1 + 0.0625C_2 + 0.0081C_3 - A\cos(\phi)/2 - A\sqrt{3}\sin(\phi)/2  = 0.0564179$$

Fazendo $D = A\cos(\phi)$ e $E = A\sin(\phi)$, obtemos:
$$C_1   +  C_2  +  C_3  + D  = 1.0$$
$$0.7C_1   - 0.5C_2   - 0.3C_3  + D\sqrt{3}/2 - E/2 = 1.2660254$$
$$0.49C_1  +  0.25C_2 +   0.09C_3 +   D/2 - E\sqrt{3}/2 = 1.2164102$$
$$0.343C_1 - 0.125C_2 - 0.027C_3  - E  = 0.7524294$$
$$0.2401C_1 + 0.0625C_2 + 0.0081C_3 - D/2 - E\sqrt{3}/2  = 0.0564179$$

Cuja solução é: $C_1  = -0.2930939$, $C_2 = 0$, $C_3 = 0.0651971$, $D = 1.2278968$ e $E = -0.8547209$, logo $A\cos(\phi) = 1.2278968$ e $E\sin(\phi) = -0.8547209$. Finalmente: $A = -1.4960877$ e $\phi = 145.16^o$.

Finalmente temos:
$$y[n] = -0.2931(0.7)^n + 0.0652(-0.3)^n - 1.4961\cos(\pi n/6 + 145,16^o)$$

*** Código Scilab para gerar o gráfico e auxiliar nos cálculos:

close; clc;
k=0:29;
x=cos(%pi*k/6);
delta = 0*x; delta(1)=1;
d=convol([1, 0.8, 0.15],[1, -0.7]);
n=[1, 0.5];
y=filter(n,d,x);
h=filter(n,d,delta);
subplot(2,2,1); bar(x,0.2); 
title('Entrada x[n] = cos(pi*k/6)')
subplot(2,2,2); bar(h,0.2);
title('Resposta ao impulso do sistema h[n]');
subplot(2,1,2); bar(y,0.2);
title('Saída y[n]');

///////////  Cálculos:
m=[(0.7)^0, (-0.5)^0, (-0.3)^0 1, 0
(0.7)^1, (-0.5)^1, (-0.3)^1 sqrt(3)/2, -0.5
(0.7)^2, (-0.5)^2, (-0.3)^2 0.5,  -sqrt(3)/2
(0.7)^3, (-0.5)^3, (-0.3)^3 0, -1
(0.7)^4, (-0.5)^4, (-0.3)^4 -0.5 -sqrt(3)/2];
v = y(1:5)';
cs = inv(m)*v;
A = -sqrt(cs(4)^2 + cs(5)^2);
phi = atan(cs(5)/A, cs(4)/A);

k=0:29;
yinf = cs(1)*(0.7).^k + cs(2)*(-0.5).^k + ...
cs(3)*(-0.3).^k + A*cos(%pi*k/6 + phi);
plot(k+1,yinf,'-o');