量子计算在增强应变梯度晶体塑性理论中的应用前景

《INTERNATIONAL JOURNAL OF PLASTICITY》：Prospect of quantum computing on Enhanced Strain Gradient Crystal Plasticity theory

　　阿米尔侯赛因·拉梅·朱伊巴里|莱昂·奇泽利
约瑟夫·斯特凡研究所，Jamova cesta 39号，卢布尔雅那，斯洛文尼亚SL-1000

**摘要**
本文介绍了一种新的增强型应变梯度晶体塑性（ESGCP）理论分支，该理论基于与每个滑移系上累积剪切应变梯度相关的二次能量贡献，并在热力学一致的框架内解释了晶体微观结构中滑移带和屈曲带的形成过程。结合最近提出的基于Nye张量的ESGCP公式，开发了一种新的微分算子来求解相应的非局部场方程（或高阶平衡方程）。在ESGCP理论的这两个分支中，高阶模量与受辐照晶体晶格在变形过程中的微观结构演变内在地耦合在一起。该框架被用来研究Hall-Petch（平均晶粒尺寸）效应，并系统地与经典应变梯度晶体塑性（CSGCP）模型进行了比较。结果表明，与CSGCP模型中晶粒尺寸效应随载荷增加而持续增强不同，ESGCP模型预测在低应变水平下晶粒尺寸敏感性增强，而在较高载荷水平下则逐渐减弱。

此外，还开发了一种基于量子傅里叶变换（QFT）的新型量子计算算法，用于在固定点迭代方案中求解经典线性动量平衡方程，而与ESGCP和CSGCP模型相关的非局部场方程则采用量子有限差分方法进行处理。实验表明，所提出的QFT方法实现了多对数级的计算加速，这对于高分辨率模拟受辐照材料具有重要意义，因为在核电站应用中，数值精度和计算效率对于可靠的结构完整性评估至关重要。

**引言**
当代计算科学正经历一场由量子计算快速发展及其与高性能和超级计算架构的整合所驱动的颠覆性变革。量子计算是一种基于非经典信息表示和数据处理的新型计算范式，它能够访问指数级大的解空间，并通过固有的并行量子操作从高度复杂的数据集中提取信息。通过利用量子力学原理（如叠加、纠缠和干涉），量子算法可以在特定类别的问题上显著提高计算效率。在经典计算中，算法加速通常通过增加硬件资源和并行化策略来实现；然而，对于许多问题，输入规模的线性增长会导致计算时间和内存需求的指数级增加。相比之下，对于某些问题类别，量子算法可以实现多项式级别的扩展，而经典算法则表现出指数级复杂性。经典计算依赖于经典比特的操作，这些比特在任何给定时间只能存在于0或1两种离散状态之一。量子计算则基于量子比特（qubit），它是量子信息的基本单位。与经典比特不同，量子比特可以存在于基态的叠加态中，在希尔伯特空间H2中表示为|0?和|1?的线性组合。在测量时，量子态根据玻恩规则不可逆地坍缩到这些基态之一，引入了量子计算特有的概率结果。

量子计算的动机可以追溯到20世纪80年代初，当时理查德·费曼（Feynman, 1982）和尤里·马宁（Manin, 1980）独立认识到经典计算机在模拟量子系统时效率极低，尤其是当相互作用粒子的数量超过几个时。从那时起，持续的理论和实验进展促使了几种具有证明优势的量子算法的发展。著名的算法包括格罗弗搜索算法（Grover, 1997）、肖尔整数分解算法（Shor, 1997）、用于求解线性方程组的哈罗-哈西迪姆-劳埃德（HHL）算法（Harrow et al., 2009）、变分量子本征值求解器（VQE）（McClean et al., 2016）以及用于哈密顿量模拟的量子算法（Low and Chuang, 2019）。格罗弗搜索算法可以通过振幅放大从大小为N的非结构化数据库中检索目标元素，只需O(N)次操作，而最优的经典算法则需要O(N)次操作。因此，格罗弗算法比经典搜索方法快两个数量级。肖尔算法可以在O(n3logn)时间内使用O(n2lognloglogn)个量子门分解n位整数N，而已知的最有效经典算法的运行时间呈亚指数级增长，为O(exp(649n(logn)23))。HHL算法可以在O(poly(logN,κ))时间内求解具有N个变量的线性方程组，其中κ表示系数矩阵的条件数，这比经典算法的O(Nκ)复杂度有了指数级的改进。HHL算法包括五个主要阶段：状态准备、量子相位估计、受控辅助旋转、逆量子相位估计和测量。变分量子本征值求解器采用混合量子-经典计算策略，通过迭代准备参数化的量子态|ψ(θi)〉，在量子硬件上评估哈密顿量的期望值，并使用经典优化算法更新参数θi直至收敛。这种方法特别适用于凝聚态和材料物理学问题，因为局部相互作用允许将哈密顿量分解为可以并行评估的算子之和。直接哈密顿量模拟的量子算法通过将时间演化分解为一系列子区间来处理薛定谔方程的数值解，在这些子区间内独立演化各个哈密顿分量。这种分解使得运行时间可以扩展为O(poly(n)t?)，其中n是量子比特的数量，t是演化时间，?是所需的模拟精度。相比之下，经典模拟的通用量子动力学的运行时间为O(exp(n)t?)，突显了量子模拟技术带来的指数级优势。

量子傅里叶变换（QFT）是量子计算中的一个基本原语，在加速量子硬件上的关键子程序中起着核心作用（Coppersmith, 2002, Sutor, 2024, Nielsen and Chuang, 2010）。QFT是几个关键量子算法（包括量子相位估计、HHL算法、肖尔分解算法和哈密顿量模拟算法）的关键构建块。特别是量子相位估计——HHL算法的核心组成部分——通常分为三个主要阶段：(i) 使用哈达玛门在时钟（控制）寄存器上准备均匀叠加态，(ii) 将本征相位信息印记到控制寄存器上的受控幺正操作，以及(iii) 应用逆量子傅里叶变换（IQFT），将累积的相位转换为二进制表示，从而使系数矩阵的特征值编码在时钟量子比特中。在肖尔算法中，QFT用于提取与模幂运算相关的周期，从而实现高效的整数分解。在量子哈密顿量模拟中，QFT用于在共轭表示之间切换（例如，位置和动量），其中变换便于在更容易实现的基态中演化哈密顿量（或其分量），随后通过IQFT返回到实空间（位置）表示。与经典快速傅里叶变换（FFT）相比，经典计算机上对N=2n个数据点应用FFT需要O(Nlog(N))时间，而量子计算机上的QFT可以在O((log(N))2)时间内实现。这种渐进式的改进表明，对于大型问题，只要数据以合适的量子编码形式提供，并且状态准备和读出成本不会主导整体复杂性，量子计算机可以比经典计算机更快地执行傅里叶变换。因此，QFT为基于傅里叶的经典数值方法和量子算法原语之间提供了自然的桥梁。

量子计算最近被探索作为固体力学中数值方法的使能技术，其潜力在于减少计算复杂性并加速某些类别的计算，超出了纯经典算法所能实现的程度。特别是，已经研究了量子算法在求解偏微分方程、边界值问题以及有限差分和有限元离散化中常见的大规模线性系统中的应用。尽管量子方法的实际性能取决于状态准备、噪声和读出开销，但其底层的算法扩展性激发了继续开发量子和量子-经典混合方法用于计算力学的动力。一些代表性的贡献展示了当前的研究进展。例如，曹等人（Cao et al., 2013）提出了一种用于泊松方程的量子算法，引入了一种将拉普拉斯算子嵌入哈密顿量并通过哈密顿量模拟技术求解的量子有限差分公式；随后，宫本和久保（Miyamoto and Kubo, 2021）报告了量子有限差分方法的进一步进展。在有限元方法（FEM）的背景下，Montanaro和Pallister（2016）开发了一种用于加速与边界值问题相关的线性系统求解的量子算法，灵感来自HHL算法，并在标准复杂性假设下展示了多项式级别的加速。将这一方向扩展到近期设备，Aroora等人（2025）使用线性和二次形状函数在噪声量子硬件上开发了一种用于稳态热方程的量子-FEM框架。除了直接的PDE求解器之外，还引入了用于结构设计的量子优化范式。例如，Sukulthanasorn等人（2025）研究了量子退火在拓扑优化和桁架及连续结构设计中的应用。还探索了量子启发的迭代方法用于FEM线性系统，如Raisuddin和De（2024）提出的量子松弛算法。在尺度桥接层面，Liu等人（2024）提出了一种量子代表性体积元素（RVE）求解器，并展示了其在泊松方程中的应用。最后，Lubasch等人（2020）报告了使用基于QFT的原语求解偏微分方程，包括不可压缩对流、热传导、各向同性声波传播和泊松方程。总体而言，这些研究表明，人们越来越感兴趣于利用量子算法原语——特别是与线性代数和傅里叶分析相关的那些——来加速与固体力学相关的计算工作流程。

在微观尺度上，金属材料由晶体晶格组成，其在外加载下的机械响应源于弹性变形和塑性变形机制的联合贡献（Tang et al., 2023, Dadhich and Alankar, 2022, Ahmadikia et al., 2021, Lindroos et al., 2022a, Agius et al., 2022, Lindroos et al., 2022b, Ren et al., 2021, Pai et al., 2022, Cao et al., 2025, White et al., 2025, Li et al., 2025, Wang et al., 2025, Zhang and Sumigawa, 2025, Toursangsaraki et al., 2025, Pai et al., 2025）。虽然弹性变形与可逆的晶格畸变相关，但晶体固体中的塑性变形主要通过位错沿特定晶体学方向和平面的滑移来实现，这些方向和平面通常被称为滑移系。晶体尺度上的塑性变形的异质性及其与晶体学取向的相互作用已在经典晶体塑性（CCP）理论框架内得到了广泛研究，该理论为模拟晶体材料中的各向异性塑性行为提供了坚实的基础（Ry? et al., 2022, You et al., 2023, Lindroos et al., 2022b, Ren et al., 2021, Pai et al., 2022, Cao et al., 2025, White et al., 2025, Li et al., 2025, Wang et al., 2025, Zhang and Sumigawa, 2025, Toursangsaraki et al., 2025, Pai et al., 2025）。塑性变形的异质性及其与晶体学取向的相互作用已在经典晶体塑性理论（CCP）框架内得到了广泛研究，该理论为模拟晶体材料中的各向异性塑性行为提供了坚实的基础（Ry? et al., 2022, You et al., 2023, Lindroos et al., 2022b, Ren et al., 2021, Pai et al., 2022, Cao et al., 2025, White et al., 2025, Li et al., 2025, Wang et al., 2025, Zhang and Sumigawa, 2025, Toursangsaraki et al., 2025, Pai et al., 2025）。然而，明确引入位错概念揭示了CCP框架的重要局限性。特别是，Nye（1953）提出的统计存储位错（SSDs）和几何必需位错（GNDs）之间的区别对于捕捉尺寸依赖的塑性效应（如Hall-Petch效应和其他尺度敏感现象）至关重要。这些观察结果激发了应变梯度晶体塑性（SGCP）理论的发展，该理论将塑性变形的梯度纳入考虑，以解释与GNDs相关的额外硬化。Gurtin（2002）在这个方向上做出了早期且有影响力的贡献，他提出了基于Nye张量的SGCP公式，建立了晶格曲率与GND密度之间的直接联系。后续的发展引入了替代的梯度度量方法，以解决基于Nye张量公式的特定缺陷。例如，Wulfinghoff和B?hlke（2012年）提出了一种基于总累积剪切应变的应变梯度晶体塑性模型，该模型在涉及应变局部化的问题中提供了更好的鲁棒性。后来，Lame Jouybari等人（2024年）扩展了这种方法，将累积剪切应变的贡献与各个滑移系结合起来，从而提高了滑移系级塑性的物理分辨率。一般来说，SGCP模型通过非局部场方程引入了一个内在的材料长度尺度或非局部参数，该参数控制着应变梯度对材料响应的影响。在大多数现有的公式中，这个长度尺度被假设在整个变形过程中是恒定的。然而，实验证据表明，内在长度尺度与不断演变的微观结构特征密切相关，不一定应该被视为一个固定的材料参数。来自微弯曲和微扭转实验的观察结果表明，应变梯度效应在变形的早期阶段最为明显，并且随着塑性流动的发展而逐渐减弱，这意味着微观结构的演变可以在更高的变形水平上减少尺寸效应。这种行为促使我们将长度尺度解释为一个反映潜在位错滑移变化的演变量。沿着这些思路，Zbib和Aifantis（1998年）为抛物线硬化和软化行为开发了依赖于梯度的本构关系，后来又将这一概念扩展到控制损伤局部化的损伤依赖长度尺度（Brepols等人，2020年）。此外，Dahlberg和Boasen（2019年）直接将内在长度尺度与位错密度联系起来，而Petryk和Stupkiewicz（2016年）将其与应变率梯度联系起来。还提出了其他微形态和梯度增强公式来调节应变局部化，包括Scherer等人（2019年）开发的基于累积剪切应变的微形态晶体塑性模型，以及Abatour和Forest（2024年）引入的具有饱和内部变量的应变梯度塑性框架。最近，Lame Jouybari等人（2025年、2026年）提出了一种基于Nye张量的增强应变梯度晶体塑性理论，在该理论中，内在长度尺度与变形过程中受辐照材料的微观结构演变明确耦合。总的来说，这些发展突显了微观结构演变在应变梯度塑性中的核心作用，并激发了在晶体塑性理论内制定物理上合理的、自适应长度尺度模型的持续努力。

自从芝加哥Pile-1反应堆在1942年12月达到临界运行状态以来，人们一直在努力提高核反应堆组件的效率和结构完整性。在反应堆运行期间，来自裂变链反应的快中子与结构内部相互作用，可能会产生辐照诱导的晶格缺陷，这些缺陷统称为辐照缺陷。这些缺陷通过创建和演变点缺陷、缺陷簇及相关损伤特征来改变微观结构，从而改变反应堆材料的机械响应和宏观行为。据报道的后果包括韧性和延展性的降低、屈服强度的增加以及加工硬化区域的收缩，所有这些都可能显著影响组件的可靠性和寿命。尽管在许多受辐照的合金中观察到明显的延展性损失，但对受辐照材料（包括不锈钢）的断裂表面分析表明，主导的断裂机制仍然可以是延展性的，尽管其特征是塑性变形的强烈微观局部化（El Shawish等人，2016年；Hesterberg等人，2019年；Lame Jouybari等人，2023年）。这种明显的悖论可以通过高度局部化的变形模式的出现来解释，即使在宏观塑性被强烈抑制的情况下，这些模式也会促进早期失效。与此解释一致的是，透射电子显微镜（TEM）观察证实了在塑性变形部分，辐照缺陷被滑移位错扫除，形成了狭窄的无缺陷通道——通常称为清晰通道（Thomas等人，2019年）。这种通道的发展将塑性滑移集中在微观结构的一小部分，从而增强了应变局部化并加速了延性损伤的积累。因此，清晰通道变形被广泛认为是受辐照合金过早延性失效的主要微观机制。

基于快速傅里叶变换（FFT）的方法用于解决异质介质中的经典平衡（柯西平衡）方程，首次在Moulinec和Suquet（1998年）的开创性工作中提出。在这种方法中，通过引入一个辅助的均匀参考弹性介质来重新表述异质材料，这使得控制方程可以重新表述为Lippmann–Schwinger积分方程。所得到的公式特别适合周期性微观结构，并且受益于傅里叶空间中卷积运算符的高效评估。由于其计算效率以及与成像获得的体素化微观结构的自然兼容性，基于FFT的方法已成为计算微观力学的典范。特别是，与标准有限元实现相比，后者的计算工作量可能随N2而增加，而基于FFT的求解器将主要成本降低到O(Nlog(N))，其中N表示离散代表域的元素或体素的数量。值得一提的是，已经提出了几种方法来降低超弹性材料有限元模拟的计算成本。值得注意的是，所谓的校正簇积分（E3C）超简化方法已被开发出来，以显著减少CPU时间（Wulfinghoff，2025年；Wulfinghoff，2026年）。此外，在微形态连续介质的背景下还探索了数据驱动的计算方法（Ulloa等人，2024年）。然而，原始FFT固定点方案的收敛行为被证明强烈依赖于辅助参考材料的选择以及异质介质中各相之间的对比度。高相位对比度和次优的参考选择可能导致收敛缓慢甚至数值发散，因此激发了后续为提高鲁棒性和效率所做的实质性努力。因此，提出了多种增强公式，包括基于极化的方案（Wicht等人，2021年）、Krylov子空间加速方法（Brisard和Dormieux，2010年）、配置方法（Zeman等人，2010年）以及傅里叶-伽辽金方案（Lucarini和Segurado，2019年）。总的来说，这些发展显著扩大了基于FFT的求解器的适用范围，特别是对于强异质材料和非线性本构行为。

本研究首次将量子傅里叶变换（QFT）方法引入经典晶体塑性框架和应变梯度晶体塑性（SGCP）理论中，以加速FFT方法求解器中出现的基于傅里叶的计算。具体来说，所提出的方法旨在使用QFT原语将经典FFT方法中傅里叶变换步骤的渐近运行时间从O(Nlog(N))降低到O((log(N))2，根据域分辨率产生多对数复杂性。这种加速对于高分辨率的多晶聚集体模拟特别有利，在三维分析中，体素数量N可能变得非常大。详细介绍了QFT方法，包括完整的算法工作流程以及所需的中间步骤和数学证明。此外，该研究还回顾了近期文献中的几个SGCP模型，这里组织为经典MicroSlip-SGCP（CMSlip）模型、经典MicroCurl-SGCP（CMCurl）模型和增强型MicroCurl-SGCP（EMCurl）模型。在这些框架的基础上，引入了一种新的ESGCP理论分支，即由每个滑移系相关的累积剪切应变驱动的增强型MicroSlip-SGCP（EMSlip）模型。此外，修改后的运算符被纳入CMCurl和EMCurl模型中；这些模型在本工作中分别称为CMCurl-Star和EMCurl-Star模型。对于每个考虑的SGCP模型，都在量子有限差分方法中推导并详细讨论了相应的QFT方法公式，确保了控制运算符的一致量子-经典映射。最后，通过二维和三维模拟受拉伸载荷的辐照多晶聚集体来演示所提出的QFT方法，从而展示了其在塑性变形中的适用性。需要强调的是，本研究中的所有量子算法都是使用Matlab（2025年）中的量子计算工具箱实现的，并在理想的（无噪声）量子模拟器上执行，为硬件开发提供了方法的基线评估。

本文的组织结构如下。第2.1节总结了用于受辐照晶体的经典晶体塑性（CCP）框架。第2.2节回顾了本工作中考虑的经典应变梯度晶体塑性（SGCP）公式，即CMSlip、CMCurl和CMCurl-Star模型。第2.3节首先回顾了EMCurl-SGCP模型，然后介绍了EMCurl-Star和EMSlip扩展。第2.5节详细介绍了所提出的量子傅里叶变换方法，包括其集成到CCP框架和每个SGCP模型中。第3节报告了使用所提出的QFT方法获得的高分辨率二维和三维多晶聚集体的数值结果，重点讨论了应变局部化（例如剪切带或清晰通道）的形成以及平均晶粒大小的影响。还提供了一整套附录以支持可重复性并记录实现细节。附录A总结了量子计算中的基本概念，包括量子比特、量子门和电路表示。附录B介绍了用于将晶体状态变量转移到量子状态的状态准备程序（幅度编码）。附录C详细介绍了QFT电路及其在量子有限差分方法中的应用。附录D提供了针对代表性本构响应（包括线性弹性、完美塑性、线性软化和非线性软化行为）的QFT方法实现的验证。最后，附录E提供了在拉伸载荷下配备SGCP模型的二维单晶模拟中滑移和屈曲带形成的额外分析。

本研究采用了晶体塑性中张量分析的典型符号，例如标量、一阶（向量）、二阶、三阶和四阶张量分别表示为A、A?、A～、A?和A≈。Levi-Civita三阶排列张量表示为??。此外，字母上的上标如tot、e、p、.、?、?和_分别表示张量的总部分、弹性部分、塑性部分、时间导数部分、傅里叶变换部分、波动部分和空间平均值部分。数学符号如?、×、?、.、:、?、?.、?×和Δ表示张量积、叉积、卷积、点积、双重收缩、梯度、散度、旋度和拉普拉斯算子。例如，二阶张量的旋度在实空间笛卡尔坐标系中计算为?×X～=?jmsXim,se?i?e?j，使用爱因斯坦求和约定。此外，二阶张量的范数计算为∥X～∥=X～:X～，二阶张量的对称部分和斜对称部分分别计算为Sym(X～)=(X～+X～T)/2和Skew(X～)=(X～?X～T)/2。最后，在量子计算中使用了狄拉克（bra-ket）符号（见附录A）。

本节介绍了本研究中使用的理论框架和计算方法。首先，在无限小变形框架中回顾了经典晶体塑性（CCP）公式。由于辐照效应对晶格的影响，如屈服应力的增加、延展性和加工硬化区域的减少以及应变局部化，这些结果报告了通过Voronoi镶嵌生成的周期性二维和三维多晶聚集体的数值结果，这些聚集体由随机取向的晶粒组成（Quey等人，2011年）。微观结构在高分辨率下离散化，能够详细评估局部变形和微观结构场。所有模拟都在相对较低的施加应变率10^-6s^-1下单轴拉伸载荷下进行，以促进准静态响应。

在本研究中，引入了增强型应变梯度晶体塑性（ESGCP）理论的一个新分支，即EMSlip模型，该模型结合了与每个滑移系上的累积剪切应变梯度相关的二次能量贡献。与基于Nye张量的ESGCP公式类似，EMSlip模型引入了一个物理上合理的长度尺度，该尺度定义为受辐照晶体的演变微观结构的函数。因此，CMRediT作者贡献声明Amirhossein Lame Jouybari：写作-审稿与编辑、写作-原始草稿、可视化、验证、软件、资源、项目管理、方法论、调查、形式分析、数据管理、概念化。Leon Cizelj：监督、项目管理、资金获取。

作者衷心感谢斯洛文尼亚研究与创新机构提供的财政支持（资助号P2-0026以及年轻研究员Amirhossein Lame Jouybari）。作者还感谢Samir El Shawish的建设性讨论。

作者声明他们没有已知的竞争性财务利益或个人关系可能会影响本文报告的工作。