基于量子算法的碎纸片拼接复原问题

The Problem of Splicing and Recovery of Pieces of Paper Based on Quantum Algorithm

摘要: 本文借助于量子算法对2013年全国大学生数学建模竞赛B题的两种碎纸片(纵切和纵横切)的拼接方法进行新的探究。运用MATLAB软件，设计出了可以快速准确地将碎纸片进行复原的量子算法程序。在第一问仅有纵切的情况下，将原问题转化为旅行商问题，并借助量子蚁群算法并进行求解。在第二问中利用改进的量子聚类算法，较完美的完成了拼接任务，极大地降低了算法的时间复杂度和拼接的错误率。

Abstract: With the help of quantum computing method, this paper makes a new research on the splicing method of two kinds of shreded paper (longitudinal and transverse) of problem B in the 2013 National College Students’ mathematical modeling competition. We use MATLAB software to design a program that can quickly and accurately restore the broken pieces of paper. In the case of the first question, which is only vertical, we transform the original question into a traveling salesman problem and design a quantum ant colony algorithm to solve it. In the second question, the improved quantum clustering algorithm is used to complete the splicing task perfectly, which greatly reduces the time complexity of the algorithm and the error rate of stitching.

1. 引言

破碎文件及图像的自动拼接复原在物证复原、文献修复等领域有很多的应用，当碎片的数量很少时，可以采用手工拼接的方法完成。但是当碎纸的数量非常大时，如果仍然依靠手工完成，可能会耗费大量的时间、人力、物力，不但会对物证造成一定的损坏，有时还会错过最佳时机，例如在破案中如不及时拼接复原会耽误破案进度。很多国家已经对破碎文件的自动修复技术进行了大量的研究。因此，把计算机视觉和模式识别与碎片复原结合起来进行研究意义重大。然而在国内，针对这方面的研究还是比较少的。图像的拼接复原技术是一个新兴课题，具有很大的实用价值。运用先进的计算机技术作为工具，借助计算机高速运算这一优点，将极大减少复原所需时间，提高正确率。本文主要研究了规则的碎纸机破碎纸片的拼接方法，我们分别针对单面仅纵切和既纵切又横切规则碎片进行了拼接研究，并借助于量子算法找到了可行的拼接方案。

量子算法与经典算法最大的区别在于它能够凭借量子计算的一些特性(例如：量子并行性)来实现指数级加速。本文中，采用量子态来表示输入信息和输出信息，通过查询量子随机存取存储器来获取我们需要的纠缠量子，算法的核心思想是利用量子态的相干性，使客观所需的结果増强，同时使非所需的结果减弱。这样在测量的时候客观所需的结果会以相当高的概率出现。

2. 单面仅纵切碎纸片拼接

2.1. 数据预处理

借用MATLAB中的imread函数先将图片数据导入MATLAB软件并进行读取。实验发现，使用经过二值化处理的数据与没处理过的数据相比复原率要高，因此，我们将数据进行二值化处理，即将图片上的灰度值置为0和255，使图片仅显现黑白两种颜色 [1] [2]。

2.2. 模型建立

求解所有图片的排列顺序，其实是一个最优化问题，于是我们将碎纸片拼接问题转化为经典的旅行商问题来进行处理 [3]。

旅行商问题是经典的组合优化中的一个NP难题，可以将它描述为：一个旅行商从指定的一个城市出发，游玩所有城市并且每个城市只能游玩一次，每个城市之间的距离是不同的，求旅行商游玩路线的最短距离。现在国内外学者提出了多种方法来解决这一问题，比如遗传算法，蚁群算法等，本文将在前人的基础上介绍一种改进的量子蚁群算法来求碎纸片拼接问题的最优解。

碎纸片之间的距离定义为 $I=\frac{a}{a+b+c}$，a为两个向量1-1配对的总数，b为两个向量不配对的总数，

c为两个向量0-0配对的总数。

将每个碎纸片看作一座城市，定义城市之间的距离D为：一条碎纸片右侧和另一条碎纸片左侧的距离。通过计算所有城市之间的距离，可以得到相关的距离矩阵(值得注意的是，城市A到城市B的距离与城市B到城市A的距离是不同的，因为它们顺序相反时对应的是碎纸片不同侧的距离)。由于文档最左侧和文档最右侧的碎纸条都为空白，所以可以形成一个哈密顿圈，由此转化为求最优解的旅行商问题。

2.3. 算法设计

解决旅行商问题的算法主要有神经网络算法、模拟退火算法和遗传算法等，但为了提高运算的速度和准确性，这里我们介绍一种量子蚁群算法 [4]。

蚁群算法是由意大利学者Dorigo、Maniezzo等人于20世纪90年代提出的 [5] [6]。他们从蚂蚁寻找食物的过程中受到了启发，蚂蚁在其经过的路径上会留下信息素，其他蚂蚁可以识别这一物质，在蚂蚁寻找食物过程中，每只蚂蚁到达食物的路径上，都会留下信息素，每条路径对应了一种“可行解”，蚂蚁将食物搬回蚁穴后还会回去继续搬食物，所以路径短的蚂蚁用的时间少，重复频率快，单位时间内往返的蚂蚁就越多，路上留下的信息素也就会越多 [7]。同时，蚂蚁会沿着信息素浓度较高的路径行走，由此形成一种正反馈机制，当蚂蚁数量非常大时，一段时间后蚂蚁会集中到最佳路径上，对应问题中的最优解。

由于量子算法具有并行性，在运算时会大大提高运算速度，下文将用量子思想对蚁群算法进行改进。我们认为蚂蚁释放的信息素是在蚂蚁当前位置留下的，蚂蚁的位置由一组量子比特的概率幅表示，蚂蚁的移动用量子旋转门实现，同时引入量子变异操作 [8]。设蚂蚁m携带的一组量子比特为

$X_m=\left(\begin{array}{ccccc}\cos {\phi }_{m1}& \cos {\phi }_{m2}& \cos {\phi }_{m3}& \cdots & \cos {\phi }_{mn}\\ \sin {\phi }_{m2}& \sin {\phi }_{m2}& \sin {\phi }_{m3}& \cdots & \sin {\phi }_{mn}\end{array}\right)$，则蚂蚁占据的位置分别为 $P_{mx}=\left(\cos {\phi }_{m1},\cdots ,\cos {\phi }_{mn}\right)$ 和

$P_{my}=\left(\sin {\phi }_{m1},\cdots ,\sin {\phi }_{mn}\right)$ [9] [10] [11]。蚂蚁的位置对应问题的可能解，我们求得的最优位置便对应着问题的最优解。虽然蚂蚁数量不变，但借助于量子并行性，搜索空间加倍了，从而提高了搜索速度 [12]。

设蚂蚁在r处的信息素为f(r)，初始时全部设置为零，可见度g(r)我们用信息的可检测性来表示，即定义为实际检测强度与原有信息素强度的比值。

(1) 定义蚂蚁m由 $r_s$ 到 $r_t$ 位置的移动规则为

$r_t=\{\begin{array}{l}\arg \max \left\{{\left[f\left(r_s\right)\right]}^{\alpha }{\left[f\left(r_t\right)\right]}^{\beta }\right\}\left(r_s\in R\right)q\le q_0\\ {r^{\prime }}_{s}q>q_0\end{array}$

这里q为介于0到1之间的随机数， $q_0$ 为0到1之间的一个常数，R为蚂蚁位置的集合， ${r^{\prime }}_s$ 为目标位置。

$P\left(r_s\right)=\frac{{\left[f\left(r_s\right)\right]}^{\alpha }{\left[f\left(r_t\right)\right]}^{\beta }}{{{\displaystyle \sum }}_{r_s,r_t\in R}{\left[f\left(r_s\right)\right]}^{\alpha }{\left[f\left(r_t\right)\right]}^{\beta }}$

其中 $\alpha$ 和 $\beta$ 代表着信息素强度和它的可检测性对移动的影响比重，可由实际情况定义。

(2) 接下来我们用量子旋转门来改变量子比特的相位，从而实现蚂蚁位置的改变。蚂蚁的量子比特为

$X_r=\left(\begin{array}{ccccc}\cos {\phi }_{r1}& \cos {\phi }_{r2}& \cos {\phi }_{r3}& \cdots & \cos {\phi }_{rn}\\ \sin {\phi }_{r2}& \sin {\phi }_{r2}& \sin {\phi }_{r3}& \cdots & \sin {\phi }_{rn}\end{array}\right)$

$X_s=\left(\begin{array}{ccccc}\cos {\phi }_{s1}& \cos {\phi }_{s2}& \cos {\phi }_{s3}& \cdots & \cos {\phi }_{sn}\\ \sin {\phi }_{s1}& \sin {\phi }_{s2}& \sin {\phi }_{s3}& \cdots & \sin {\phi }_{sn}\end{array}\right)$

量子旋转门转角的大小其实是十分重要的，我们要尽可能使 $X_r$ 旋转后趋近于 $X_s$。这里我们取0.001π到0.05π之间。记

$A_i=\left|\begin{array}{cc}\cos {\phi }_{ri}& \cos {\phi }_{si}\\ \sin {\phi }_{ri}& \sin {\phi }_{si}\end{array}\right|$

转角由下式确定：

$\Delta {\theta }_i=-\mathrm{sgn}\left(A_i\right){\theta }_0{\text{e}}^{-t}$

其中 ${\theta }_0$ 为迭代初值，t为优化步数，对蚂蚁进行量子旋转门操作后，蚂蚁移动到 $r_t$ 位置。

(3) 为了使蚂蚁有新的空间位置，我们用量子非门对蚂蚁进行变异处理，方法如下

$\left[\begin{array}{cc}0& 1\\ 1& 0\end{array}\right]\left[\begin{array}{c}{\alpha }_{si}\\ {\beta }_{si}\end{array}\right]=\left[\begin{array}{c}{\beta }_{si}\\ {\alpha }_{si}\end{array}\right]$

(4) 设蚂蚁的前一个位置为 $r_q$，当前位置为 $r_s$，移动后的位置为 $r_t$，可得出更新规则为

$f\left(r_t\right)=f\left(r_s\right)+\mathrm{sgn}\left(\Delta fit\right)\times {\left|\Delta fit\right|}^{\alpha }$

$\Delta fit=fit\left(r_t\right)-fit(\; r\; s\; )$

$g\left(r_t\right)=g\left(r_s\right)+\mathrm{sgn}\left(\Delta \partial fit\right)\times {\left|\Delta \partial fit\right|}^{\alpha }$

$\Delta \partial fit=\underset{1\le i\le n}{\max }\frac{\partial fit}{\partial r_{ti}}-\underset{1\le i\le n}{\max }\frac{\partial fit}{\partial r_{si}}$

$fit$ 为适应度函数，定义为 $fit\left(x\right)=A_{\max }-f\left(x\right)$，其中 $A_{\max }$ 可为输入的恰当值，也可以为优化过程的最大数值。

当所有蚂蚁进行一次循环后，按照下面的公式更新信息素

$f\left(r_u\right)=\{\begin{array}{l}\left(1-\rho \right)f{r}_u+\rho fit\left(r_u\right)r_u=r^{\prime }\\ \left(1-\rho \right)f{r}_u{r}_u\ne r^{\prime }\end{array}$

其中， $\left(1-\rho \right)$ 为信息素挥发系数，大小在0到1之间， $r^{\prime }$ 为最优解。

(5) 算法实现步骤

步骤1：将每一只蚂蚁的信息素强度和可见度进行初始化，赋以相同的值。变异概率设为 $p_n$，

最大迭代次数为Max，当前迭代次数为0，随机生成初始群体为

$p_i=\left(\begin{array}{ccccc}\cos {\phi }_{i1}& \cos {\phi }_{i2}& \cos {\phi }_{i3}& \cdots & \cos {\phi }_{in}\\ \sin {\phi }_{i1}& \sin {\phi }_{i2}& \sin {\phi }_{i3}& \cdots & \sin {\phi }_{in}\end{array}\right)$

其中 ${\phi }_{ij}=2\pi \times q$， $i=1,2,3,\cdots ,m$ ； $j=1,2,3,\cdots ,n$ ；m是种群规模，n是量子位数。

步骤2：对种群中的每个蚂蚁按其移动规则和转移概率用量子旋转门进行移动，随机选取蚂蚁用量子非门进行变异操作。

步骤3：按照更新规则对信息素强度和可见度进行更新。

步骤4：迭代次数加一，当t < Max时，重复b；当t > Max时，输出结果。

2.4. 仿真结果及分析

借用量子改进蚁群算法，按照以上步骤编写算法程序，并在计算机上进行仿真运算。结果显示：成功拼接了碎纸片，完美地还原了原图片，且无需人工干预，得到了碎纸片的拼接序列(如表1所示)。

Table 1. Result sequence table of Appendix 1

表1. 附件1结果顺序表

还原结果图：(如图1所示)。

Figure 1. The restored picture of Annex 1

图1. 附件一还原的图片

量子蚁群算法由于利用了量子并行性，其运算速度得到了很大提高，本例中只有19条碎纸片，当碎纸片数量达到上千上万甚至更大后，量子算法的优越性将会更好的体现。考虑到改进的量子蚁群算法具有运算速度快，正确率高的优点，还可以把该算法运用到其他优化问题的求解中。例如银行在破碎现金的复原中，多数拼接方法因准确率问题大都需要人工辅助完成，而改进的量子蚁群算法经过测试后具有较高的准确性，可直接完成拼接任务。

3. 纵横切纸片的拼接

为了不影响原图片且方便图像处理，我们将原图片转化为二值图像，每个像素点都用0~255之间的一个灰度值表示，每个图片就可以分别表示为一个矩阵。我们首先对碎纸片进行聚类，找到碎纸片所在的行，然后在行内借助于问题一中的方法对碎纸片进行拼接。

设a为文字上方距离上边缘的距离，b为文字最左侧的位置与左边缘的距离，c为文字最下方的位置与下边缘的距离，d为文字最下方与下边缘的距离，如果在边缘检测到了像素点，那就令对应的距离为0，设矩阵 $e={\left[a,b,c,d\right]}^{\text{T}}$。

首先，把每个样品看作一类，并将每两个样品间的相似性转换为碎纸片的距离。其中 $c_{ij}$ 表示变量i和变量j之间的相关系数，看作类与类之间的距离，然后将距离最近的两类合成新的一类，每次减少一类，重新进行最近种类的合并，直至所有的变量合并成一类。

3.1. 算法设计

我们引入一种新型的聚类方法，是在传统的K-means算法基础上改进得来的 [13]。

(1) 数据的量子比特的表示

空间中的任意两个点a和b，可以把点a和b表示成向量，接下来将向量a和b用量子态表示，

$a=\left|a\right||a\rangle$， $b=\left|b\right||b\rangle$

接下来我们讨论二维情况

$|a\rangle =x_1|0\rangle +x_2|1\rangle$， $|b\rangle =y_1|0\rangle +y_2|1\rangle$

其中 $|0\rangle$， $|1\rangle$ 为计算基， $|a\rangle$ 、 $|b\rangle$ 表示计算后的概率，并且对它们进行了归一化。

(2) 用量子方法，描述碎纸片间的距离，定义两个态之间的距离为：

$C=\left|a-b\right|=\sqrt{\left(\left|a\right|\langle a|-\left|b\right|\langle b|\right)\left(\left|a\right||a\rangle -\left|b\right||b\rangle \right)}$

添加一个辅助量子比特来构建分类点a到聚类中心b的量子纠缠态，并测量这个纠缠态，从而得到概率和欧氏距离的关系。

3.2. 实现步骤

步骤1：根据最左边碎纸片边距要大的规律选择出最左边的11个碎纸片

$E=\left\{{\left(\left|b_1\right||b_1\rangle \right)}^0,{\left(\left|b_2\right||b_2\rangle \right)}^0,\cdots ,{\left(\left|b_n\right||b_n\rangle \right)}^0\right\}$

步骤2：对剩余的碎纸片归类到已有的11类中，其中归类方法是按最小距离原则进行，从而使得距离最近的碎纸片被分为一类。

步骤3：更新聚类中心，聚类中点为

${\left(\left|b_i\right||b_i\rangle \right)}^{t+1}=\frac{{\displaystyle \underset{j=1}{\overset{m_i}{\sum }}\left|a_j\right||a_j\rangle }}{m_i}$

其中 $m_i$ 表示中点数， $a_j$ 为距离元素。

步骤4：若 ${\left(\left|b_i\right||b_i\rangle \right)}^{t+1}={\left(\left|b_i\right||b_i\rangle \right)}^t$，结束并输出 ${\displaystyle \underset{j=1}{\overset{m_i}{\sum }}\left|a_j\right||a_j\rangle }$，否则转至步骤2。

这里我们选择初始聚类中心时，没有用K-means算法的传统方法随机选取，而是通过分析数据的特点即页边距得出了初始聚类中心，大大地减少了迭代次数，从而使运算速度和正确率都有了很大程度地提升。

3.3. 结果分析

按照上述算法的四个步骤，利用MATLAB软件进行编程，通过详细运算成功的还原了原图像。碎纸片的聚类结果和排列顺序：(如表2所示)。

Table 2. Result sequence table of Appendix 3

表2. 附件3结果顺序表

还原出的结果图：(如图2所示)。

Figure 2. The restored picture of Annex 3

图2. 附件三还原的结果图

事实证明，改进的量子K-means算法的聚类速度和准确率都得到了大幅提升，多次试验后，需要手工辅助的次数也比传统方法少了很多 [14]。此次实验仅有209个元素，当碎纸片被分割的更细小时，量子算法的优越性将会更好的展现。

4. 未来展望

由于量子算法并行性的特点，相对于传统算法具有运算速度快、准确率高等优势 [15]。现如今量子算法只能在传统计算机上局部实现，未来如果能研制出量子计算机，量子算法的速度将会大大提升，而且量子算法在计算量较大时，它的优势体现的更明显。在当今大数据时代，它的重要性不言而喻，量子技术和量子计算的应用有着十分广阔的前景 [16]。

基金项目

2021年山东省大学生创新训练项目(项目编号202110429213)，2020年山东省大学生创新训练项目(项目编号S202010429197)和2020年山东省本科教学改革研究重点课题(项目编号Z2020045)和资助。

NOTES

^*通讯作者。