problem in creating function file using ode45

2023 3 10
2 个回答
回答已采纳
更新时间：2023 3 12
3 次查看（30 天）

0 个投票

I am trying to solve a second order coupled ode using ode45, however I am facing issues in making it function file and providing variables as input. How to convert symbolic input to numerical in main file using ode45.

请先登录，再回答此问题。

 采纳的回答

Torsten
Torsten 2023-3-11
编辑：Torsten 2023-3-11
Ran in:

0 个投票

Ran in:
Be careful with the order of the solution vector Y. It is (y, Dy, x, Dx). You have to supply the vector of initial values "ic" in the same order !
rng("default")
clc
clear all
close all
O=rand;
a=rand;
g=9.81;
L=rand;
[y]=coupled_ode(O, a, g, L);
VF = 
Subs = 
ftotal = function_handle with value:
@(t,Y)[Y(2);Y(1).*(-7.725211393067794e+1)-sin(9.057919370756192e-1).*Y(4).*1.629447372786358;Y(4);Y(3).*(-7.725211393067794e+1)+sin(9.057919370756192e-1).*Y(2).*1.629447372786358]
Disp=y;
function [y] = coupled_ode(Onum, anum, gnum, Lnum)
syms O a g L x(t) y(t) t Y
dx = diff(x);
d2x = diff(x,2);
dy = diff(y);
d2y = diff(y,2);
Eq1 = d2x == 2*O*sin(a)*dy - (g/L)*x(t);
Eq2 = d2y == -2*O*sin(a)*dx - (g/L)*y(t);
[VF,Subs] = odeToVectorField(Eq1, Eq2)
VF = subs(VF,[O a g L],[Onum anum gnum Lnum]);
ftotal = matlabFunction(VF,'Vars',{t,Y})
tspan = [0 25]; % Choose Appropriate Simulation Time
ic = [0 1 0 1]; % Choose Appropriate Initial Conditions
[t,y] = ode45(ftotal, tspan, ic);
figure
plot(t, y)
grid
legend(string(Subs))
end

2 个评论
显示 无 隐藏 无

Susmita Panda
Susmita Panda 2023-3-12
Thank you @ Torsten sir. rng(default) helped me.This is exactly what I was searching.I didnt have any idea about this command.
Torsten
Torsten 2023-3-12
编辑：Torsten 2023-3-12
Ran in:
Ran in:
It's always a good idea if programming parts that serve different tasks have their own subfunctions. But in your case, I don't understand why and how you divided your code in the two different .m files.
E.g. to have control over what is done in the function coupled_ode, you will have to pass tspan, the initial conditions for the variables and the parameter values for O, a, g and L. The values returned should be t and y. The plotting part also should be done in the main program because it makes no sense to do the plotting in coupled_ode while passing the results as output parameters back to the calling program.
I'd arrange the code as
rng("default")
clc
clear all
close all
O = rand;
a = rand;
g = 9.81;
L = rand;
f = function_to_integrate(O, a, g, L)
f = function_handle with value:
@(t,Y)[Y(2);Y(1).*(-7.725211393067794e+1)-sin(9.057919370756192e-1).*Y(4).*1.629447372786358;Y(4);Y(3).*(-7.725211393067794e+1)+sin(9.057919370756192e-1).*Y(2).*1.629447372786358]
tstart = 0.0;
tend = 25.0;
tspan = linspace(tstart,tend,100); % Choose Appropriate Simulation Time
ic = [0 1 0 1]; % Choose Appropriate Initial Conditions
[t,y] = ode45(f, tspan, ic);
figure
plot(t, y)
grid
legend(string(Subs))
Unrecognized function or variable 'Subs'.
function f = function_to_integrate(Onum, anum, gnum, Lnum)
syms O a g L x(t) y(t) t Y
dx = diff(x);
d2x = diff(x,2);
dy = diff(y);
d2y = diff(y,2);
Eq1 = d2x == 2*O*sin(a)*dy - (g/L)*x(t);
Eq2 = d2y == -2*O*sin(a)*dx - (g/L)*y(t);
[VF,Subs] = odeToVectorField(Eq1, Eq2);
VF = subs(VF,[O a g L],[Onum anum gnum Lnum]);
f = matlabFunction(VF,'Vars',{t,Y});
end

请先登录，再进行评论。

更多回答（1 个）

Alan Stevens
Alan Stevens 2023-3-10
Ran in:

0 个投票

Ran in:
Are you looking for something like this? (I might have misinterpreted your equations!)
O=1;
a=pi/2;
g=9.81;
L=100;
tspan = [0 15]; % Choose Appropriate Simulation Time
ic = [0 1 0 1]; % Choose Appropriate Initial Conditions
[t,y] = ode45(@(t,y) ftotal(t,y,O,a,g,L), tspan, ic);
figure
plot(t, y)
grid
legend('y1','dy1/dt','y2','dy2/dt')
function dydt=ftotal(~,y,O, a, g, L)
y1 = y(1); v1 = y(2); y2 = y(3); v2 = y(4);
dydt = [ v1;
2*O*sin(a)*y2 - g/L;
v2;
-2*O*sin(a) - g/L*y1];
end

6 个评论
显示 4更早的评论 隐藏 4更早的评论

Susmita Panda
Susmita Panda 2023-3-11
编辑：Susmita Panda 2023-3-11
@Alan sir, this is a demo question, I have two ode equations coupled in such a such that the acceleration of one variable x also contains y also, e.g xddot=yddot+A+B ydot= xdot+M+N
I have solved successfully using the procedure however i am facing issue in providing variables like in this question.
Is there any way we can defined acceleration of first variable with the second variable using your procedure?
Alan Stevens
Alan Stevens 2023-3-11
Hmm! You need an equation for yddot. What are your actual equations? Try uploading them exactly in mathematical notation (using an image if necessary).
Susmita Panda
Susmita Panda 2023-3-11
Sir please find the attached file containing second order coupled odes. Is there any simple way to do using ode45 because the method I tried having problem of defining variables in another function file.
Torsten
Torsten 2023-3-11
编辑：Torsten 2023-3-11
Ran in:
Ran in:
Call ode45 with the function
fun = @(t,y) [y(2); xddot; y(4); yddot]
where xddot and yddot are computed as
syms a b c d e f g h i t x y xdot ydot xddot yddot
eqn1 = (1-e)*xddot + e*yddot == -2*a*xdot - x - sin(b*t) + c*(x-y) + d*(xdot-ydot)
eqn1 = 
eqn2 = (g-h)*xddot - (i+g)*yddot == -f*(x-y)
eqn2 = 
solve([eqn1 eqn2],[xddot yddot])
ans = struct with fields:
xddot: -(g*x + i*x + g*sin(b*t) + i*sin(b*t) + 2*a*g*xdot - c*g*x + 2*a*i*xdot + c*g*y + e*f*x - d*g*xdot - c*i*x - e*f*y + d*g*ydot + c*i*y - d*i*xdot + d*i*ydot)/(g + i - e*h - e*i) yddot: (f*x - f*y - g*x + h*x - g*sin(b*t) + h*sin(b*t) - 2*a*g*xdot + 2*a*h*xdot + c*g*x - c*g*y - c*h*x - e*f*x + d*g*xdot + c*h*y + e*f*y - d*g*ydot - d*h*xdot + d*h*ydot)/(g + i - e*h - e*i)
Susmita Panda
Susmita Panda 2023-3-11
Thank you @Torsten sir for this method. I will try this.
but is there any way we can call function "coupled_ode" in "Main_file" which I posted without error. This method is convenient to me and produces good result when they variables are called in same file.I am getting error of this kind when I am making it function file and providing variables in "Main_file":
Error using superiorfloat
Inputs must be floats, namely single or double.
Error in odearguments (line 114)
dataType = superiorfloat(t0,y0,f0);
Error in ode45 (line 107)
odearguments(odeIsFuncHandle,odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);
Error in coupled_ode (line 18)
[t,y] = ode45(@(t,y) ftotal(t,y,O,a,g,L), tspan, ic);
Error in Main_file (line 9)
[y]=coupled_ode(O, a, g, L);
Walter Roberson
Walter Roberson 2023-3-12
After you get the xddot and yddot using solve(), call odeFunction to create something to pass to ode45() .
I recommend following the workflow shown in the first example to odeFunction

请先登录，再进行评论。

类别

在 帮助中心 和 File Exchange 中查找有关 Numeric Solvers 的更多信息

产品

版本

R2022b

标签

Community Treasure Hunt

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

Start Hunting!

Translated by