Effacer les filtres
Effacer les filtres

Please help me solve this differential equation

2 vues (au cours des 30 derniers jours)
Armel Kapso
Armel Kapso le 18 Juil 2022
Commenté : Torsten le 29 Juil 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 commentaires
John D'Errico
John D'Errico le 18 Juil 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 le 19 Juil 2022
Please take a2=3 a3=2, b2=5 and b3=2

Connectez-vous pour commenter.

Réponse acceptée

Sam Chak
Sam Chak le 19 Juil 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 commentaires
Armel Kapso
Armel Kapso le 29 Juil 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 le 29 Juil 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))

Connectez-vous pour commenter.

Plus de réponses (1)

Chunru
Chunru le 19 Juil 2022
Modifié(e) : Chunru le 19 Juil 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 commentaires
Armel Kapso
Armel Kapso le 19 Juil 2022
Chunru
Chunru le 19 Juil 2022
You need to save the code in a file, e.g., testdiff.m. Then run the program testdiff.m

Connectez-vous pour commenter.

Community Treasure Hunt

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

Start Hunting!

Translated by