solving ode equation using ode45

I have an ode which I need to solve with respect to the time variable.

I have 2 seperate scripts, one contains the function and the other contains the command.

my function is,

===================================

function dy_dt=func(t,y)
dy_dt=0.003445*((y-0.072)/y^3)^0.5;
end

===================================

In the other script I have,

===================================

clc, clear, close all

time=0:0.05:8.5;
[t,y]=ode45(@myfun1,time,0.073);
plot(t,y(:,1));

==================================

When the initial value for y=0.072 I get the following graph;

because y=0.072 leads to an erroneous graph I decided to give y=0.073 as the initial value for which I got the following graph;

This is also erroneous as this is not the result I expect. The graph I expect is attached below.

Can someone assist me on what errors I have done in my code? Thank you very much.

[darova]

Please attach original equation

[nadeeshani Haggalla]

[nadeeshani Haggalla]

y0 is the initial value which is 0.072

[darova]

Your code is correct. Check your initial conditions

[James Tursa]

The derivative expression y-y0 will evaluate to be exactly 0 initially and thus there is nothing to drive the solution away from the initial condition by even a tiny amount. An exactly horizontal line at y0 would be the expected result in this case. As a result, your differential equation doesn't really make any sense.

[nadeeshani Haggalla]

@James Tursa @darova Thank you for your response. Yes, which is why I changed the initial condition and ran the code again. Eventhough I'm getting a graph, it is not the result I expect. That's what I'm struggling with.

[James Tursa]

It you pick something other than exactly y0 for the initial condition, what does the y0 in the differential equation even mean anymore? It isn't the initial condition.

Looked at another way, suppose you start the system at some arbitrary y1 value at some future time t1 and then ran the system backwards in time until you got to t=0. You would end up at a value other than y0 ... i.e., a value other than your supposed initial condition.

Again, this is why this differential equation as written doesn't make sense to me. At least not if the y0 in the differential equation is supposed to be the initial value of y. If y0 is the initial value of y, then the only solution is going to be the constant solution

y(t) = y0

[David Goodmanson]

Hello NH/JT

The eqn is

y'(t) = f(y,t)

and in this case f is a function of y only, with no explicit dependence on t. That means that if y(t) is a solution, then y(t-t0) is also a solution, and you have translational invariance. So you can put t = 0 basically anywhere you want along the curve as the initial condition (except where y' = 0).

[James Tursa]

Translational invariance isn't the issue. The issue is that the differential equation, as written, has only one solution which is the constant solution y(t) = y0 if we are to assume that the y0 in the derivative function is the initial value of y. I don't care what the value of t is at that initial condition point. That doesn't affect the solution.

I would like to see where this differential equation comes from ...

[David Goodmanson]

Why are you assuming that y0 is the initial value of y? y0 appears to be just a constant in a particular function, just like the .003445. Why couldn't, say, .08 be set as the initial value of y?

[James Tursa]

Because OP wrote:

"... initial value for y=0.072 ..."

Some clarification from OP would help here, as well as some info on where this differential equation comes from.

[David Goodmanson]

I agree about clarification from the OP, who wrote "When the initial value for y=0.072 I get the following graph" but who also later wrote "I decided to give y=0.073 as the initial value" so it did not appear to me that the initial value = y0 was set in stone.

Answer 1

[David Goodmanson]

Hi nadeeshani,

UPDATED ANSWER

This is a speculative answer, but perhaps you did not take the square root when you calculated the constant. At any rate, plugging in sqrt(.003445) = 0.0587 in place of .003445 (NOTE: and using initial value 0.073) gives a plot that is very close to the expected one.