gnuplot 绘制球谐函数图

因为对形变感兴趣,最近在看 P. Ring & P. Schuck 的书。在形变很小((eta << 1))时,椭球体的表面可以用低阶球谐函数来近似描述。所以,不妨用gnuplot画两个图出来看看,比较直观,也许以后展示的时候也用得着。

1. 球谐函数

使用 Griffiths 的量子力学教材上的相位约定,球谐函数定义为

[Y^m_l ( heta, phi) = epsilon sqrt{ frac{ (2l-1)(l-|m|)! }{ 4 pi ( l+ |m|)! } } e^{imphi} P^{|m|}_l(cos heta), ]

(m geq 0) 时,(epsilon = (-1)^m), (mleq 0) 时有(epsilon=1)。Griffiths的书上写的是(P^m_l(cos heta)),应该是笔误,因为连带勒让德的(m)必须是非负整数,所以我更正为(P^{|m|}_l(cos heta))。Peter Ring 他们的书应该用的相同的相位约定,因为他们提到(Y^*_{lambda mu} = (-1)^mu Y_{lambda -mu}),与这个定义一致。

2. gnuplot 绘制三维参数图

我先尝试了python的绘图,从网上抄了一段代码试了试,发现还没有我熟悉的gnuplot好使,所以决定用gnuplot。

2.1 python 绘制三维图

下面是那段python代码

from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt, numpy as np
fig = plt.figure()
ax = fig.add_subplot(111, projection = '3d' )
u = np.linspace( 0, 2*np.pi, 1000)
v = np.linspace( 0, 2*np.pi, 1000)
x = np.outer( np.cos(u), np.sin(v) )
y = np.outer( np.sin(u), np.sin(v) )
z = np.outer( np.ones(np.size(u)), np.cos(v))
ax.plot_surface( x, y, z, cmap = 'rainbow' )

得到图片:
image
我有以下抱怨:
①三维图不能旋转视角,这点真远不如 gnuplot
②图片小,不够清楚
③代码细节太多,需要手动地定义所有离散的 x, y, z 坐标

2.2 gnuplot 绘制三维图

下面是我写的 gnuplot 脚本,用来绘制同样的三维球面

set param       # 参数图模式
# set pm3d depthorder   # 前面的色块挡住后面的,形成立体感
set urange [0:pi]       # 设置 u, v 的取值范围
set vrange [0:2*pi]
set iso 30              # 值越大取点越多
set view equal xyz      # 显示时 x, y, z 三个坐标标度一致
unset tics              # 不显示坐标轴
unset border            # 不显示边框
splot sin(u)*cos(v), sin(u)*sin(v), cos(u) 

画出来的球体:
image
如果想画彩图,打上灯光,弄得美美的,可以如下:

set param       # 参数图模式
set pm3d depthorder     # 前面的色块挡住后面的,形成立体感
set pm3d lighting primary 0.5 specular 0.6
set palette rgbformulae 7, 5, 15
set urange [0:pi]       # 设置 u, v 的取值范围
set vrange [0:2*pi]
set iso 100             # 值越大取点越多
set view equal xyz      # 显示时 x, y, z 三个坐标标度一致
unset tics              # 不显示坐标轴
unset border            # 不显示边框
splot sin(u)*cos(v), sin(u)*sin(v), cos(u) w pm3d

效果如下:
image

3. gnuplot 绘制低阶球谐函数

我们要对低阶球谐函数,制作 (R( heta, phi) = Y^m_l( heta, phi)) 的图。以下是 Griffiths 的书上的表格。
image
显然,(Y^0_0)就是一个球体。下面挨个画,苦哈哈。注意,为了方便,我略去了所有归一化因子。

3.1 (Y^m_1)

gnuplot 脚本如下,不仅做球谐函数的图,也做了一个球形网格,作为相对位置参考。

set param       # 参数图模式
set pm3d depthorder     # 前面的色块挡住后面的,形成立体感
set pm3d lighting primary 0.5 specular 0.6
set palette rgbformulae 7, 5, 15
set urange [0:pi]       # 设置 u, v 的取值范围
set vrange [0:2*pi]
set iso 50              # 值越大取点越多
set view equal xyz      # 显示时 x, y, z 三个坐标标度一致
unset tics              # 不显示坐标轴
unset border            # 不显示边框
set dummy u,v
set multiplot layout 1,2
# subplot
r(u) = cos(u)
splot r(u)*sin(u)*cos(v), r(u)*sin(u)*sin(v), r(u)*cos(u) w pm3d t "R({/Symbol q },{/Symbol j})=Y^0_1 ( {/Symbol q },{/Symbol j} )", sin(u)*cos(v), sin(u)*sin(v), cos(u) t "Sphere for reference"
# subplot
set urange [0:pi]
r(u) = sqrt(0.5)*sin(u)
splot r(u)*sin(u)*cos(v), r(u)*sin(u)*sin(v), r(u)*cos(u) w pm3d t "R({/Symbol q },{/Symbol j})=Y^1_1 ( {/Symbol q },{/Symbol j} )", sin(u)*cos(v), sin(u)*sin(v), cos(u) t "Sphere for reference"
unset multiplot

彩图如下:
image
可见,(Y^0_1)是个球体,偏离球心一段距离,这可以推导出来,因为(Y^0_1 propto cos( heta) e^{imphi}),所以

[x = cos( heta) sin( heta) cos(phi) = frac{1}{2} sin(2 heta) cos(phi), \ y = cos( heta) cos( heta) sin(phi) = frac{1}{2} sin(2 heta) sin(phi), \ z = cos( heta) cos( heta) = frac{1}{2} cos(2 heta) + frac{1}{2}. ]

确实表征一个球体(每个点被描绘2次)沿 (z) 方向向上偏离。
但是 (|Y^{pm 1}_1 |) 的形状是一个 Donut,为何 Peter Ring 说也表征球体平移,我就不懂了。

3.2 (Y^{m}_2)

set param       # 参数图模式
set pm3d depthorder     # 前面的色块挡住后面的,形成立体感
set pm3d lighting primary 0.5 specular 0.6
set palette rgbformulae 7, 5, 15
set urange [0:pi]       # 设置 u, v 的取值范围
set vrange [0:2*pi]
set iso 50              # 值越大取点越多
set view equal xyz      # 显示时 x, y, z 三个坐标标度一致
unset tics              # 不显示坐标轴
unset border            # 不显示边框
set dummy u,v
set multiplot layout 1,3
# subplot
r(u) = 0.5*(3*cos(u)*cos(u)-1)
splot r(u)*sin(u)*cos(v), r(u)*sin(u)*sin(v), r(u)*cos(u) w pm3d t "R({/Symbol q },{/Symbol j}): Y^0_2 ( {/Symbol q },{/Symbol j} )", sin(u)*cos(v), sin(u)*sin(v), cos(u) t "Sphere for reference"
# subplot
r(u) = 0.5*3*sin(u)*cos(u)
splot r(u)*sin(u)*cos(v), r(u)*sin(u)*sin(v), r(u)*cos(u) w pm3d t "R({/Symbol q },{/Symbol j}): Y^1_2 ( {/Symbol q },{/Symbol j} )", sin(u)*cos(v), sin(u)*sin(v), cos(u) t "Sphere for reference"
# subplot
r(u) = 0.25*3*sin(u)*sin(u)
splot r(u)*sin(u)*cos(v), r(u)*sin(u)*sin(v), r(u)*cos(u) w pm3d t "R({/Symbol q },{/Symbol j}): Y^2_2 ( {/Symbol q },{/Symbol j} )", sin(u)*cos(v), sin(u)*sin(v), cos(u) t "Sphere for reference"
unset multiplot

得到图:
image
没弄明白怎么和椭球体联系起来。

原文地址:https://www.cnblogs.com/luyi07/p/14713231.html