Filter löschen
Filter löschen

Please help me solve this differential equation

2 Ansichten (letzte 30 Tage)
Armel Kapso
Armel Kapso am 18 Jul. 2022
Kommentiert: Torsten am 29 Jul. 2022
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 Kommentare
John D'Errico
John D'Errico am 18 Jul. 2022
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 am 19 Jul. 2022
Please take a2=3 a3=2, b2=5 and b3=2

Melden Sie sich an, um zu kommentieren.

Akzeptierte Antwort

Sam Chak
Sam Chak am 19 Jul. 2022
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 Kommentare
Armel Kapso
Armel Kapso am 29 Jul. 2022
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 am 29 Jul. 2022
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))

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (1)

Chunru
Chunru am 19 Jul. 2022
Bearbeitet: Chunru am 19 Jul. 2022
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 Kommentare
Chunru
Chunru am 19 Jul. 2022
You need to save the code in a file, e.g., testdiff.m. Then run the program testdiff.m

Melden Sie sich an, um zu kommentieren.

Community Treasure Hunt

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

Start Hunting!

Translated by