计算化学程序的实现：粒子数表象下波函数的表示

占据数矢量的表示

我们现在来考虑程序的实现。需要说明的是，这种计算并不是很轻松，所以一些简单好用的编程语言是不适合拿来编写量子化学的计算程序的。基本上，可以选择的只有 Fortran 以及 C/C++。在此文中使用 C++ 来实现，并且我会用到一些不常见的手段来提高运行速度。如果你觉得阅读这部分有困难，那么可以读完此部分之后再去读一些关于 C/C++ 语言的位运算之类的文章。

程序的第一步从构建波函数开始，而波函数则是由一个或多个占据数矢量构成。既然每个占据数矢量都可以拆分成 $\alpha$ 和 $\beta$ 两部分，每部分都是一个由“占据”或者“非占据”两种状态组成的串，串上的每一个位置表示一个自旋轨道。那么很自然地，这种形式可以使用由 0 和 1 组成的序列来描述。但是应该使用怎样的数据类型呢？

最直接的想法是，既然每个位置的状态只有 0 和 1 两种，那么使用一个 bit（位） 来表示一个位置是最好的选择。不过计算机本身进行操作的最小内存单位是 byte（字节），等于 8 个 bit。计算机中所有的数据类型都至少是 1 byte 大小，即使是布尔类型也不例外。并且，任意长度的 bit 串是不可行的，这里存在内存对齐的问题。目前直接实现的 bit 串主要有以下两个：

• std::bitset<n>，需要在编译时就指定其大小。除非你想一遍遍地对不同的体系重新修改和编译程序，否则这并不是一个好的解决办法。
• std::vector<bool>，对于某些实现而言，这个类型是 bit 数组而非布尔类型的数组。但是此类型的访问速度很慢，无论是什么情况都不应该使用（除非内存占用确实很重要，重要到可以让你放弃性能）。

使用 bit 串存在的另一个问题是如何对其进行循环遍历，或者，至少需要一种方法能够高效地判断给定位置的状态。

另一个方法是使用数组，事实上很多现有的计算程序就是采用这种方式。数组形式的好处是可以方便直观地进行循环遍历。考虑到 int 类型的大小是 4 byte（32 位整数，对于 x86 以及 x86_64 处理器而言这是默认的情况），而实际需要表示的只有 0 和 1 两个数，所以即使真的要使用数组，也应该使用 int8_t（8 位整数）或者 char 这种最小的整数类型。但无论使用什么数据类型，数组形式都会浪费大量的内存空间，不适用于行列式数目非常多的情况。

这里我们使用一种折衷的实现方式：使用 64 位无符号整数 uint64_t 类型来表示一个长度为 64 的 0/1 序列。这种实现方式与文章 https://arxiv.org/abs/1311.6244 中的实现方式类似。对于有 norb 个轨道（此处指的是空间轨道，即 alpha 自旋轨道与 beta 自旋轨道都各有 norb 个）组成的体系，如果 norb 不大于 64，则需要 2 个 uint64_t 整数来表示这个占据数矢量，即 alpha 与 beta 各一个；否则，若 norb 大于 64，则需要多个 uint64_t 整数。简单来说，因为每 64 个轨道为一组，我们需要 nset 组 alpha 自旋轨道才能完全表示总共 norb 个 alpha 自旋轨道，并且有

$$nset = (norb + 63) / 64$$
而总的 uint64_t 整数的个数为 $2 \times nset$。这样我们就可以定义如下的类：

class SlaterDet
{
public:
    // 省略的代码 ...
private:
    static int _nset; // 此处为静态成员变量
    std::uint64_t* _strs; // 注意此处数组的长度为 2*_nset
};

注意这里只使用了一个总长度为 2 * nset 的数组，其中前 nset 个是 alpha 序列，后 nset 个是 beta 序列。每个数都有 64 个 bit，表示 64 个自旋轨道，并且 1 表示有电子占据而 0 表示无电子占据。考虑到 x86 处理器是 little-endian，这里的每一个 uint64_t 数也采用从低位到高位逐位递进的方式。注意此处的“逐位”，我们的操作对象不再是 byte 而是 bit。例如，对序列 2ud0，其中 2 表示 alpha/beta 双占据，u 表示 alpha（或 spin up）占据，d 表示 beta（或 spin down）占据，0 表示此空间轨道无电子占据，这样其 alpha 序列为 1100 而 beta 序列为 1010，在内存中其形式为

alpha: ...0011 # 索引为 [63, 62, ..., 3, 2, 1, 0]，... 表示 60 个零
beta: ...0101 # 索引为 [63, 62, ..., 3, 2, 1, 0]，... 表示 60 个零

简单来说就是与 2ud0 的书写方式左右颠倒，索引值是从左到右递减的，最右侧的轨道是第一个轨道。对于 $nset > 1$ 的情况也是类似的，只不过将整个序列切成了多个分块。例如，对于 $norb = 150$，有 $nset = 3$，此时索引值为

alpha[63, ..., 1, 0], alpha[127, ..., 64], alpha[191, ...,