Custom Euler's Method for Second Order ODE

4 views (last 30 days)
Hello, I am trying to develop a way to solve a specific differential equation using Euler's method. I have some code so far, but it is not giving the correct approximation and I cannot see where I made an error. The ODE represents a second order RLC circuit.
The ODE I am trying to approximate: v'' + 5v' + 1000v = u(t-1) - u(t-3) where u(t) is the unit step function.
The approximation is supposed to look something like this:
Instead, my code gives me something like this:
Here is my code:
% Euler's Method Code:
v(1,1) = 0; v(2,1) = 0; time(1)=0; %Initial conditions
h = .0001; %step size
N = 5/h; %Set up recursions
%Set up A matrix
A = [0,1;-5,-1000];
for k = 1:N
time(k+1) = time(k)+h; %concatenates new time position to time vector
v(:,k+1)=v(:,k) + h * A * v(:,k) + ...
h * [0; heaviside(time(k)-1)-heaviside(time(k)-3)];
end
figure('Name','Part 1 Question 11')
plot(time,v(1,:));
xlabel('time (s)')
ylabel('Capacitor Voltage')
xlim([0,5]);
I have a feeling it may be something to do with the driving function. The driving function is made up of the difference of the two unit step/ heaviside functions.
I have tried using changing the matrix inside the for loop from this:
[0; heaviside(time(k)-1)-heaviside(time(k)-3)]
To this:
[0; (time(k)>1)-(time(k)>3)]
But then the graphed solution is 0 for the duration of the graph. I've been staring at this for hours trying to figure out where my solution is messed up.
Thank you in advance for any help or suggestions.

Accepted Answer

James Tursa
James Tursa on 24 Apr 2020
Edited: James Tursa on 24 Apr 2020
Your derivative matrix is wrong. It should be this (the -5 and -1000 are switched):
A = [0,1;-1000,-5];
  2 Comments
Jacob McElwain
Jacob McElwain on 24 Apr 2020
Thank you so much! I can't believe I overlooked such a minor error. I fixed it and it gave me the correct solution.

Sign in to comment.

More Answers (0)

Categories

Find more on Numerical Integration and Differential Equations in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!