WaveFunctionCollapse

什么是 WFC？

WFC 中文全称为波函数坍塌。这是是个量子力学中的概念，只不过我们在计算中借用了这一概念。相信我们听说过“薛定谔的猫”的理论，一个物体在观测前拥有多种状态，但在观测后便固定成为一种状态。 通常这种算法用来实现程序化建筑或者地图生成。

以数独来举例，我们如何生成一个完备的数独表格，也可以用 WFC 的方式看待。 该算法中的约束，便是数字的选取规则。这里进行一下简单复习，正确的数独表格最终形成每一行、每一列以及每三乘三网格区域中都必须包含 1 到 9 的所有数字，同时相同的数字不会出现在同一行、同一列以及同一个三乘三的网格区域。 当我们观察一个空白的数独表格时，每个单元格都拥有全部的可能性，我们也可以称它们都处于“叠加态”，也就是说当不存在任何数时每个单元格都可以在 1-9 中取任意值。

当任意的单元格被“观测”，或者说确定取值“坍塌”成唯一解，根据约束将会影响到该行、该列以及该网格的取值范围，“传导”至受影响的单元格，把该数字从“叠加态”中移除。

在确定所有的初始限制条件后，“坍塌”和“传导”将会清理很多的“叠加态”。

然后我们优先找到叠加可能性最少的单元格，也就是找到“熵”最小的单元格，并使其坍塌成唯一解，并把对应的影响传递至其他单元格。优先处理最小熵单元格的原因是，只需要考虑更小的可能性，做出错误的可能也会越低。持续不断的坍塌和传导，直到所有单元格都坍塌成唯一解。实际上也会出现无法坍塌的错误情况，只不过概率较低，碰到后重新来一遍就好。

信息熵

熵可以简单的将其理解成情况的复杂程度。如果并不想了解更多的理论，以下知识点可以不用深入理解。

理解信息熵之前先要理解信息量。对时间的度量是时间量，单位是秒，那么对信息的度量是信息量单位是比特。当我们考虑一个离散的随机变量 $x$ 时，关于这个变量的概率事件被称为信息（如果不理解没关系，就像时间，我们知道什么叫时间增加减少，但是仍然不知道如何理解时间这个词语），信息的变化度量称作信息量。

当我们说某个事件携带了一定的信息量，我们实际上是在说这个事件的发生减少了我们的不确定性。信息量就是对这种减少的量化。例如，如果我们知道某个非常不可能发生的事件发生了，这个信息会极大地改变我们对该事件的看法，因此它携带了大量的信息。

信息的大小跟随机事件的概率有关。越小概率的事情发生了产生的信息量越大。 假设两个不相关事件，同时发生的信息和表示为

$$h(x+y) = h(x) + h(y).$$

由于是两个不相关事件，概率计算为

$$p(x+y) = p(x) \times p(y).$$

概率和信息量之间的关系和 $\log$ 在数学中的计算公式相似，推导得出

$$h(x) = -\log_{n}{p(x)}.$$

1. 负号是为了保证计算值仍然为正（概率的值通常为 [0,1]，$\log$ 值小于等于 0）
2. 底数 n 是可以任意取值的。（物理上常用 e，数学计算中常用 10，计算机相关常用 2）

信息量度量的是一个具体事件发生了所带来的信息，而熵则是考虑该随机变量的所有可能取值，即所有可能发生事件所带来的信息量的期望。

$$H(X) = -\sum_{i=1}^{n} P(x_i) \log_2 P(x_i)$$

该公式对比期望值公式发现很相似（概率和值的相乘后累加得到的值通常被称为“期望值”）。期望值理解为概率的长期平均结果。那我们也可以将熵理解为信息量的长期平均结果，也就是信息量的随机平均取值。

2D 实现

波函数坍塌算法是针对每一个单元格的求解过程，每个单元格都包含了所有“可能状态”的叠加态，我们称之为波函数。每种可能状态都有自己的邻接规则，比如图形 A 只允许上部分邻接 BCD，左部分邻接 ABE...，这便是该算法的约束。坍塌单元格并以相邻单元格进行传导，使用邻接规则将受影响的状态从叠加态中移除，也就是缩小邻居的值域。选取拥有最少值域的单元格坍塌成唯一值，然后传递影响不断重复该过程，直到所有单元格拥有唯一解。

这有个很好的 Demo 可供使用帮助理解。 Wave Function Collapse - Mixed Initiative Demo by Martin Donald (itch.io)

Martin Donald 提出一种“插槽系统”，使用某种标记标识所有的边界，通过标识的匹配来去除无效的邻居。

该网站使用边沿上固定三个点的颜色标记很好的展示了这一方法。 Wave - by Oskar Stålberg (oskarstalberg.com) 而它们的邻接关系是由这三个固定点的颜色是否匹配进行判定的。

无限生成城市的案例也是对该方案的支持。 使用波函数坍塌算法的无限程序生成城市 |玛丽安的博客 (marian42.de)

为了验证单元之间彼此彼邻接关系的有效性，每个模块都包含了 6 个邻接列表

3D 实现

从二维扩展到三维空间的实现是很容易的，相比于二维平面三维空间下只不过是将单元格变成了单元快，4 个邻接单元格增加到了 6 个邻接块。

Tile[][] neighbors = new Tile[6][]  
{  
    new Tile[] { tile1, tile2, tile3 .. },  //px
    new Tile[] { tile1, tile2, tile3 .. },  //nx
    new Tile[] { tile1, tile2, tile3 .. },  //py
    new Tile[] { tile1, tile2, tile3 .. },  //ny
    new Tile[] { tile1, tile2, tile3 .. },  //pz
    new Tile[] { tile1, tile2, tile3 .. },  //nz
};

然后通过事先检输入好的简单规则判断模块之间的连接是否符合要求，最后确定哪些可以拼接在一起。

首先要做的是循环遍历集合中的每个模块，并存储 tile 在六个边界上的每个顶点信息存储并记录这些边界配置通过一个字典来组织标记以及对应信息。还需为“对称”和“不对称”的插槽做一些特殊标记。给“对称”的插槽标记上“s”，无论方向如何“对称”的插槽始终与自身匹配。而不对称的插槽只有放置在自身的镜像旁才有效。通过镜像每个模块后，检测它们仍然相同，对于不同的我们将存储镜像版本和非镜像版本这两种插槽，镜像相同的其中一个将会标记为“f”。底部和顶部插槽使用 0、1、2、3 来标记它们的旋转方向。在垂直方向上，如果相邻插槽的标记以及旋转方向都能对应，就判定连接有效。不存在顶点的面上我们标记插槽为“-1”。

“原型”是指我们所用模块的元数据。包含了所用模型、模型旋转以及六个方向的有效邻接列表。模型旋转很好的解决了为每个模型的旋转方向都单独导出模型的问题。

public class Prototype
{
	public string model = "mesh.obj";
	public int rotation = 0;
	public string[] sockets = new string[6]
	{
		"0", //posX
		"1s", //negX
		"1s", //posY
		"0f", //negY
		"-1", //posZ
		"v0_0" //negZ
	};
}

构建

首先我们我们会为每个模块创建四个原型信息，每一个对应一个旋转方向。这个阶段原型只包含 3D 模型信息（可以时模型名称也可以是游戏对象）、旋转方向以及对应的插槽标识。最终的原型将基于其中的信息进行旋转。有一个特殊的原型需要注意，它不包含 3D 模型信息，六个面的插槽都对应“-1”，被称作“空白原型”。

比较原型

下一步我们将每个原型与其他原型进行 6 次比较，每个方向一次，使用我们的邻接规则来判断对应方向的插槽，如果符合则将对应原型作为有效邻居添加到列表中。

for (int i = 0; i < prototypes.Count; i++)
{
	var current = prototypes[i];
	for (int j = 0; j < prototypes.Count; j++)
	{
		var contrast = prototypes[j];
		for (int s = 0; s < 6; s++)
		{
			if (MatchCheck(current.sockets[s], contrast.sockets[matchingSurface[s]], s,matchingSurface[s]))
			{
				if (s == 0)
				{
					current.posX.Add(j);
				}else if (s == 1)
				{
					current.negX.Add(j);
				}else if (s == 2)
				{
					current.posY.Add(j);
				}else if (s == 3)
				{
					current.negY.Add(j);
				}else if (s == 4)
				{
					current.posZ.Add(j);
				}else if (s == 5)
				{
					current.negZ.Add(j);
				}
			}
		}
	}
}

public class Prototype
{
	public string model = "mesh.obj";
	public int rotation = 0;
	public string[] sockets = new string[6]
	{
		"0", //posX
		"1s", //negX
		"1s", //posY
		"0f", //negY
		"-1", //posZ
		"v0_0" //negZ
	};

	public string[] posX = new string[]{...};
	public string[] negX = new string[]{...};
	public string[] posY = new string[]{...};
	public string[] negY = new string[]{...};
	public string[] posZ = new string[]{...};
	public string[] negZ = new string[]{...};
}

插槽为 2 和 3 代表其是上下方向，垂直方向的判断只需断定插槽命名是否完全一样。分对称判别字符后缀"f"，对称以"s"结尾。需要对 Blank 空白模型特殊处理，仅让其水平方向可匹配“-1”插槽。

public static bool MatchCheck(string socketA, string socketB, int socketAI, int socketBI)
{
	//vertical
	if (socketAI == 2 | socketAI == 3)
	{
		if (string.Equals(socketA, socketB))
		{
			return true;
		}
	}
	//asymmetrical
	if (string.Equals(socketA,socketB+"f") | string.Equals(socketA + "f",socketB) )
	{
		return true;
	}
	//symmetrical
	if (socketA.EndsWith('s') & socketB.EndsWith('s'))
	{
		if (string.Equals(socketA,socketB))
		{
			return true;
		}
	}
	//empty
	if ((socketAI != 2 & socketAI != 3) & string.Equals(socketA,"-1") & string.Equals(socketB, "-1"))
	{
		return true;
	}
	
	return false;
}

加载原型数据

读取配置与预处理好的数据文件。

public static List<Prototype> LoadJson(string path)
{
	var json = System.IO.File.ReadAllText(path);
	List<Prototype> prototypes = JsonConvert.DeserializeObject<List<Prototype>>(json);
	return prototypes;
}

初始化波函数

使用原型数据填充每个网格单元，初始状态下都拥有全叠加态，在不断的坍塌和传递叠加态会陆续减少。

collapsedCount = size.z * size.y * size.x;
wave = new List<int>[collapsedCount];
int[] prototypeIndexes = Enumerable.Range(0, allPrototypes.Count).ToArray();

for (int i = 0; i < wave.Length; i++)
{
	wave[i] = prototypeIndexes.ToList();
}

主函数中创建 While 循环直到波函数坍塌完成再跳出循环。

while (!IsCollapsed())
{
	if (!Iterate())
	{
		...
	}
}

IsCollapsed 所有坍塌完成的状态监测方式最终是知道所有网格单元是否只具备唯一解。我们可以通过遍历网格查看它们可能的邻居数量得知，或者使用更节约性能的方法，在每次坍塌完成时修改未坍塌数量标记，直到为 0。

坍塌未完成程序迭代不止。迭代的过程是从网格中获取熵值最低的单元格，利用随机算法得到唯一解，向四周传递。

public bool Iterate()
{
	int coord = GetEntropyCoord();
	CollapseTo(coord);
	if (!Propagate(coord))
	{
		return false;
	}
	return true;
}

计算当前坍塌后方向上规则正确的可能原型，然后对比该方向上原型列表，从中剔除不该存在原型，并再次传递。现在该过程与树的深度遍历较为相似。

合理中的不合理

有的时候按照插槽标准生成的块规则上是合理的但在视觉上我们无法认同，即插槽命名规则并不是完美的。

如何规避这样不合理情况，我们有两种办法

使用过滤列表

手动的为每个原型配置增加 filter 过滤数据，该数据准确记录了六个方向上被过滤原型，从而在比较原型数据时将其排除在外。坍塌与传递时该数据便会在该原型的指定方向上消失。

{
        "model": "A",
        "rotation": 0,
        "sockets": [
            "1f",
            "-1",
            "-1",
            "v2_0",
            "-1",
            "1"
        ],
        "filter":[
            [...],
            [...],
            [...],
            [...],
            [...],
            [...]
        ],
        "posX": [...],
        "negX": [...],
        "posY": [...],
        "negY": [...],
        "posZ": [...],
        "negZ": [...]
}

使用权重熵

相比于过滤列表的使用，权重的配置更加简单。我们希望根据配置较低权重来以较小的概率生成不合理的情况。所以我们将配置文件字段设置如下，增加 weight 字段，初始为 0。

{
        "model": "A",
        "rotation": 0,
        "sockets": [
            "1f",
            "-1",
            "-1",
            "v2_0",
            "-1",
            "1"
        ],
        "weight":0,
        "posX": [...],
        "negX": [...],
        "posY": [...],
        "negY": [...],
        "posZ": [...],
        "negZ": [...]
}

在配置处理阶段，根据匹配出现频率获取权重，然而预先设置为 1 的权重我们不进行更改，让其在权重列表中保持最低值，以达到手动设置的目的。

public static void PrototypePreprocess(List<Prototype> prototypes)
{
	int[] weight = new int[prototypes.Count];
	for (int i = 0; i < prototypes.Count; i++)
	{
		var current = prototypes[i];
...
		if (MatchCheck(current.sockets[s], contrast.sockets[matchingSurface[s]], s,matchingSurface[s]))
		{
			weight[j]++;
...
		}
	}

	for (int i = 0; i < prototypes.Count; i++)
	{
		if (prototypes[i].weight != 1)
		{
			prototypes[i].weight = weight[i];
		}
	}
}

坍塌时根据权重计算随机值，有效的保证叠加状态越多时，较低权重的“态”越不可能被选中（当然选中机率仍然存在）。

/// <summary>
/// 坍塌唯一解
/// </summary>
/// <param name="coord">坍塌网格快</param>
public void CollapseTo(int coord)
{
	List<int> possiblePrototypes = wave[coord];
	for (int i = 0; i < weightCache.Length; i++)
	{
		weightCache[i] = 0;
	}
	float sumOfWeights = 0;
	for (int i = 0; i < possiblePrototypes.Count; i++)
	{
		int index = possiblePrototypes[i];
		weightCache[index] = allPrototypes[index].weight;
		sumOfWeights += weightCache[index];
	}
	int prototypeIndex = 0;
	float range = Random.Range(0f, sumOfWeights);
	for (int i = 0; i < weightCache.Length; i++)
	{
		range -= weightCache[i];
		if (range <= 0)
		{
			prototypeIndex = i;
			break;
		}
	}
	
	// int range = Random.Range(0, possiblePrototypes.Count);
	//坍塌完成
	possiblePrototypes.Clear();
	possiblePrototypes.Add(prototypeIndex);
	entropies[coord] = 0;
	collapsedCount--;
}

Random.Range(0, possiblePrototypes.Count) 代码展示的是完全随机的坍塌方式，使用权重的随机算法，影响原型选取。

熵计算：在信息论中，熵是一个衡量不确定性的指标。对于概率分布，熵的公式是 $H(X) = -\sum_{i=1}^{n} P(x_i) \log_2 P(x_i)$，其中 $P_{i}$ 是事件发生的概率。熵可以用来评估模型的不确定性或者信息量。

权重除了影响坍塌时的随机解，同时也可用于坍塌网格单元时选取策略。某原型的 权重/权重和 计算值代表该单元某原型的出现概率。$w$ 表示单元网格内可能原型的数组，使用 $w_i$ 表示下标为 $i$ 的原型权重，$S(w)$ 表示原型权重累加，即 $S(w) =\sum_{i=1}^{n}\limits{w_{i}}$。计算该单元格的熵公式如下:

$$\begin{array}{} \text{Entropy} &= -\sum_{i=1}^n \limits{\frac{w_{i}}{S(w)}log_{2}{\frac{w_{i}}{S(w)}}} \\ &= -\sum_{i=1}^n \limits{\frac{w_{i}}{S(w)} (\log_{2}{w_{i}} - \log_{2}{S(w)})}\\ &= \sum_{i=1}^{n}\limits{\frac{w_{i}}{S(w)} \log_{2}{S(w)}} - \sum_{i=1}^{n}\limits{\frac{w_{i}}{S(w)} \log_{2}{w_{i}}}\\ \end{array}$$

计算单元网格熵时，该单元 $S(w)$ 值始终不变，这意味着

$$\sum_{i=1}^{n}\limits{\frac{w_{i}}{S(w)} \log_{2}{S(w)}} = \frac{\log_{2}{S(w)}}{S(w)}\sum_{i=1}^{n}\limits{w_{i}} = \log_{2}{S(w)}.$$

所以得:

$$\text{Entropy} = \log_{2}{S(w)} - \frac{1}{S(w)}\sum_{i=1}^{n}\limits{w_{i}\log_{2}{w_{i}}}.$$

另 $L(w)$ 表示权重对数累加，即 $L(w)=\sum_{i=0}^{n}\limits{w_{i}\log_{2}{w_{i}}}$，最终公式表示为

$$E(w) = \log_{2}{S(w)} - \frac{1}{S(w)}L(w).$$

当波函数坍塌以及传递时，使用原先值的记录进行剔除从而避免重新累加计算。另 $w'$ 是原 $w$ 权重数组剔除 $\Delta{w}$ 元素后的子数组。

$$w' = w - \Delta{w}.$$

那么可得:

$$\begin{align} E(w') &= \log_{2}{S(w')} - \frac{1}{S(w')}L(w') \\ &=\log_{2}{(S(w)-S(\Delta{w}))} - \frac{1}{S(w) - S(\Delta{w})}(L(w) - L(\Delta{w})). \end{align}$$

首次坍塌前记录初始值 $S_{o},w_{o},L_{o}$

private int[] weights; //权重列表
private int sumOfWeights; //权重总和
private double[] weightLogWeights; //权重对数
private double sumOfWeightLogWeights; //权重对数总和
private double startingEntropy;//初始熵
private double sumOfLogWeght; //权重对数总和
private double sumOfLogWeights; //权重对数总和列表

private int[] sumsOfWeights; //单元格权重总和列表
private double[] sumsOfWeightLogWeights; //单元格权重对数总和列表
private double[] entropies; //单元格熵列表

weights = new int[allPrototypes.Count];
weightLogWeights = new double[allPrototypes.Count];
sumOfWeights = 0;
sumOfWeightLogWeights = 0.0;
for (int i = 0; i < weights.Length; i++)
{
	weights[i] = allPrototypes[i].weight;
	weightLogWeights[i] = weights[i] * Math.Log(weights[i]);
	sumOfWeights += weights[i];
	sumOfWeightLogWeights += weightLogWeights[i];
}
startingEntropy = Math.Log(sumOfWeights) - sumOfWeightLogWeights / sumOfWeights;

sumsOfWeights = new int[collapsedCount];
sumsOfWeightLogWeights = new double[collapsedCount];
entropies = new double[collapsedCount];

for (int i = 0; i < collapsedCount; i++)
{
	sumsOfWeights[i] = sumOfWeights;
	sumsOfWeightLogWeights[i] = sumOfWeightLogWeights;
	entropies[i] = startingEntropy;
}

当叠加态变化时，计算改变后的熵值。

/// <summary>
/// 熵变
/// </summary>
/// <param name="coord">变化单元网格</param>
/// <param name="prototypeIndex">剔除原型</param>
private void EntropieChange(int coord, int prototypeIndex)
{
	double sum = sumsOfWeights[coord] -= weights[prototypeIndex]; //计算权重总和
	sumsOfWeightLogWeights[coord] -= weightLogWeights[prototypeIndex]; 
	entropies[coord] = Math.Log(sum) - sumsOfWeightLogWeights[coord] / sum;
}

参考

mxgmn/WaveFunctionCollapse: Bitmap & tilemap generation from a single example with the help of ideas from quantum mechanics (github.com)

Superpositions, Sudoku, the Wave Function Collapse algorithm. - YouTube

使用波函数坍塌算法的无限程序生成城市 |玛丽安的博客 (marian42.de)