Free Rotation of an Asymmetric Top
Free Rotation of an Asymmetric Top
a=4;b=2;c=1;q1=0.01;q2=3;q3=1;m=1000;T=5;kk=0.0005;
I1=m(b*b+c*c)/12;I2=m(c*c+a*a)/12;I3=m(a*a+b*b)/12;
EE=0.5(I1*q1*q1+I2*q2*q2+I3*q3*q3);MM=Sqrt[I1*I1*q1*q1+I2*I2*q2*q2+I3*I3*q3*q3];
range=kkMax[Sqrt[2EE*I1],Sqrt[2EE*I2],N@Sqrt[2EE*I3]];
sgnQ1=Sign[q1];sgnQ3=Sign[q3];
Ω[t_]:={Ω1[t],Ω2[t],Ω3[t]};AngularMomentum[t_]:=rotationMatrix[t].{I1Ω1[t],I2Ω2[t],I3Ω3[t]};
If2EE*I2<,Block{ch,ss,k2,AA,Δt},ss=q2;k2=;ch=;AA=-NIntegrate,{s,0,ss};Δt=AAch;Period=4chNIntegrate,{s,0,1};Ω1[t_]=sgnQ1JacobiCN[(t-Δt)/ch,k2];Ω2[t_]=JacobiSN[(t-Δt)/ch,k2];Ω3[t_]=sgnQ3JacobiDN[(t-Δt)/ch,k2];,Block{ch,ss,k2,AA,Δt},ss=q2;k2=;ch=;AA=-NIntegrate,{s,0,ss};Δt=AAch;Period=4chNIntegrate,{s,0,1};Ω3[t_]=sgnQ3JacobiCN[(t-Δt)/ch,k2];Ω2[t_]=JacobiSN[(t-Δt)/ch,k2];Ω1[t_]=sgnQ1JacobiDN[(t-Δt)/ch,k2];;
2
MM
I2(I3-I2)
2EEI3-MM*MM
(I2-I1)(2EEI3-MM*MM)
(I3-I2)(MM*MM-2EEI1)
I1I2I3
(I3-I2)(MM*MM-2EEI1)
1
(1-s*s)(1-k2*s*s)
1
(1-s*s)(1-k2*s*s)
2EEI3-MM*MM
I1(I3-I1)
2EEI3-MM*MM
I2(I3-I2)
MM*MM-2EEI1
I3(I3-I1)
I2(I1-I2)
2EEI1-MM*MM
(I2-I3)(2EEI1-MM*MM)
(I1-I2)(MM*MM-2EEI3)
I1I2I3
(I1-I2)(MM*MM-2EEI3)
1
(1-s*s)(1-k2*s*s)
1
(1-s*s)(1-k2*s*s)
2EEI1-MM*MM
I3(I1-I3)
2EEI1-MM*MM
I2(I1-I2)
MM*MM-2EEI3
I1(I1-I3)
rotationMatrix[t_]=
;
a11[t] | a12[t] | a13[t] |
a21[t] | a22[t] | a23[t] |
a31[t] | a32[t] | a33[t] |
solution=NDSolveFlatten@Thread/@ThreadrotationMatrix[t]rotationMatrix[t].
,Flatten@(Thread/@Thread[rotationMatrix[0]IdentityMatrix[3]]),{a11,a12,a13,a21,a22,a23,a31,a32,a33},{t,0,T},MaxStepsInfinity;
∂
t
0 | -Ω3[t] | Ω2[t] |
Ω3[t] | 0 | -Ω1[t] |
-Ω2[t] | Ω1[t] | 0 |
aM=First[{a11,a12,a13,a21,a22,a23,a31,a32,a33}/.solution];
rotationMatrix[t_]:={{aM[[1]][t],aM[[2]][t],aM[[3]][t]},{aM[[4]][t],aM[[5]][t],aM[[6]][t]},{aM[[7]][t],aM[[8]][t],aM[[9]][t]}};
points1=Table[{t,Ω1[t]},{t,0,T,0.1}];points2=Table[{t,Ω2[t]},{t,0,T,0.1}];points3=Table[{t,Ω3[t]},{t,0,T,0.1}];ΩI1=Interpolation[points1];ΩI2=Interpolation[points2];ΩI3=Interpolation[points3];
Brick[a_,b_,c_]:={Green,GraphicsComplex[{{-a,-b,-c},{a,-b,-c},{a,b,-c},{-a,b,-c},{-a,-b,c},{a,-b,c},{a,b,c},{-a,b,c}},Polygon[{{1,4,3,2},{3,4,8,7},{1,5,8,4},{2,3,7,6},{5,6,7,8},{1,2,6,5}}],VertexColors{Blue,Blue,Blue,Blue,Green,Green,Green,Green}]};
Arrow3D[{x0_,y0_,z0_},{x1_,y1_,z1_},{HeadSizeH_,HeadSizeW_,ShaftSize_,color_}]:=Module[{Vertex1,Vertex2,Vertex3,Vertex4,C,D1,D2,D3,D4,k1,k2,ch,p},k1=HeadSizeW;k2=HeadSizeH*Norm[{x1-x0,y1-y0,z1-z0}];p=1;ch=0;If[z1z0,If[y1y0,If[x1x0,p=2,ch=1],ch=2],ch=3];C=k1{x1-x0,y1-y0,z1-z0}+{x0,y0,z0};D1={0,0,0};D2={0,0,0};D3={0,0,0};D4={0,0,0};If[ch3,D1[[1]]=1;D1[[2]]=0;D1[[3]]=((-x0+x1)(1-x0-k1(-x0+x1))+(-y0+y1)(-y0-k1(-y0+y1)))/(z0-z1)+z0+k1(-z0+z1);];If[ch2,D1[[1]]=1;D1[[3]]=0;D1[[2]]=((-x0+x1)(1-x0-k1(-x0+x1))+(-z0+z1)(-z0-k1(-z0+z1)))/(y0-y1)+y0+k1(-y0+y1);];If[ch1,D1[[2]]=1;D1[[3]]=0;D1[[1]]=((-y0+y1)(1-y0-k1(-y0+y1))+(-z0+z1)(-z0-k1(-z0+z1)))/(x0-x1)+x0+k1(-x0+x1);];D1=D1-C;D2=Cross[{x1-x0,y1-y0,z1-z0},D1];D3=k2Normalize[D1];D4=k2Normalize[D2];D1=-D3;D2=-D4;D3=D3+C;D4=D4+C;D1=D1+C;D2=D2+C;Vertex1=D3;Vertex2=D4;Vertex3=D1;Vertex4=D2;{{Thickness[ShaftSize],color,Line[{{x0,y0,z0},{x1,y1,z1}}],GraphicsComplex[{{x1,y1,z1},Vertex1,Vertex2,Vertex3,Vertex4},Polygon[{{1,2,3},{1,3,4},{1,4,5},{1,5,2}}]]},{color,Point[{x0,y0,z0}]}}[[p]]];