Quit

eqs={m1 x1''[t]==-T1 Sin[θ1[t]]+T2 Sin[θ2[t]],
m1 y1''[t]==T1 Cos[θ1[t]]-T2 Cos[θ2[t]]-m1 g,
m2 x2''[t]==-T2 Sin[θ2[t]],
m2 y2''[t]==T2 Cos[θ2[t]]-m2 g}

x1[t_]=l1 Sin[θ1[t]]
y1[t_]=-l1 Cos[θ1[t]]
x2[t_]=x1[t]+l2 Sin[θ2[t]]
y2[t_]=y1[t]-l2 Cos[θ2[t]]

eqs

eqs//MatrixForm

Solve[{eqs[[2]],eqs[[4]]},{T1,T2}]//FullSimplify

sistema={eqs[[1]],eqs[[3]]}/.{T1->Sec[θ1[t]] (l2 m2 Cos[θ2[t]] (((θ2)'[t])^2)+(m1+m2) (g+l1 Cos[θ1[t]] (((θ1)'[t])^2)+l1 Sin[θ1[t]] (θ1)''[t])+l2 m2 Sin[θ2[t]] (θ2)''[t]),T2->m2 Sec[θ2[t]] (g+l1 Cos[θ1[t]] (((θ1)'[t])^2)+l2 Cos[θ2[t]] (((θ2)'[t])^2)+l1 Sin[θ1[t]] (θ1)''[t]+l2 Sin[θ2[t]] (θ2)''[t])}//FullSimplify

TraditionalForm[MatrixForm[sistema]]

condin={θ1[0]==θinicial,θ1'[0]==0,θ2[0]==θinicial,θ2'[0]==0}

entrada=Join[sistema,condin]

m1=10;m2=1;l1=1;l2=0.5;g=10;
θinicial=(π/2)

entrada

sol=NDSolveValue[entrada,{θ1,θ2},{t,0,30}]
θ1sol=sol[[1]];
θ2sol=sol[[2]];

Plot[{θ1sol[t],θ2sol[t]},{t,0,30}]

pos1[t_]={x1[t],y1[t]}/.{θ1->θ1sol,θ2->θ2sol};
pos2[t_]={x2[t],y2[t]}/.{θ1->θ1sol,θ2->θ2sol};

ParametricPlot[pos1[t],{t,0,30}]

ParametricPlot[pos2[t],{t,0,30}]

Graphics[{Point[{0,0}],Point[pos1[3]],Point[pos2[3]]}]

bolinhas[t_]:=Graphics[{PointSize[0.02],Red,Point[{0,0}],Point[pos1[t]],Point[pos2[t]]}]

hastes[t_]:=Graphics[Line[{{0,0},pos1[t],pos2[t]}]]

rastro[t_]:=ParametricPlot[pos2[τ],{τ,Max[t-0.5,0],t}]

Manipulate[Show[bolinhas[t],hastes[t],rastro[t],PlotRange->1.2 (l1+l2)],{t,0.000001,30}]

Δt=0.1;     (*  intervalo de tempo entre frames  *)
fps=(1/(Δt));  (*  frames por segundo  *)

tab=Table[Show[bolinhas[t],hastes[t],rastro[t],PlotLabel->"(Subscript[m, 1]) dez vezes maior que (Subscript[m, 2])",PlotRange->1.2 (l1+l2)],{t,0.000001,30,Δt}];

ListAnimate[tab,fps,AnimationRunning->False]

Export["penduloduplo.gif",tab,"DisplayDurations"->(1/fps),AnimationRepetitions->∞,ImageSize->300]