how to do implement difference equation in matlab

Hi i am stuck with this question
Write a MATLAB program to simulate the following difference equation 8y[n] - 2y[n-1] - y[n-2] = x[n] + x[n-1] for an input, x[n] = 2n u[n] and initial conditions: y[-1] = 0 and y[0] = 1
(a) Find values of x[n], the input signal and y[n], the output signal and plot these signals over the range, -1 = n = 10.
The book has told to user filter command or filtic
my code is down kindly guide me about initial conditions

6 commentaires

moonman
moonman le 14 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
I have done in this way
n=0:10;
a = [8 -2 -1]; % left hand side of difference equation
b = [1 1 0]; % right hand side of difference equation
x=2.^n;
k=filter(b,a,x)
plot(n,k)
Is it correct How to cater for initial conditions y[-1]=0 and y[0]=1
moonman
moonman le 14 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
I have to do
Find values of x[n], the input signal and y[n], the output signal and plot these signals over the range, -1 = n = 10.
moonman
moonman le 14 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
can somebody guide me
Fangjun Jiang
Fangjun Jiang le 15 Nov 2011
Modifié(e) : Walter Roberson le 5 Mar 2021
My advice:
2. Use proper punctuation mark, e.g. comma, period, question mark.
3. Read your own question again after posting, e.g. what is x[n] = 2n u[n]?
4. Update your question with comments, not answers.
Ahmed ElTahan
Ahmed ElTahan le 25 Mar 2016
Modifié(e) : Ahmed ElTahan le 25 Mar 2016
Here is a function I have built to calculate it with added example. https://www.mathworks.com/matlabcentral/fileexchange/56142-output-estimation-difference-equation
Micah
Micah le 5 Mai 2025
@Fangjun Jiang You are annoying

Connectez-vous pour commenter.

Réponses (5)

I think your b is incorrect, it should be [1 1] instead. To accommodate for the initial value of y and x, you need to translate them into the corresponding filter state. filter command is an implementation of direct form II transpose, so you can use filtic to convert y and x to the state.
Here is an example, where n runs from 1 to 10. Based on you example, x[0] is 1.
n = 1:10;
a = [8 -2 -1];
b = [1 1];
yi = [1 0];
xi = 1;
zi = filtic(b,a,yi,xi)
y = filter(b,a,2.^n,zi)
BTW, I doubt if the input is really 2.^n as this becomes unbounded very quickly. Are you sure it's not 2.^(-n)?
HTH

8 commentaires

moonman
moonman le 14 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
Thanks Honglei this is the question of book. I hope it will be correct.U can read it again
*Write a MATLAB program to simulate the following difference equation: 8y[n] - 2y[n-1] - y[n-2] = x[n] + x[n-1] for an input, x[n] = 2n u[n] and initial conditions: y[-1] = 0 and y[0] = 1. (Hint : Try the commands ‘filter’ or ‘filtic’). a. Find the values of x[n], the input signal and y[n], the output signal and plot these signals over the range, -1 = n = 10. (All plots must be properly labeled). * So can u help me for this
moonman
moonman le 14 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
Is this solution Ok
n = -1:10;
a = [8 -2 -1]; % left hand side of difference equation
b = [1 1 0]; % right hand side of difference equation
yi = [1 0];
xi = 0;
x=2*n;
zi = filtic(b,a,yi,xi)
y = filter(b,a,2*n,zi)
subplot(2,1,1)
plot(n,y)
xlabel('--------n----------->')
ylabel('-------Values------->')
title(' Plot of output signal y over the range from n=-1 to n=10')
subplot(2,1,2)
plot(n,x)
xlabel('--------n----------->')
ylabel('-------Values------->')
title(' Plot of input signal x over the range from n=-1 to n=10')
Honglei Chen
Honglei Chen le 14 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
You are on the right track in general, but there are two things I want to point out:
  1. If your x is 2*n, then x(0) is 0.
  2. You already have y[0], y[-1], x[0], x[-1], so you don't need to compute them again.
Basically you should update your code to use
n = 1:10
and
xi = 0
This way, the y you get is y(1) through y(10). You can then concatenate y(0) and y(-1) to y to form the total vector. Similar things can be done for x too.
HTH
moonman
moonman le 15 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
Any more help will be appreciated Can anybody else confirm me that my code is correct
Walter Roberson
Walter Roberson le 15 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
You are working on a homework question. You do the best you can and if *you* cannot find an errors in your work, then you submit your answer and take your chances.
Have you considered working the results out manually and comparing them to the computed results?
Honglei Chen
Honglei Chen le 16 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
I think what Walter suggests is a great way to proceed. Or you can do a program, as Fangjun suggested to output the result using a program and see if the two approaches agree. Debug and verification is just as important as programming.
BTW I think I mention in my post that you have issue with your coefficient and that is not fixed yet.
Divya K M
Divya K M le 9 Mai 2020
Déplacé(e) : DGM le 26 Fév 2023
output of difference equation should be discrete but the output waveform is coming continuous
Aditya Kumar Singh
Aditya Kumar Singh le 17 Nov 2020
Déplacé(e) : DGM le 26 Fév 2023
use stem instead of plot... you should get the correct waveform after that

Connectez-vous pour commenter.

%%DSP LAb Task 4
% Difference equation implementation in matlab
%
clc
clear all
close all
% using filter function
n=-5:1:10;
index=find(n==0);
x=zeros(1,length(n));
x(index)=1;
subplot(2,2,1)
stem(n,x)
grid on
axis tight
b=[1 0];
a=[1 -2];
y=filter(b,a,x);
subplot(2,2,2)
stem(n,y,'filled','r')
grid on
axis tight
% Now without filter function
y1=zeros(1,length(n));
for i=1:length(n)
if(n(i)<0)
y1(i)=0;
end
if (n(i)>=0)
y1(i)=2*y1(i-1)+x(i);
end
end
subplot(2,2,3)
stem(n,y1,'filled','k')
grid on
axis tight

2 commentaires

Is this code for given question if so means can I explain

Daniel
Daniel le 27 Mar 2024
Why do you use negative samples in the n vector? I am trying to use only positve samples and it doesnt work

Connectez-vous pour commenter.

Fangjun Jiang
Fangjun Jiang le 14 Nov 2011

0 votes

It might be a filter. But I thought all the assignment was asking you to do is to write a for-loop to generate the y series data based on the equation and the initial conditions.

1 commentaire

moonman
moonman le 14 Nov 2011
Déplacé(e) : DGM le 26 Fév 2023
no book has told to use filter is my code correct?? How to cater for initial conditions y[-1]=0 and y[0]=1

Connectez-vous pour commenter.

BHOOMIKA MS
BHOOMIKA MS le 4 Déc 2024

0 votes

% Define the symbolic variable syms z n
% Define the Z-transform of the right-hand side 2^n rhs_z = 1 / (1 - 2*z^(-1));
% Define the left-hand side: Y(z) - 2z^(-1)Y(z) + z^(-2)Y(z) lhs_z = (1 - 2*z^(-1) + z^(-2)) * sym('Y(z)');
% Set up the equation for Y(z) eq = lhs_z == rhs_z;
% Solve for Y(z) Y_z = solve(eq, 'Y(z)');
% Simplify the expression for Y(z) Y_z_simplified = simplify(Y_z);
% Perform partial fraction decomposition on Y(z) Y_z_decomp = partfrac(Y_z_simplified, z);
% Display the decomposed Y(z) disp('Decomposed Y(z):'); disp(Y_z_decomp);
% Now take the inverse Z-transform for each term y_n = iztrans(Y_z_decomp);
% Display the time-domain solution y(n) disp('The time-domain solution y(n) is:'); disp(y_n);
% Create a numerical sequence for plotting % Define the range for n (e.g., n from 0 to 10) n_values = 0:10;
% Evaluate y(n) for each n using subs (substitute n into the expression) y_values = double(subs(y_n, n, n_values));
% Plot the solution figure;
% Plot y(n) stem(n_values, y_values, 'filled', 'LineWidth', 2); title('Time-domain solution y(n)'); xlabel('n'); ylabel('y(n)'); grid on;
% Add labels to the graph for clarity text(0, y_values(1), ['y(0) = ', num2str(y_values(1))], 'VerticalAlignment', 'bottom', 'HorizontalAlignment', 'right');
BHOOMIKA MS
BHOOMIKA MS le 4 Déc 2024

0 votes

% Define the symbolic variable syms z n
% Define the Z-transform of the right-hand side 2^n rhs_z = 1 / (1 - 2*z^(-1));
% Define the left-hand side: Y(z) - 2z^(-1)Y(z) + z^(-2)Y(z) lhs_z = (1 - 2*z^(-1) + z^(-2)) * sym('Y(z)');
% Set up the equation for Y(z) eq = lhs_z == rhs_z;
% Solve for Y(z) Y_z = solve(eq, 'Y(z)');
% Simplify the expression for Y(z) Y_z_simplified = simplify(Y_z);
% Perform partial fraction decomposition on Y(z) Y_z_decomp = partfrac(Y_z_simplified, z);
% Display the decomposed Y(z) disp('Decomposed Y(z):'); disp(Y_z_decomp);
% Now take the inverse Z-transform for each term y_n = iztrans(Y_z_decomp);
% Display the time-domain solution y(n) disp('The time-domain solution y(n) is:'); disp(y_n);
% Create a numerical sequence for plotting % Define the range for n (e.g., n from 0 to 10) n_values = 0:10;
% Evaluate y(n) for each n using subs (substitute n into the expression) y_values = double(subs(y_n, n, n_values));
% Plot the solution figure;
% Plot y(n) stem(n_values, y_values, 'filled', 'LineWidth', 2); title('Time-domain solution y(n)'); xlabel('n'); ylabel('y(n)'); grid on;
% Add labels to the graph for clarity text(0, y_values(1), ['y(0) = ', num2str(y_values(1))], 'VerticalAlignment', 'bottom', 'HorizontalAlignment', 'right');

Question posée :

le 14 Nov 2011

Commenté :

le 5 Mai 2025

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by