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 Kommentare

moonman
moonman am 14 Nov. 2011
Verschoben: DGM am 26 Feb. 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 am 14 Nov. 2011
Verschoben: DGM am 26 Feb. 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 am 14 Nov. 2011
Verschoben: DGM am 26 Feb. 2023
can somebody guide me
Fangjun Jiang
Fangjun Jiang am 15 Nov. 2011
Bearbeitet: Walter Roberson am 5 Mär. 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 am 25 Mär. 2016
Bearbeitet: Ahmed ElTahan am 25 Mär. 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 am 5 Mai 2025
@Fangjun Jiang You are annoying

Melden Sie sich an, um zu kommentieren.

Antworten (5)

Honglei Chen
Honglei Chen am 14 Nov. 2011

2 Stimmen

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 Kommentare

moonman
moonman am 14 Nov. 2011
Verschoben: DGM am 26 Feb. 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 am 14 Nov. 2011
Verschoben: DGM am 26 Feb. 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 am 14 Nov. 2011
Verschoben: DGM am 26 Feb. 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 am 15 Nov. 2011
Verschoben: DGM am 26 Feb. 2023
Any more help will be appreciated Can anybody else confirm me that my code is correct
Walter Roberson
Walter Roberson am 15 Nov. 2011
Verschoben: DGM am 26 Feb. 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 am 16 Nov. 2011
Verschoben: DGM am 26 Feb. 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 am 9 Mai 2020
Verschoben: DGM am 26 Feb. 2023
output of difference equation should be discrete but the output waveform is coming continuous
Aditya Kumar Singh
Aditya Kumar Singh am 17 Nov. 2020
Verschoben: DGM am 26 Feb. 2023
use stem instead of plot... you should get the correct waveform after that

Melden Sie sich an, um zu kommentieren.

Mohazam Awan
Mohazam Awan am 10 Okt. 2017

2 Stimmen

%%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 Kommentare

Divya K M
Divya K M am 9 Mai 2020

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

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

Melden Sie sich an, um zu kommentieren.

Fangjun Jiang
Fangjun Jiang am 14 Nov. 2011

0 Stimmen

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 Kommentar

moonman
moonman am 14 Nov. 2011
Verschoben: DGM am 26 Feb. 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

Melden Sie sich an, um zu kommentieren.

BHOOMIKA MS
BHOOMIKA MS am 4 Dez. 2024

0 Stimmen

% 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 am 4 Dez. 2024

0 Stimmen

% 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');

Gefragt:

am 14 Nov. 2011

Kommentiert:

am 5 Mai 2025

Community Treasure Hunt

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

Start Hunting!

Translated by