5 views (last 30 days)

How to get two lines in one graph for

y''-10y'-11y=0 with initial conditions y(0)=1,y'(0)=-1

and with the change in the first intial condition y(0)=1+epsilon ,y'(0) =-1

I run Heun's method but the solutions are slightly different though it should show the stability for the first line and unstability for the 2nd line.

John D'Errico
on 8 Oct 2020

Your question about how to get two lines in one graph is a simple one. Use the command "hold on" between plots. That is...

plot(1:5,rand(1,5),'r-o')

hold on

plot(1:5,rand(1,5),'b-x')

hold off

So now we have one plot, with two curves drawn on it.

Now look at the problem you are trying to solve. Though I won't do your homework for you, writing a Heun's method solver to produce several results plotted together, we can look at the problem and explain what should happen. It was:

y''-10y'-11y=0 with initial conditions y(0)=1,y'(0)=-1

What is happening there? Look at this curve, in the form of a symbolic solution.

syms y(t)

dy = diff(y,t);

ysol = dsolve(diff(y,t,2) - 10*dy - 11*y == 0,y(0) == 1,dy(0) == -1)

As you should see, this solution is indeed a nice, simple thing. In the analytical solution, we have a second term in there multiplying exp(11*t), but in the analytical solution for the unperturbed problem the multipler of exp(11*t) was exactly ZERO. The full analytical solution could have been written

ysol = exp(-t) + 0*exp(11*t)

given that set of initial conditions. That is theory. Theory is always so nice, right? But then practice gets in the way.

What happens now if we change the initial conditions just by a tiny amount? A tiny error in the nether-bits is enough to nudge the solver to give me a second non-zero term in there.

ysol2 = dsolve(diff(y,t,2) - 10*dy - 11*y == 0,y(0) == 1.0000000000001,dy(0) == -1)

vpa(expand(ysol2),20)

So what happens with this set of initial conditions is as the solver moves along, you need to remember you are using a numerical ODE solver. Is the solution produced by any numerical solver exact? No. In fact, you should know it has an error term. You should expect some tiny amount of slop in the solution. But now once we wander off the true colution, even by the tiniest amount, that second part comes into the problem, and it introduces itself with a vengeance. That second term gets HUGE, and it does so very rapidly. So as soon as any tiny amount of crap gets into the problem, the ODE solver can start to chase along the bad solution.

It will get worse if you start away from the stable part of the solution by just a tiny amount. Now the solver is already looking at something that will get huge quickly. It does not even need to introduce numerical errors into the problem, to then be massively amplified in later steps.

Opportunities for recent engineering grads.

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

Start Hunting!
## 2 Comments

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/607816-how-to-find-unstable-solution#comment_1044761

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/607816-how-to-find-unstable-solution#comment_1044761

## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/607816-how-to-find-unstable-solution#comment_1045146

⋮## Direct link to this comment

https://fr.mathworks.com/matlabcentral/answers/607816-how-to-find-unstable-solution#comment_1045146

Sign in to comment.