二つの重りがバネに繋がれている場合の連成振動

二つ以上の振動子がお互いに作用を及ぼし合いながら振動する場合、これを連成振動という．下図のような場合の運動方程式をとき、振動運動をアニメーションにしてデモンストレートする．
Mathematicaはドローイング用のプログラムではないので次のインプットの例で明らかなように絵を書く場合にはやや面倒である．しかし後にわかるように、すばらしいアニメーションを作り出すことができる．

注意：うまくアニメーションが作れないとき

パラメターの値が固定されている場合、うまくアニメーションが作れないときがあります．パラメーターをリセットしてもう一度、トライしてください．

ClearAll["Global`*"]

パラメター設定

a1=0;a2=2;Show[Graphics[{Hue[0.8],PointSize[0.09],Point[{a2,0}]}],Graphics[{Hue[1],PointSize[0.09],Point[{a1,0}]}],Graphics[{Hue[0.3],Line[{{a1,0},{((a2-a1)/7)+a1,(1/2)}}],Line[{{a2-((a2-a1)/7),-(1/2)},{a2,0}}],Table[Line[{{(1/7) (n-1) (a2-a1)+a1,Mod[n-1,2]-(1/2)},{a1+(1/7) (a2-a1) n,Mod[n,2]-(1/2)}}],{n,2,6}]}],Graphics[{Hue[0.5],Line[{{-2,0},{((a1+2)/11)-2,1}}],Line[{{a1-((a1+2)/11),-1},{a1,0}}],Table[Line[{{(1/11) (a1+2) (n-1)-2,2 Mod[n-1,2]-1},{(1/11) (a1+2) n-2,2 Mod[n,2]-1}}],{n,2,10}],Line[{{a2,0},{((4-a2)/11)+a2,1}}],Line[{{4-((4-a2)/11),-1},{4,0}}],Table[Line[{{(1/11) (n-1) (4-a2)+a2,2 Mod[n-1,2]-1},{a2+(1/11) (4-a2) n,2 Mod[n,2]-1}}],{n,2,10}]}],Graphics[{Text["k",{-1,1.5}],Text["λ",{1,1.5}],Text["k",{3,1.5}],Text[(Subscript[m, 1]),{0,-1}],Text[(Subscript[m, 2]),{2,-1}]}],PlotRange->{{-2,4},{-1.7,2.2}}]

以下のように下付けを定義しておくとSubscript[x, 1]をｘの関数と扱わないので便利である．

<<Notation`

Symbolize[(NotationTemplateTag[(Subscript[expr_, 1])])]

Symbolize[(NotationTemplateTag[(Subscript[expr_, 2])])]

運動方程式は各々の重りに対してニュートンの第二法則を使って導いてもよいが、ここではLagrangeの方法を使う．

K=(1/2) (Subscript[m, 1]) (((Subscript[x, 1]^′)[t])^2)+(1/2) (Subscript[m, 2]) (((Subscript[x, 2]^′)[t])^2);U=(1/2) (Subscript[k, 1]) (((Subscript[x, 1])[t])^2)+(1/2) (Subscript[k, 2]) (((Subscript[x, 2])[t])^2)+(1/2) λ ((((Subscript[x, 2])[t]-(Subscript[x, 1])[t]))^2);
Lagrange= K-U;

微分方程式はDSolveで解ける

SetAttributes[{(Subscript[m, 1]),(Subscript[m, 2])},Constant]

解くべき微分法方程式は

{Dt[(∂_((Subscript[x, 1]^′)[t]))Lagrange,t]-(∂_((Subscript[x, 1])[t]))Lagrange==0,Dt[(∂_((Subscript[x, 2]^′)[t]))Lagrange,t]-(∂_((Subscript[x, 2])[t]))Lagrange==0}

resultD=DSolve[{(Subscript[k, 1]) (Subscript[x, 1])[t]-λ (-(Subscript[x, 1])[t]+(Subscript[x, 2])[t])+(Subscript[m, 1]) (Subscript[x, 1])''[t]==0,(Subscript[k, 2]) (Subscript[x, 2])[t]+λ (-(Subscript[x, 1])[t]+(Subscript[x, 2])[t])+(Subscript[m, 2]) (Subscript[x, 2])''[t]==0},{(Subscript[x, 1])[t],(Subscript[x, 2])[t]},t]

微分方程式は確かにDSolveでも解けるが、解は複雑すぎるので、ここでは数値的に解き、アニメーションを作ることにする．

アニメーションを作ってみよう simulation
ここでは数値的に解くことにする．

(Subscript[m, 1])=1;(Subscript[m, 2])=1;λ=.05;(Subscript[k, 1])=1;(Subscript[k, 2])=1;
result=NDSolve[{Dt[(∂_((Subscript[x, 1]^′)[t]))Lagrange,t]-(∂_((Subscript[x, 1])[t]))Lagrange==0,Dt[(∂_((Subscript[x, 2]^′)[t]))Lagrange,t]-(∂_((Subscript[x, 2])[t]))Lagrange==0,(Subscript[x, 1]^′)[0]==0, (Subscript[x, 2]^′)[0]==0,(Subscript[x, 1])[0]==-1,(Subscript[x, 2])[0]==0},{(Subscript[x, 1])[t],(Subscript[x, 2])[t]},{t,0,130}]

この計算の結果を横軸を時間に縦軸を振幅として図に描く．

Show[Plot[(Subscript[x, 1])[t]/.result, {t, 0, 130}, PlotStyle->Hue[1]],Plot[(Subscript[x, 2])[t]/.result, {t, 0, 130},PlotStyle->Hue[.8]], AxesLabel->{"time", "amplitude"}]

上図より運動のエネルギーが一方の質量から他方に漸次移り、また戻ってくるという様子が窺える．実際、以下をを入力すれば沢山の図がでてきてアニメーションが楽しめる．

a1=(Subscript[x, 1])[t]/.result[[1]];a2=2+(Subscript[x, 2])[t]/.result[[1]];Manipulate[Evaluate[Show[Graphics[{Hue[0.8],PointSize[0.09],Point[{a2,0}]}],Graphics[{Hue[1],PointSize[0.09],Point[{a1,0}]}],Graphics[{Hue[0.3],Line[{{a1,0},{((a2-a1)/7)+a1,(1/2)}}],Line[{{a2-((a2-a1)/7),-(1/2)},{a2,0}}],Table[Line[{{(1/7) (n-1) (a2-a1)+a1,Mod[n-1,2]-(1/2)},{a1+(1/7) (a2-a1) n,Mod[n,2]-(1/2)}}],{n,2,6}]}],Graphics[{Hue[0.5],Line[{{-2,0},{((a1+2)/11)-2,1}}],Line[{{a1-((a1+2)/11),-1},{a1,0}}],Table[Line[{{(1/11) (a1+2) (n-1)-2,2 Mod[n-1,2]-1},{(1/11) (a1+2) n-2,2 Mod[n,2]-1}}],{n,2,10}],Line[{{a2,0},{((4-a2)/11)+a2,1}}],Line[{{4-((4-a2)/11),-1},{4,0}}],Table[Line[{{(1/11) (n-1) (4-a2)+a2,2 Mod[n-1,2]-1},{a2+(1/11) (4-a2) n,2 Mod[n,2]-1}}],{n,2,10}]}],PlotRange->{{-2,4},{-1.7,2.2}}]], {t, 0, 130, .5}]

list=Table[Evaluate[Show[Graphics[{Hue[0.8],PointSize[0.09],Point[{a2,0}]}],Graphics[{Hue[1],PointSize[0.09],Point[{a1,0}]}],Graphics[{Hue[0.3],Line[{{a1,0},{((a2-a1)/7)+a1,(1/2)}}],Line[{{a2-((a2-a1)/7),-(1/2)},{a2,0}}],Table[Line[{{(1/7) (n-1) (a2-a1)+a1,Mod[n-1,2]-(1/2)},{a1+(1/7) (a2-a1) n,Mod[n,2]-(1/2)}}],{n,2,6}]}],Graphics[{Hue[0.5],Line[{{-2,0},{((a1+2)/11)-2,1}}],Line[{{a1-((a1+2)/11),-1},{a1,0}}],Table[Line[{{(1/11) (a1+2) (n-1)-2,2 Mod[n-1,2]-1},{(1/11) (a1+2) n-2,2 Mod[n,2]-1}}],{n,2,10}],Line[{{a2,0},{((4-a2)/11)+a2,1}}],Line[{{4-((4-a2)/11),-1},{4,0}}],Table[Line[{{(1/11) (n-1) (4-a2)+a2,2 Mod[n-1,2]-1},{a2+(1/11) (4-a2) n,2 Mod[n,2]-1}}],{n,2,10}]}],PlotRange->{{-2,4},{-1.7,2.2}}]], {t, 0, 130, .5}];

Export["TwoCoupledHarmomnicOscillators.gif",list]