Please help me with the code

1 次查看(过去 30 天)
Nicia Nanami
Nicia Nanami 2017-8-27
编辑: Stephen23 2017-8-28
I have an equation:
dy/dx = 8*y^k - 3*y
The given numbers:
k=[0.1:0.01:0.2]
y0=0.6
x=[0 10]
y1=2
How can I use ODE to solve the equation and use that solution to find x1 (x1 is at y1) for every value of k. CANNOT use arrayfun or find.
  6 个评论
Nicia Nanami
Nicia Nanami 2017-8-27
@Image Analyst: yes this is my homework. I'm so sorry for not tagging it as homework since this is my first time using this website. Thank you for reminding me.
Jan
Jan 2017-8-27
编辑:Stephen23 2017-8-28
@Nicia: This looks confused:
[ dydx ] = func( x,y)
k=0.1
dydx = (8.*(y.^k))-(3.*y);
[x,y]=ode45(@func,[0:0.5:10],0.6)
Better:
k = 0.1
dy = @(x,y) 8 * (y.^k) - 3 * y;
[x,y] = ode45(@dy, 0:0.5:10, 0.6)
Remark: Unnecessary square brackets should be avoided (see Why not use square brackets). I prefer the use of optional round parentheses, if they improve the clarity of the code. A multiplication has a higher precedence than the subtraction, so here I would omit the parentheses, but e.g. -2^k could be confused with (-2)^k, so I write -(2^k). It's a question of taste.

请先登录,再进行评论。

回答(1 个)

Jan
Jan 2017-8-27
If this is an initial value problem, all you need is a loop over the elements of k:
kList = 0.1:0.01:0.2;
Now pick a specific k and:
k = kList(1);
dy = @(t, y) 8*y^k - 3*y;
[x, y] = ode45(dy, [0, 10], 0.6);
Is the last element of y the searched value? Then all you need to do is to embed this in a loop and collect the results.
  5 个评论
Jan
Jan 2017-8-27
编辑:Jan 2017-8-27
Yes, the "t" could be a "x" also. But as long as both do not appear in the calculations, this does not matter. < must be the 2nd input, that's all.
As said, you can get the x value, at which y=2 by using an event function for the integrator. This is defined by odeset. This function can stop the integration when y=2 is reached and then you get the corresponding x value as output also:
opt = odeset('Events', @myEventsFcn);
[x, y] = ode45(dy, [0, 10], 0.6, opt);
function [value, isterminal, direction] = myEventsFcn(x, y)
...
Now check for y==2 and set the isterminal flag to stop the integration.
See: doc ode45 and search for "ODE Event Location".
Nicia Nanami
Nicia Nanami 2017-8-28
Thank you so much for your help. I greatly appreciate it.

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by