連成振動

n個の原子からなる一次元chainを考えよう．i番目の原子の座標をaiとするとニュートンの第二法則より

m (((∂^2) (Subscript[a, i]))/(∂(t^2)))==c ((Subscript[a, i - 1]) - (Subscript[a, i])) - c ((Subscript[a, i]) - (Subscript[a, i + 1]))

という運動方程式を解けばよい．ここで、

(((∂^2) (Subscript[a, i]))/(∂(t^2))) 　==-(ω^2) (Subscript[a, i])

であることに注意．従って、k == c/m として、次の行列の固有値を求めれば角振動数ωが求まる．ここでは境界では自由端と考える（nが大きくなると自由端であろうが固定端であろうが結果にあまり影響がないが）．

a[n_] := 
      Table[Which[i == j, 2 k, i == j + 1, -k, i == j - 1, -k, True, 
          0], {i, n}, {j, n}]

たとえばn=10とすると、

MatrixForm[a[10]]

という行列になる．さて、この行列の固有値を次のように定義する．（便宜上k =1としておく．）

b[n_]  :=  N[Eigenvalues[a[n] /. k -> 1]]

求めた固有値を小さい順から並べてプロットする．

g[n_] := 
      ListPlot[Sort[Sqrt[b[n]]],  
        PlotRange -> {{0,  n + 1},  {0,  2.5}},  
        AxesLabel -> {"P",  "\:30a8\:30cd\:30eb\:30ae\:30fc"},  
        PlotStyle -> {Hue[1/12 n],  PointSize[ .02]}]

例えばn=48とすると、

g[48]

また、n= 1から12のでのプロットを一緒に描くと、

Show[Table[g[n],  {n, 12}],  
      PlotRange -> {{0,  12},  {0,  2}}]

となる．

ここで横軸は波数(== 2 π/λ)に対応する．

課題（optional)
 2種類の原子がn 組、交互に１次元的に chainを形成しているとする．原子間距離をa 
 として、上と同じように、波数とエネルギーの関係（分散）を計算してみよ．