Asked by Takey Asaad
on 14 Sep 2018 at 23:29

my function is

dy/dt=k*y*exp(450/y)

k is constant and y(0)=40 and y(15)=95 solve this equation by using ode45 can someone pleaseeeeeeeeeee check the code and make it work .

tspan = [0 300]; y0 = 40;y15=95 [t,y] = ode45(@(t,y) 'k'*y*exp(450/y), tspan, y0,y15); plot(t,y,'-o')

Answer by Stephan Jung
on 16 Sep 2018 at 20:39

Edited by Stephan Jung
on 17 Sep 2018 at 10:05

Accepted Answer

Hi,

if you search a value for k that complies with the boundary conditions, you could use ** fzero** to solve the problem numerically:

k = fzero(@calculate_k, 0.001); disp(['k = ', num2str(k)]) [t,y] = ode45(@(t,y) k*y*exp(450/y), [0 300], 40); plot(t,y,'r') hold on scatter(15,95,'ob') text(23,96,'y(t=15) = 95','Color','b') hold off

function k_value = calculate_k(x) tspan = [0 15]; y0 = 40; [~,y] = ode45(@(t,y) x*y*exp(450/y), tspan, y0); k_value = y(end) - 95; end

This will give you as result:

k = 0.00010435

or if you need it more accurate:

k =

1.043495007807761e-04

and the plot which belongs to this result looks like the boundary conditions are met (To be precise, it keeps the boundary conditions. This is because ** fzero** is a powerful tool and also because you can think of a very good initial value for x0 with just a few estimates...):

.

Best regards

Stephan

Takey Asaad
on 17 Sep 2018 at 13:03

Thanks so much

Stephan Jung
on 17 Sep 2018 at 13:16

No problem - please be patient with your questions. Dont ask the same question a second time, if it needs more time than you expected to get an answer.

Better try to give as much as and as good as possible informations that are needed to let the contributers give you a useful answer.

Sign in to comment.

Answer by James Tursa
on 15 Sep 2018 at 0:29

Remove the quotes from 'k', and be sure to define k before you call ode45. Also, ode45 is an initial value problem solver, so the y15 variable is not applicable (remove it from the call).

Takey Asaad
on 15 Sep 2018 at 1:53

Takey Asaad
on 15 Sep 2018 at 1:56

Takey Asaad
on 15 Sep 2018 at 2:02

Sign in to comment.

Answer by Takey Asaad
on 15 Sep 2018 at 2:04

if true function bvp4 xlow=0;xhigh=300; solinit=bvpinit(linspace(xlow,xhigh,300),[40 95]); sol=bvp4c(@bvp4ode,@bvp4bc,solinit); xinit=linspace(xlow,xhigh,20); sxinit=deval(sol,xint); plot(xinit,sxinit(1,:); function dydx=bvp4ode(x,y) dydx=[y(2) ,a*y*exp(450/y)]; function res=bvp4bc(ya,yb) res=[ya(1)-1,yb(1)] end

Sign in to comment.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.