Question

怎么用matlab解这个公式

A MATLAB script that solves the differential equation for a linearised and the more accurate non-linear pendulum consisting of a mass on the end of a very light string of length l=24.8cm ( g is the acceleration due to gravity 9.8ms -2.)
A diagram and explanation of the symbols.
A plot of the motion of the pendulum showing linear and non-linear pendulums for a starting angle (in radians) of q = 0.5
The error made by a clock using such a pendulum in a day.
谢谢大家!
spritecoca 我发现我用你写的东西算出来不一样的
C2*exp((t*(-L*g)^(1/2))/L) + C3/exp((t*(-L*g)^(1/2))/L)

solve(int((2^(1/2)*(-1/(C20 - (g*cos(y))/L))^(1/2))/2, y) = C22 + t, y)
solve(int(-(2^(1/2)*(-1/(C20 - (g*cos(y))/L))^(1/2))/2, y) = C22 + t, y)
wacs5 能不能详述一下如何找x(1)和x(2)呀

Answer 1

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function dydx=zhidao_septemaugust_lin(t,x)
%%θ=x(1), dθ/dt=x(2)
theta=x(1);
Dtheta=x(2);
dydx=[x(2);-10/0.248*theta];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function dydx=zhidao_septemaugust_nolin(t,x)
%%θ=x(1), dθ/dt=x(2)
theta=x(1);
Dtheta=x(2);
dydx=[x(2);-10/0.248*sin(theta)];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

t_start=0;
t_end=5.1;
theta0=[0.5;0]; %初值
[t_lin,theta_lin]=ode45('zhidao_septemaugust_lin',[0,t_end],theta0);
[t_nolin,theta_nolin]=ode45('zhidao_septemaugust_nolin',[0,t_end],theta0);

plot(t_lin,theta_lin(:,1),t_nolin,theta_nolin(:,1));
legend('linear','non-linear',3);
xlabel('t');

%%%%%%%%%%%%%%%%%
%%zoom the figure, linearised used 4.94s and non-linear system used 5.03s in the five cycle.
%%In other word,5.03s have error 0.09s as a day have error 1500s, because the initial angle is too large.

Answer 2

为了方便书写，暂且把thata写作y,MATLAB程序写作为：
dsolve('D2y+g/L*y=0','t')
dsolve('D2y+g/L*sin(y)=0','t')
结果分别为：
1. C1*sin(1/L^(1/2)*g^(1/2)*t)+C2*cos(1/L^(1/2)*g^(1/2)*t)
2. Int(L/(L*(2*g*cos(_a)+C1*L))^(1/2),_a = .. y)-t-C2 = 0
Int(-L/(L*(2*g*cos(_a)+C1*L))^(1/2),_a = .. y)-t-C2 = 0