Please help me solve this differential equation

2022 7 18
2 个回答
回答已采纳
更新时间：2022 7 29
12 次查看（30 天）

1 个投票

please i will like to solve this eqution for w=188 rad/s
Please take a2=3 a3=2, b2=5 and b3=2

3 个评论
显示 1更早的评论 隐藏 1更早的评论

Torsten
Torsten 2022-7-18
编辑：Torsten 2022-7-18
Solve for x_double_dot, then substitute x = y1, x_dot = y2 and write your equation as a system of two first-order differential equations:
y1' = y2
y2' = x_double_dot = ...
Finally, use ODE45 to solve.
John D'Errico
John D'Errico 2022-7-18
You tell us what the value of w is, but not all of the other important parameters. They are just as important. Thus what are a2, a3, b2 b3? If you don't have them, then no numerical solution will be possible. My expectation is no analytical solution will exist in any case, so a numerical solution is your only option.
Armel Kapso
Armel Kapso 2022-7-19
Please take a2=3 a3=2, b2=5 and b3=2

请先登录，再进行评论。

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

 采纳的回答

Sam Chak
Sam Chak 2022-7-19
Ran in:

1 个投票

Ran in:
Hi @Armel Kapso
Some basic code from the ode45 documentation. Probably looks like this, try manupulating the parameters yourself.
tspan = [0 1]; % simulation time interval
initv = [1 0]; % initial values
[t, x] = ode45(@odefcn, tspan, initv);
plot(t, x(:, 1), 'linewidth', 1.5)
grid on
xlabel({'$t$'}, 'Interpreter', 'latex')
ylabel({'$x(t)$'}, 'Interpreter', 'latex')
title('System Response')
function dxdt = odefcn(t, x)
dxdt = zeros(2, 1);
a2 = 3;
a3 = 2;
b2 = 5;
b3 = 2;
omega = 188;
num = 2*a3*omega^3*sin(x(1) - omega*t)*x(2) + 2*a3*omega*sin(x(1) - omega*t)*x(2)^3 - b3*sin(x(1)) - b3*sin(omega*t);
den = 2*a2 - 2*a3*omega*(2*x(2) + omega)*cos(x(1) - omega*t);
dxdt(1) = x(2);
dxdt(2) = num/den;
end

9 个评论
显示 7更早的评论 隐藏 7更早的评论

Armel Kapso
Armel Kapso 2022-7-26
please this code is not going when i try to run it with my MATLAB R2018b
please help me
Sam Chak
Sam Chak 2022-7-26
@Armel Kapso, I don't know what happened because I don't have R2018b. Moreover, the ode45() was introduced before R2006a.
I indistinctly guess that you copied the code and directly pasted on the Workspace, and pressed [ENTER].
Armel Kapso
Armel Kapso 2022-7-28
@Sam Chak No i didn't it's the fucntion you defined that i copied but i did not see that it was long i mean 'num' and so it was crossing the space allocated for the new function i wanted to add
that's why it was not working
Thanks a Lot Sir @Sam Chak
Sam Chak
Sam Chak 2022-7-28
I'm glad to hear that you have sorted it out, @Armel Kapso.
Armel Kapso
Armel Kapso 2022-7-28
please @Sam Chak help plot the phase space i'm having difficulty in doing it
Sam Chak
Sam Chak 2022-7-28
Ran in:
Ran in:
That's the way to plot the phase.
tspan = [0 1];
initv = [1 0];
[t, x] = ode45(@odefcn, tspan, initv);
plot(x(:, 1), x(:, 2), 'linewidth', 1.5)
Armel Kapso
Armel Kapso 2022-7-29
Hi @Sam Chak please the system i'm studying is actually a chotic system. my challenge now is to plot the Lyapunov exponent. i went through some of the contribution you have shared in that same domain and i have tried to emplement this one but it doesn't work. please help me again.
a2 = 3;
a3 = 2;
b2 = 5;
b3 = 2;
omega = 210;
f = @(t, x) [(2*a3*omega^3*sin(x(1) - omega*t)*x(2) +...
2*a3*omega*sin(x(1) - omega*t)*x(2)^3 +...
- b3*sin(x(1)) - b3*sin(omega*t)/(2*a2 - 2*a3*omega*(2*x(2) + omega)*cos(x(1) - omega*t))];
tspan = linspace(0, 88, 8801);
init = [1 1 1];
[t, x] = ode45(f, tspan, init);
plot3(x(:,1), x(:,2), x(:,3)), view(45, 30)
Torsten
Torsten 2022-7-29
Ran in:
Ran in:
Corrected.
a2 = 3;
a3 = 2;
b2 = 5;
b3 = 2;
omega = 210;
f = @(t,x) [x(2);(2*a3*omega^3*sin(x(1)-omega*t)*x(2)+...
2*a3*omega*sin(x(1)-omega*t)*x(2)^3-...
b3*sin(x(1))-b3*sin(omega*t))/(2*a2-2*a3*omega*(2*x(2)+omega)*cos(x(1)-omega*t))];
tspan = linspace(0, 88, 8801);
init = [1 1];
options = odeset('RelTol',1e-8,'AbsTol',1e-8);
[t, x] = ode15s(f, tspan, init, options);
plot(x(:,1), x(:,2))

请先登录，再进行评论。

更多回答（1 个）

Chunru
Chunru 2022-7-19
编辑：Chunru 2022-7-19
Ran in:

1 个投票

Ran in:
opts = odeset('RelTol', 1e-2, 'AbsTol', 1e-4);
[t, x] = ode45(@diffeqn, [0 1], [0.1; .1], opts);
whos
Name Size Bytes Class Attributes cmdout 1x33 66 char opts 1x1 3712 struct t 4369x1 34952 double x 4369x2 69904 double
plot(t, x(:,2))
function dxdt = diffeqn(t, x)
a2=3; a3=2; b2=5; b3=2; w=188;
% Check this out
dxdt(1, 1) = x(2);
dxdt(2, 1) = ( (2*a3*w^3*sin(x(2)-w*t))*x(1) + ...
(2*a3*w*sin(x(2)-w*t))*x(1).^3 + ...
(-b2*sin(w*t) - b3*sin(x(2))) ) ./ ...
( 2*a2 -2*a3*w*(2*x(1)+w)*cos(x(2)-w*t) );
end

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

Armel Kapso
Armel Kapso 2022-7-19
Yes actually it's a chaotic system. a small change in the parameter produces significant effects
i intend to plot the motion curve of the system and also its space phase trajectories for different values of the parameters
my main problem is that i don't have strong knowledge in matlab
Armel Kapso
Armel Kapso 2022-7-19
i've tested this code on my MATLA R2018b is not going
Chunru
Chunru 2022-7-19
Can you show what is the problem in 2018b, such as error message?
Chunru
Chunru 2022-7-19
I did a correction in the code above.
Chunru
Chunru 2022-7-19
You need to save the code in a file, e.g., testdiff.m. Then run the program testdiff.m

请先登录，再进行评论。

类别

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

产品

标签

Community Treasure Hunt

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

Start Hunting!

Translated by