pythonでMatplotlib使おう！

では、今度は
f(x,y)=exp(-X*X-Y*Y)*(2*X)*(4*Y*Y-2)
をプロットしてみましょう。

----python---script---

import math
from pylab import *
def gaussfunc(x,y,a,b):
   return exp(-a*(x)*(x)-b*(y)*(y))

dx, dy = 0.05, 0.05
x = arange(-3.0, 3.0, dx)
y = arange(-3.0, 3.0, dy)
X,Y = meshgrid(x, y)
Z = gaussfunc(X, Y,1.0,0.5)*(2*X)*(4*Y*Y-2)
im = imshow(Z,origin='lower' ,alpha=.9)
colorbar(im)

cset = contour(Z, arange(-2.0,2.0,0.2), origin='lower'  )

clabel(cset,
inline=1,
fmt='%1.1f',
fontsize=10)

hot()
show()

----python script----


ちょっと不思議な絵になりましたね？


一応物理的な解釈をこの絵に持ち込むことも可能です。
量子力学での調和振動子の波動関数の解は基本的に（エルミート多項式）＊ガウス関数とかけます。任意の多項式はエルミート多項式に分解出来る事を使うと、上の濃淡プロットは波動関数をプロットしていて、波動関数はn=3までの状態の重ね合わせで描けているというように。
エルミート多項式に関しては
http://homepage2.nifty.com/eman/quantum/oscillator.html
に詳しく書いてあります。
では、どういう重ね合わせなのかというと、
H0(x)=1
H1(x)=2x
H2(x)=4x^2-2
f0(x)=exp(-x^2)
f1(x)=2xexp(-x^2)
f2(x)=(4x^2-2)exp(-x^2)
とかけるので、
(2x)(4y^2-2)=H1(x)*H2(y)
よって、波動関数f(x,y)=f1(x)*f2(y)となります。
これは何を意味しているかというと、0や1や2は節の数を表しています。
x方向に１個の節で、y方向には２個の節をもった波動関数を描いて
振動している事が分かります。