Clear Filters
Clear Filters

Please help me solve this differential equation

2 views (last 30 days)
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 Comments
John D'Errico
John D'Errico on 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.

Sign in to comment.

Accepted Answer

Sam Chak
Sam Chak on 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 Comments
Armel Kapso
Armel Kapso on 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 on 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))

Sign in to comment.

More Answers (1)

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

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!