필터 지우기
필터 지우기

Please help me solve this differential equation

조회 수: 2 (최근 30일)
Armel Kapso
Armel Kapso 2022년 7월 18일
댓글: Torsten 2022년 7월 29일
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
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일
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
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일
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일
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
Chunru
Chunru 2022년 7월 19일
You need to save the code in a file, e.g., testdiff.m. Then run the program testdiff.m

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Ordinary Differential Equations에 대해 자세히 알아보기

태그

Community Treasure Hunt

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

Start Hunting!

Translated by