Please help me with the code

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 个评论

here x=1/(8y^(k-1)-3), Now put the K values 0.1000 0.1100 0.1200 0.1300 0.1400 0.1500 0.1600 0.1700 0.1800 0.1900 0.2000... After that question is not clear. y0=0.6??
Why can't you use arrayfun or find()? Usually these kinds of weird restrictions are only for homework problems. But you didn't tag it as homework. Is it? http://www.mathworks.com/matlabcentral/answers/8626-how-do-i-get-help-on-homework-questions-on-matlab-answers
Jan
Jan 2017-8-27
编辑:Jan 2017-8-27
@Nicia: "x1 is at y1" is not clear to me. Actually all you need is a loop over the elements of k and a call of ODE45. But the "y1=2" is confusing. Is this an initial or boundary value problem?
If this is a homework, I will not post the complete solution to give you a chance to submit your solution.
@KALAN, @Jan Simon: Thank you for replying my post. This question asks to find x1, when y is equal 2, it's called y1. At that time, x is equal x1. I understand how to use ODE45 to find the solution of the equation. This is what I do. I create an m.file
[ dydx ] = func( x,y)
k=0.1
dydx=(8.*(y.^k))-(3.*y);
[x,y]=ode45(@func,[0:0.5:10],0.6)
The output will be y over x. However, I don't know how to find x1 without using find or arrayfun function. Do you know any functions that I can use in this situation? And since k is a range that means I have to create many m.file for each k value. Is there any better ways for me, instead? I just want to be given a hint not a complete solution. Thank you in advance
@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.
@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 个)

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 个评论

Thank you for helping me. The question asks to find x1( x is equal x1 when y reaches 2 called y1). I have found the solution of the equation. The output is y over x. But I'm having difficult time to find x1 without using find or arrayfunc function. Please guide me what I should do next. Thank you in advance.
Does this mean that you want to stop the integration when y=2 is reached? Then an event function will help: See doc odeset. fzero might be useful also. I have no idea, why arrayfun could be used here.
Do not create different M-files for the elements of k, but a loop:
kList = 0.1:0.01:0.2;
for ik = 1:length(kList)
k = kList(ik);
dy = @(t, y) 8*y^k - 3*y;
...
end
Is this a homework for numerics and you have heard something about boundary value problems with free end time?
I haven't heard anything about boundary value problems with free and time. What I intend to do is like drawing out x1 from the solution got by using ode45 to solve the equation. And x1 satisfies a condition that y=2. I don't if any function can help me with that. I have been searching around but I couldn't find it. That's why I come here and ask for help. However, there maybe another way to solve this problem better than mine.
In your code: dy=@(t,y).. what is 't', does you mean 'x'?
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".
Thank you so much for your help. I greatly appreciate it.

请先登录,再进行评论。

类别

帮助中心File Exchange 中查找有关 Programming 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by