Siesta 这个免费的第一性原理计算程序总体说是不错的，但手册做得挺简陋。关于能带计算，只给了一个面心立方的例子，网上其他资源包括教程里面也全都是面心立方。问题是，一旦体系不是面心立方的，该怎么计算能带呢？除了面心立方，别的都不会了怎么行。

不妨考虑一个简单正交的原胞，三个方向上的晶格常数分别为 Lx, Ly, Lz。Siesta 要计算能带的时候必须设置一个 LatticeConstant 参数。如果不设这个参数，也能计算，但结果中就不会有能带结构，这一点真挺坑人的。我建议对于非立方的晶系，为要计算能带，一律设置：

LatticeConstant 1 Ang

K点的选取也很微妙。官方说有两种 K 点的格式：

BandLinesScale pi/a

或者

BandLinesScale ReciprocalLatticeVectors

其中前者是默认的，后者也可选。但我对非立方晶系用后者算出的结果不对，网上也有人说后者“never works”。好吧，这样不如就用默认的 pi/a。

到底什么是 pi/a 呢？pi 就是圆周率，但 a 是什么呢？对于正交晶系，有三个晶格常数，到底哪个是 a 呢？实践证明，a 既不是 Lx 也不是 Ly 也不是 Lz。a 就是你设置的 LatticeConstant，在我们这个例子中就是 1 埃。

然后就是开动脑力，计算想要的 K 点在倒格子中的绝对坐标。以简单正交的 X 点为例。X 方向倒格子基矢的长度是 2pi/Lx。由于 -X $\rightarrow$ Gamma $\rightarrow$ X 就构成了一个周期，X 点的绝对坐标为（pi/Lx, 0, 0）。然而我们约定在 Siesta 输入中以 pi/a 为单位，就是以 pi 为单位，所以计算 X $\rightarrow " style="font-family:'comic sans ms';$ Gamma 方向上的能带结构，就写为：

%block BandLines

1      1/Lx  0.0  0.0    X

100  0.0    0.0  0.0    Gamma

%endblock BandLines

其中 1/Lx 要换成具体数值。从这里也可以看出 LatticeConstant 设成 1 Ang 的好处，就是说 pi/a 简单地等于 pi，运算方便。

最后有必要跟 Siesta 官方手册上的面心立方格子的例子做对比。在那个例子中，X 点的坐标是（2, 0, 0），为什么不是（1, 0, 0）呢？为什么 X 的坐标就可以达到 2pi/a 呢？岂不是 X 方向布里渊区跨度达到 4pi/a？这个问题应该从以下两个角度理解。首先，在计算面心立方晶体的能带时，我们的 a 设成了结晶学原胞的晶格常数，但输入的却是一个固体物理学原胞，其体积是结晶学原胞的 1/4。因此，布里渊区体积理应是预想中的 4 倍大才对。其次，面心立方的布里渊区是个截角八面体，其体积只是相同跨度立方体的一半。最终结论是，每个方向上的布里渊区的尺度都扩大为 2 倍，体积扩大为 8 倍，但因截角八面体的形状问题，布里渊区的体积的确只是预想中的 4 倍。相比较而言，简单正交原胞的布里渊区是个长方体，其体积也就实打实地等于三边乘积，所以 X 点坐标必然是 pi/Lx 而不是 2pi/Lx。