张量网络的自动微分（文献翻译）

痛定思痛。 2023-01-02 15:27 151阅读 0赞

本文为文献 Differentiable Programming Tensor Networks 的原创译文，转载请注明出处。  
Author：小岚岚

### 目录 ###

*   *  I.引言
     *  II．一般理论
     *   *  A．自动微分
         *  B．张量网络收缩的计算图
         *   *  1.张量重整化群
             *  2.角转移矩阵重整化群
     *  III．技术成分
     *   *  A . 通过代数操作稳定后向传播
         *   *  1.对称eigensolver
             *  2. 奇异值分解
             *  3.QR因子分解
         *  B．运用checkpointing函数节省内存的反向自动微分
         *  C． 不动点迭代的反向传播
         *  D． 高阶微分
     *  IV．应用
     *   *  A . 自由能高阶微分
         *  B． iPEPS基于梯度的优化
     *  V . 讨论

## I.引言 ##

张量网络是研究经典统计物理和量子多体物理问题的重要方法\[1-3\]。近年来，它的应用迅速扩展到不同的领域，包括量子电路的仿真与设计\[4-7\]、量子误差校正\[8,9\]、机器学习\[10-14\]、语言建模\[15,16\]、量子场论\[17-20\]、nolographic二象性\[21,22\]。

张量网络在一般情况下的优化是众多研究方向的核心问题之一。尽管一维矩阵积状态的优化方案非常成功\[23-28\]，但二维或高维张量网络的优化一直是一个具有挑战性的课题。这个困难部分原因是张量收缩的计算成本高，还要部分原因是在高维情况下缺乏有效的优化方案。对于无限平动不变量子系统，优化张量网络状态的挑战尤其紧迫。在这些情况下，同一个张量以多种方式影响变分能量，从而导致高度非线性的优化问题。基于近似虚时算符的优化方案\[29-33\]努力解决目标函数中的非局部依赖性。文献\[34,35\]采用基于梯度的优化，结果明显改善。  
然而，即使对于一个简单的物理哈密顿量，从解析的角度推导张量网络状态梯度也很麻烦，而且容易出错。这一事实限制了基于梯度的张量网络状态优化在更复杂系统中的广泛应用。替代方法如使用数值导数计算梯度，精度和效率有限;因此，它们只适用于变化参数较少的情况\[36,37\]。使用链式法则手动推导梯度只适用于精心设计的简单张量网络结构\[38\]。

可微编程通过以完全可微的方式组合整个张量网络程序并使用自动微分来优化它们，为这些问题提供了优雅而有效的解决方案。在本文中，我们提出基本的自动微分技术，计算张量网络程序的(高阶)导数有效地数值精度。这一进展为一般情况下基于梯度(甚至基于hessian)的张量网络状态优化打开了大门。此外，计算张量网络算法输出的(高阶)梯度提供了一个简单的方法来计算物理量，如比热和磁化率。可微规划方法不知道详细的格点几何、哈密顿量和张量网络收缩方案。因此，该方法的通用性足以支持广泛的张量网络应用。我们的重点是涉及二维无限张量网络的应用，与传统方法相比，可微规划技术提供了显著的好处。我们证明，在通过奇异值分解(SVD)解决数值稳定微分和不动点迭代的内存高效实现等重大技术挑战后，可以获得张量网络状态变分优化的最先进结果。  
本文的组织结构如下。在第二节中，我们在张量网络算法的背景下引入自动微分，并将张量网络收缩表示为计算图。在第三节中，我们给出了张量网络算法的稳定和可扩展可微规划的关键技术。在第四节中，我们证明了该方法的能力，并将其应用于无限正方晶格上的经典Ising模型和量子海森堡模型。最后，我们在第五节对未来的研究方向进行了展望。我们的代码实现是公开的(参见参考\[39\])。

## II．一般理论 ##

### A．自动微分 ###

通过计算图自动微分是一个统一的框架，涵盖了机器学习训练神经网络，优化量子物理张量网络等等。我们首先回顾自动微分的核心思想，然后解释它在以计算图形式表述的各种张量网络收缩算法中的应用。  
自动微分机械地计算用计算机程序\[40\]表示的计算过程的导数。与数值微分不同，自动微分计算的导数值精确到机器精度。其性能有、普适的理论保证，不超过原程序的算法复杂度\[41,42\]。  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2xpdHRsZWxhbmxhbg_size_16_color_FFFFFF_t_70]

*图1所示。计算图的反向自动微分。黑色箭头表示从输入到输出的正向评价函数。红色箭头表示伴随反向传播的后退步骤。(a)链图。(b)一个更一般的计算图。在向后传递过程中，给定节点的伴随点根据公式(2)计算。*

自动微分是现代深度学习应用的计算引擎\[43,44\]。在量子最优控制\[45\]和量子化学计算中也有应用，如计算力\[46\]，优化基参数\[47\]等。自动微分的中心是计算图的概念。计算图是由基本计算步骤组成的有向无环图。图的节点表示数据，数据可以是标量、向量，矩阵，或者张量\[48\]。图连通性表示计算过程中数据流的依赖性。最简单的计算图是图1(a)所示的链。例如，从输入参数向量theta开始，可以计算一系列中间结果，直到最终的输出L，我们假设它是一个标量。所谓的正向求值就是按顺序遍历链图，即 θ → T 1 → . . . → T n → L \\theta \\rightarrow T^1\\rightarrow...\\rightarrow T^n \\rightarrow \\mathcal\{L\} θ→T1→...→Tn→L。  
为了计算目标函数相对于输入参数的梯度，我们可以利用链式法则  
 ∂ L ∂ θ = ∂ L ∂ T n ∂ T n ∂ T n − 1 . . . ∂ T 2 ∂ T 1 ∂ T 1 ∂ θ . ( 1 ) \\frac\{\\partial \\mathcal\{L\} \}\{\\partial \\theta\}= \\frac\{\\partial \\mathcal\{L\} \}\{\\partial T^n\} \\frac\{\\partial T^n \}\{\\partial T^\{n-1\}\}... \\frac\{\\partial T^2 \}\{\\partial T^1\} \\frac\{\\partial T^1\}\{\\partial \\theta\}. (1) ∂θ∂L=∂Tn∂L∂Tn−1∂Tn...∂T1∂T2∂θ∂T1.(1)  
由于我们考虑的情况是，输入维数大于输出维数，在Eq.(1)中通过使用一系列雅可比向量乘积从左到右相乘来计算梯度是更有效的。根据图1(a)所示的计算图，将图向后遍历，并将梯度信号从输出端传播回输入端。用这种方法计算导数叫做反向自动微分。这种方法通常称为反向传播算法\[43\]，被认为是训练深度神经网络最成功的方法。引入伴随变量 T ˉ = ∂ L / ∂ T \\bar T=\\partial \\mathcal\{L\} /\\partial T Tˉ=∂L/∂T 来表示最终输出L相对于变量T的梯度。我们可以看到，反向自动微分将伴随变量从 T n ‾ = L ˉ ∂ L / ∂ T n \\overline \{T^n\} =\\bar \\mathcal\{L\} \\partial \\mathcal\{L\}/\\partial T^n Tn=Lˉ∂L/∂Tn， L ˉ = 1 \\bar \\mathcal\{L\}=1 Lˉ=1  
一直传播到 T i ‾ = T i + 1 ‾ ∂ T i + 1 / ∂ T i , i = n − 1 , . . . 1 \\overline\{T^i\}=\\overline\{T^\{i+1\}\}\\partial T^\{i+1\}/\\partial\{T^i\},i=n-1,...1 Ti=Ti\+1∂Ti\+1/∂Ti,i=n−1,...1  
最后计算 θ ‾ = T 1 ‾ ∂ T 1 / ∂ θ \\overline \\theta=\\overline \{T^1\}\\partial\{T^1\}/\\partial \\theta θ=T1∂T1/∂θ 。在每一步中，通过局部向量-雅可比矩阵乘积向后传播伴随矩阵。伴随反向传播图可以很好地推广到更复杂的计算图。例如图1(b)中的数据节点Tl通过两条不同的下游计算路径影响最终输出。在向后传递过程中，需要为伴随节点把子节点的所有贡献求和。总的来说反向传播满足  
 T i ‾ = ∑ j : c h i l d   o f   i T j ‾ ∂ T j ∂ T i ( 2 ) \\overline\{T^i\}=\\sum\_\{j:child\\space of\\space i\}\\overline\{T^j\}\\frac\{\\partial T^j\}\{\\partial T^i\}(2) Ti=j:child of i∑Tj∂Ti∂Tj(2)  
反向自动微分算法可以理解为计算图上的一个消息传递过程。通过前向传播定义的计算图的拓扑排序后（？），从伴随 L ‾ = 1 \\overline \\mathcal\{L\}=1 L=1的输出节点向后来看图。每个节点从它的子节点收集信息来计算它自己的伴随节点，然后将这个信息传递给它的父节点。因此，我们可以在一个前向和一个后向的过程中计算关于所有参数的梯度。通常，在前向传递中缓存必要的信息，以便在后向传递中有效地计算向量-雅可比矩阵乘积。  
可微程序的组成部分称为基元。基元可以是基本操作，如加法、乘法和数学函数\[50\]。除了普通的前向函数外，每个基元都有一个相关的反向函数。反向函数根据向量雅可比矩阵乘积Eq.(2)来反向传播伴随矩阵。注意，不需要显式实例化或存储完整的雅可比矩阵。此外，人们可以把许多基本计算步骤组合在一起作为一个基元。例如，像矩阵乘法这样的线性代数运算可以被视为一个基元。通过这种方式，这些自定义基元的前向传递可以作为黑盒操作。以这样一种模块化的方式设计微分程序，可以控制自动微分的间隔级别。为特定领域的问题定制基元有几个优点：首先，这可以减少计算和内存成本。其次，在某些情况下，将几个步骤组合在一起在数值上更稳定。第三，可以将对外部库的函数调用包装成基元，而不需要跟踪每个单独的计算步骤。  
现代机器学习框架通过各种机制支持自动微分。例如，TensorFlow\[51\]使用定制语言显式地构造计算图;autograd\[52\]、PyTorch\[53\]、Jax\[54\]在运行时跟踪程序的执行顺序;Zygote\[55\]在编译器级执行源代码转换。这些框架允许通过函数调用、控制流、循环和许多其他计算指令来进行微分。在这些基础结构的基础上，可微编程正在成为一种新的编程范式，它强调组装可微组件，并通过端到端优化\[44\]学习整个程序。让机器机械地处理自动判别问题，大大减少了人类的费力工作，并将人类的注意力重新分配到设计更复杂、更有创造性的深度学习算法上。  
最后，我们注意到Eq.(1)也可以从右到左求值，这对应于前向模式自动微分。前向模自动微分的运算类似于微扰理论。我们可以在一个单一的向前传递过程计算目标函数和梯度，而不存储任何中间结果。然而，前向自动判别对于输入维数远大于输出维数\[44\]的计算图是不利的。因此，深度学习的大部分工作都采用了反向自动微分。出于同样的原因，我们将重点讨论张量网络程序的反向自动微分。

### B．张量网络收缩的计算图 ###

张量网络程序将输入张量映射到输出，我们假设输出是一个标量。根据本文，输入量可以表示经典的配分函数或量子波函数，输出可以是各种感兴趣的物理量。在概念上，将自动微分应用于张量网络程序是很简单的，方法是将计算，特别是将张量网络收缩表示为计算图。  
作为一个教学的例子，考虑无限一维Ising模型的配分函数,  
 Z = l i m N → ∞ T r ( T N ) , T = ( e β e − β e − β e β ) Z = lim\_\{N\\rightarrow \\infty\}Tr(T^N),T=\\begin\{pmatrix\} e^\{\\beta\} & e^\{-\\beta\} \\\\ e^\{-\\beta\} & e^\{\\beta\} \\end\{pmatrix\} Z=limN→∞Tr(TN),T=(eβe−βe−βeβ)  
其中T 是二阶的表示玻尔兹曼分布的张量(矩阵)。通过反复平方矩阵T并求迹，可以在热力学极限下数值得到配分函数。计算图为图1(a)所示的简单链结构。对于这种计算过程的微分涉及到通过矩阵跟踪和乘法运算反向传播伴随矩阵，这很简单。  
在这一点上，也值得区分精确和近似的张量网络收缩方案。在树形计算图中张量网络可以被精确高效地收缩。  
然而，其他精确的方法，如用于计数和模拟量子计算的方法\[56,57\]，通常会随着问题的尺寸扩大呈指数级增长。然而，对这些精确的算法应用自动微分是很简单的，因为它们主要涉及张量的重整和收缩。  
我们将重点放在通过近似张量收缩求导的不那么trivial的情况上，这通常涉及截断张量因式分解或变分逼近。它们涵盖了重要的张量网络应用，与其他数值方法相比具有很大的优势\[1-3\]。特别地，我们对收缩无限张量网络感兴趣，其中的基本数据结构是体张量。压缩方案大致分为三类:基于粗粒化变换的压缩方案\[58-64\]，基于角转移矩阵的压缩方案\[65-67\]，以及基于矩阵乘积状态的压缩方案\[2,24,68\]。由于最后两种压缩方案与\[2\]密切相关，下面我们将重点研究张量重正化群(第II b1)和角转移矩阵重正化群方法(第II b2)的自动微分。

#### 1.张量重整化群 ####

张量重整化群(TRG)方法通过对体张量的分解和迭代截断\[58\]来收缩张量网络。图2(a)显示了作为计算图的训练迭代的一个步骤，其中包括(1)使用SVD以两种方式分割体张量，其中-我们将奇异值和向量截断为规定的维度x，(2)将在最后一步中生成的四个3leg张量组装成一个4 leg张量。在这个收缩之后，我们得到了下一个迭代的一个新的张量。  
![在这里插入图片描述][20210106145844778.png]

TRG方法使晶格尺寸以指数速度增长。所以，在几十次迭代后很快就会达到热力学极限。注意，为了数值稳定性，我们需要在每次迭代后对张量元素重新进行尺度变换。TRG方法的计算代价尺度为x的六次方，存储代价尺度为x的四次方。展开迭代之后，TGR方法的计算图类似于图1(a)所示的简单链图。在每个迭代步骤中，基本操作是张量指标置换、截断SVD和张量收缩。由于这些操作都是可微的，我们可以通过TRG程序反向传播来计算下游目标函数对输入张量的导数。

#### 2.角转移矩阵重整化群 ####

角转移矩阵重正化群(CTMRG)\[65\]的计算图拓扑结构更有趣。CTMRG计算的目标是获得代表体张量的环境自由度(d.o.f.)的收敛角张量和边张量。  
在大部分情况下张量具有正方晶格的全部对称性，CTMRG一个迭代的步骤图2所示(b):(1)合同大部分张量与角落和边缘张量形成4-leg张量,(2)执行截断奇异值分解4-leg张量,保持单一维度截止x(保持截断奇异矩阵等角投影);(3)从第一步开始，将等距法应用于4-leg张量，求出新的角张量;(4)为下一步应用相同的等距来找到一个新的边张量。迭代这个过程直到收敛。我们可以看到，在CTMRG迭代的每一步中，维度为d的体张量都是相同的。由于这个原因，收敛环境张量将以一种复杂的方式依赖于体张量。  
与\[58\]训练方法不同，CTMRG方法使系统的大小线性增长。所以，你可能需要迭代更多的步骤才能在CTMRG中达到一致。另一方面，CTMRG的计算复杂度O(dx3)和存储代价O(dx2)在截止键维数上都小于TRG。  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2xpdHRsZWxhbmxhbg_size_16_color_FFFFFF_t_70 1]

*图2所示。(a)TRG训练的迭代步骤。(b) CTMRG的迭代步骤。每个张量是计算图中的一个节点。计算图中的基本函数是SVD和张量收缩。*

## III．技术成分 ##

利用反向自动微分来计算张量网络程序的梯度，需要追踪原函数的组成，并将伴随信息在计算图上传回。现代可微编程框架\[51-55\]已经自动地处理了它们的基本数据结构——可微张量的跟踪和反向传播。  
我们需要关注的是识别出合适的张量网络程序的基元，并为反向传播定义它们的向量-雅可比矩阵乘积。张量网络算法的关键成分是矩阵和张量代数。通过这些操作，已有关于反向传播的既定结果\[69-71\]。首先，用定制的后向函数将所有BLAS例程包装为基元是很简单的。接下来，尽管不那么简单，但也可以为许多LAPACK例程(如特征解算器、SVD和QR分解)推导后向规则\[69\]。把这些线性代数运算当作基本运算，就可以组成一个可微程序，并有效地实现矩阵库。然而，可微张量网络程序的稳定和可扩展的实现存在一些实际障碍。首先，特征解和奇异值分解的后向传递可能面临特征值或奇异值简并的数值不稳定。其次，反向自动微分可能会导致大量内存消耗，从而使我们无法达到与普通张量网络程序相同的维度。我们将在下面给出这些问题的解决方案。

### A . 通过代数操作稳定后向传播 ###

关于矩阵导数我们提出几个关键的结果，涉及线性代数运算与相关的张量网络算法。由于反向自动微分的模块化特性，只需指定本地的后向函数就可以将这些组件集成到可微分程序中。我们评述了它们与物理文献的联系，并特别关注稳定的数值实现\[39\]。更多信息，可以参考参考资料。

#### 1.对称eigensolver ####

正向传递读取 A = U Λ U T A=U\\Lambda U^T A=UΛUT，其中  Λ \\Lambda Λ 是“特征值的对角矩阵;正交矩阵U的每一列都是一个对应的特征向量。在计算图中，节点A有两个子节点U和 λ。  
在反向传递中，已知伴随U和A，我们有\[69,71\]  
 A ‾ = U \[ Λ ‾ + F ⊙ ( U T U ‾ − U ‾ T U ) / 2 \] U T , \\overline A=U\[\\overline \\Lambda+F\\odot(U^T\\overline U-\\overline U^T U)/2\]U^T, A=U\[Λ\+F⊙(UTU−UTU)/2\]UT,(3)  
where F i j = ( λ j − λ i ) − 1 F\_\{ij\}=(\\lambda\_j-\\lambda\_i)^\{-1\} Fij=(λj−λi)−1if  i ≠ j i\\ne j i=jand zero otherwise.

符号 ⊙ \\odot ⊙表示元素的阿达玛乘积。我们可以很容易地证明梯度也是一个对称三矩阵。  
式(3)可视为“反向”扰动理论。当下游计算不依赖于本征态时，即 U ‾ = 0 \\overline U=0 U=0，则反向方程与著名的赫尔曼-费曼定理\[721\]有关，该定理将微扰与哈密顿量及其本征值联系起来。参考文献\[73\]将公式(3)应用于基于能谱的反向哈密顿设计的特殊情况。  
特征值差出现在F的分母中是一阶非简并微扰理论的提示。参考文献\[47\]通过特征求解器考虑了反向传播的稳定性，因此转向效率较低的前向自动微分来实现Hartree-Fock基的变分优化。在许多物理问题中，最终的目标函数依赖于本征值和本征态，与规范无关，（？）例如，占据态的二次形式。在这种情况下，F的分母中只会出现占据态与未占据态之间的特征值差异，这是线性响应理论中常见的模式\[74\]。因此，简并特征值并不一定会对这些物理应用造成问题\[75\]。在实践中，我们发现，采用洛伦兹展宽 1 / x → x / ( x 2 + ϵ ) , ϵ = 1 0 − 12 1/x \\rightarrow x/(x^2+\\epsilon),\\epsilon =10^\{-12\} 1/x→x/(x2\+ϵ),ϵ=10−12 可以使计算稳定，但代价是引入一个小误差的梯度;参见参考文献\[71\]。

#### 2. 奇异值分解 ####

张量网络算法中一个普遍存在的运算是矩阵SVD，它用于张量网络的规范化和因式分解\[1,2,25\]。前向传递形式为  A = U S V T A=USV^T A=USVT，其中A的尺寸为(m, n)，  
 U , V T U,V^T U,VT的尺寸分别为(m, k)和(k, n)， k= min(m, n)。注意，S是一个包含奇异值s的对角矩阵。在反向自动微分中，给定伴随矩阵 U ‾ , S ‾ , V ‾ \\overline U,\\overline S,\\overline V U,S,V，可以得到\[70\]  
![在这里插入图片描述][20210106151530861.png]

否则为0。  
为了防止简并奇异值情况下的数值问题，我们对第一项使用与第III节A 1相同的洛伦兹展宽，这在我们的经验中很有效。在实际应用中，从随机张量开始的变分张量网络计算，得到精确简并特征值的机会很小。即使出现这种情况，四舍五入也是一个合理的解决方案。然而，对于由于内在原因引起的简并情况\[60,62,76,77\]，只要端到端梯度有明确的定义，仍然可以得到正确的梯度。最后，对式(4)中消失的奇异值求逆并不需要考虑，因为对应的空间通常会被截断。

#### 3.QR因子分解 ####

QR因子分解常被用于张量网络的正则化。在前向传递中，将A=QR因式分解，其中  Q T Q = I Q^TQ=I QTQ=I，R为上三角矩阵\[78\]。根据矩阵A的维数(m, n)，反向函数有两种情况。  
对于矩阵A的输入形状 m ≥ n m\\ge n m≥n ，R是一个nxn矩阵。向后传递为\[71\]  
 A ‾ = Q ‾ + Q c o p y l t u ( M ) \] R − T ,   ( 5 ) \\overline A=\\overline Q+Qcopyltu(M)\]R^\{-T\},\\space (5) A=Q\+Qcopyltu(M)\]R−T, (5)  
其中  M = R R ‾ T − Q ‾ T Q M=R\\overline R^T-\\overline Q^T Q M=RRT−QTQ，copyltu函数通过将输入矩阵的下三角形复制到上三角形来生成对称矩阵： \[ c o p y l t u ( M ) \] i j = M m a x ( i , j ) , m i n ( i , j ) \[copyltu(M)\]\_\{ij\}=M\_\{max(i,j),min(i,j)\} \[copyltu(M)\]ij=Mmax(i,j),min(i,j) 。  
我们可以通过求解一个三角系数矩阵的线性方程组来处理 R − T R^\{-T\} R−T的乘法。  
对于m <n的情况，Q的尺寸为(m, m)， R是一个mxn矩阵。我们记A=(X, Y)和R= (U, V)，其中X和U是大小为(m, m)的满秩方阵。分解可以分为两个步骤:首先，X=QU唯一地决定Q和U，然后计算 V = Q T Y V=Q^TY V=QTY 。应用链式法则，可以得到逆向法则  
![在这里插入图片描述][20210106152023203.png]

### B．运用checkpointing函数节省内存的反向自动微分 ###

张量网络的反向自动微分的直接实现有很大的内存开销，因为需要在前向传递中存储中间结果，以便在后向传递中计算向量雅可比矩阵乘积。存储变量的数量与实现自动区分时的间隔级别有关。在任何情况下，反向模式自动区分的内存消耗将与计算图的深度成比例。这种开销对于具有大bond维度和大量重正化群(RG)迭代的张量网络尤其令人担忧。解决反向自动微分内存问题的方法是一种著名的技术，称为checkpointing\[50,79\]。其思想是用计算时间交换节省的内存。以式(1)中的链计算图为例，在正向过程中每隔几步就可以存储张量。在向后传递过程中，只要有需要，就通过向前运行一小段计算图来重新计算中间张量。通过这种方式，可以极大地减少内存使用，而不需要超过两倍的计算工作量。从另一个角度来看，checkpointing相当于定义定制的基元，它封装了计算图的很大一部分。这些基元有它们自己的特殊向后规则，在本地再次运行向前传递，然后反向传播伴随函数。因此，在正向传递中，不需要缓存这些checkpointing基元的内部状态。  
checkpointing是应用于任意拓扑结构的计算图。当应用于张量网络算法时，我们可以自然地考虑图2中的重整化步骤作为checkpointing基元。这样，就可以避免在正向传递中存储一些较大的中间张量。

### C． 不动点迭代的反向传播 ###

不动点迭代是张量网络算法中的一种循环迭代模式。例如，迭代函数 T i + 1 = f ( T i , θ ) T^\{i+1\}=f(T^i,\\theta) Ti\+1=f(Ti,θ) 直到达到一个收敛张量T\*，并将其用于下游计算。要计算关于参数theta的梯度，当然可以将迭代展开到一个深度计算图，并直接应用反向自动微分。然而，这种方法的缺点是，如果需要很长时间的迭代才能找到不动点，就会消耗大量内存。我们可以利用不动点方程上的隐函数定理来解决这个问题\[80\]。对T\*两边求导，得到  
![在这里插入图片描述][20210106152114342.png]

第二行将方括号中的矩阵逆展开为几何级数。因此，要通过一个不动点迭代进行反向传播，基本操作仅涉及向量-雅可比矩阵与单步迭代函数的乘积。在后向函数中，执行迭代来累积伴随函数直到达到收敛。几何级数在不动点的收敛速度与前向迭代相同\[80\]许多张量网络收缩方案，包括在第二节B 2中回顾的CTMRG方法。属于定点迭代的框架。因此，可以使用Eq.(7)进行CTMRG计算，迭代在图2(b)所示的RG步骤上。我们注意到在参考文献中导出的无限张量网络收缩的解析梯度。\[34,35\]包含了类似的模式，它是几何级数的求和。  
与第III B节的checkpointing技术类似，Eq.(7)也减少了反向自动微分的内存使用，因为在前向迭代中不需要存储长链的中间结果。此外，由于下游目标函数与不动点张量是如何获得的无关，人们也可以利用前向传播\[82\]中的加速迭代方式\[81\]。然而，当将式(7)应用于张量网络RG的微分算法时，有一个警告。人们可能需要特别注意RG迭代中的冗余全局规范，以确保不动点方程确实成立。

### D． 高阶微分 ###

由于梯度也可以用计算图来表示，所以可以通过对图进行自动微分来计算二阶导数。这样，在原则上，可以计算任意高阶导数的程序使用自动微分\[50\]。深度学习框架\[51-55\]对计算高阶导数提供开箱即用的支持计算高阶导数的能力支持基于hessian的张量网络状态优化，如牛顿法\[83\]。然而，显式地计算和反演完整的Hessian矩阵可能是非常昂贵和不必要的。我们可以通过  
![在这里插入图片描述][20210106152130815.png]  
而无需显式构造Hessian矩阵\[84\]。这一过程对于用于牛顿法的迭代线性方程求解器是充分的。

## IV．应用 ##

我们提出两个应用来证明微分编程张量网络方法在统计物理和量子多体问题中的通用性。公开代码实现\[39\]采用PyTorch\[85\]和一个定制的线性代数自动分化图书馆提高数值稳定性(见讨论第三秒。)。然而,我们注意到可以通过深度学习框架如autograd \[86\], TensorFlow87号，Jax\[88\]，Zygote\[89\]轻松得出结果。

### A . 自由能高阶微分 ###

考虑温度倒数为beta的正方格点Ising模型，配分函数可以表示为 bond维度为2的二位张量网络  
![在这里插入图片描述][20210106152157932.png]

体张量为  
![在这里插入图片描述][20210106152207231.png]

其中  λ u = e β + ( − 1 ) u e − β \\lambda\_u=e^\{\\beta\}+(-1)^ue^\{-\\beta\} λu=eβ\+(−1)ue−β。（？）我们使用第IIB 1节中讨论的TRG方法收缩无限张量网络。我们使用截断bond维度x=30，并迭代30个TRG步骤。最后，我们通过追踪体积张量得到配分函数Eq.(8)和自由能。  
接下来，我们通过自动微分直接对自由能求导来计算物理可观测值，如能量密度和比热，如图3所示。注意，能量密度显示了一个扭结。比热在临界温度  β c = l n ( 1 + 2 ) / 2 ≈ 0.44068679 \\beta\_c=ln(1+\\sqrt 2)/2 \\approx0.44068679 βc=ln(1\+2)/2≈0.44068679附近出现峰值。与数值微分不同，这些结果不受有限差分误差\[61,91\]的影响。精确计算张量网络算法的高阶导数将有助于研究热相变和量子相变。我们注意到，原则上可以通过直接计算能量方差来获得比热\[35,92\]，但这涉及到用张量网络表示的几何级数的繁琐求和。  
有其他的方法来用自动微分计算比热。例如，可以直接利用杂质张量计算能量，然后求一阶导数得到比热。另外，也可以使用前向模式自动微分，因为只有一个输入参数beta需要微分。我们特意选择了目前的方法，通过反向自动微分技术来突出可微分编程的强大功能。通过反向传播整个TRG过程，特别是SVD，可以使用高阶导数计算物理可观测值。值得注意的是，考虑到Ising模型\[47\]的Z2对称性导致的许多简并奇异值，这种方法是有效的。为了获得正确的物理结果，使用第IIIA 2节提到的数值稳定方法的SVD后向函数是很关键的。

### B． iPEPS基于梯度的优化 ###

我们考虑正方晶格的反铁磁海森堡模型的变分研究，哈密顿量为  
![在这里插入图片描述][20210106152343565.png]

我们考虑一个iPEPS作为变分波函数，变分参数则是iPEPS中的元素  
![在这里插入图片描述][2021010615234816.png]

S 代表物理指标，其他的 u、l、d、r 则是bond维度D的虚拟 dof（环境自由度）。我们用随即高斯分布初始化元素。iPEPS的重叠部分组成一个张量网络，其中体张量是bond维度为D2的双层张量  
![在这里插入图片描述][2021010615235389.png]

为了收缩由这个体张量形成的无限张量网络，我们使用在第二节B 2中复习过的CTMRG方法。我们通过部分追踪体张量的leg来初始化角张量和边张量；然后进行CTMRG迭代，直到我们在角落和边缘张量上达到收敛。收缩后，我们可以评估预期能量  ⟨ Ψ ∣ H ∣ Ψ ⟩ / ⟨ Ψ ∣ Ψ ⟩ \\langle \\Psi\\vert H\\vert\\Psi\\rangle/\\langle\\Psi\\vert\\Psi\\rangle ⟨Ψ∣H∣Ψ⟩/⟨Ψ∣Ψ⟩。由于问题的平动不变性，考虑键上的期望能量就足够了。  
![在这里插入图片描述][20210106152402929.png]

式(13)中的黑矩形为作用于键上的哈密顿算符。我们在哈密顿函数中做了基旋转这样基态就会有一个单位置的单元（？）。我们用截止bond维度X= 30,50,80,100,144,160来分别表示D = 2.3. …7。由于期望能量随截止维数的减小而减小\[36,94\]，近似的CTMRG收缩给出了基态能量的一个变分上限的期望能量。(13)通过角张量和边张量对式(11)中的变分参数具有显式和隐式的依赖性。  
我们用自动微分法计算式(13)对单层张量式(11)的梯度，自动解决计算图中复杂的内部依赖关系。梯度计算所花费的时间与预期能量的前向估计相当。然后，我们使用带有自动计算梯度的准牛顿L-BFGS算法\[83\]优化iPEPS。经过几百次的函数和梯度计算后，可以迅速达到最优。图4(a)显示了不同bond维度下与外推量子蒙特卡罗(QMC)结果\[95\]相比的能量相对误差。基态能量的精度可与最先进的结果相媲美\[34,35\]，它被证明比基于假想时间投影的简单和完全更新算法更精确\[29-33\]\[96\]。请注意，Ref.\[35\]中的ansatz和我们的ansatz只包含Ref.\[34\]中的一半的变分参数，所以在D = 2,3时，能量结果略高于Ref.\[34\]。然而，对于更大的键维数D= 4,5,6,7，我们的计算达到了无限方格海森堡模型的最低变能。图4(b)为优化状态下测得的交错磁化，在更大的bond维度达到了外推QMC的结果。  
![图4所示。(a)二维S=1/2反铁磁海森堡模型与以前变分结果的相对误差\[33,35\]。精确度是相对于“外推的QMC结果\[95\]”进行测量的。(b)交错磁化的比较，虚线为外推的QMC结果\[95\]。简单而完整的更新引用数据也来自Ref.\[34\]。][4_a_S_1_2_33_35_QMC_95_b_QMC_95_Ref._34]

*图4所示。(a)二维S=1/2反铁磁海森堡模型与以前变分结果的相对误差\[33,35\]。精确度是相对于“外推的QMC结果\[95\]”进行测量的。(b)交错磁化的比较，虚线为外推的QMC结果\[95\]。简单而完整的更新引用数据也来自Ref.\[34\]。*

为了获得bond维度D> 4的结果，我们需要使用第三节B中的checkpointing技术或第三节C中的不动点迭代，以保持足够低的内存预算，以容纳一个12G内存的Nvidia P100 GPU卡。令人鼓舞的是，只要付出适度的努力，就可以在不断优化iPEPS时达到最先进的性能\[34,35\]。这一成功也证明了我们确实可以稳定出现在科学计算\[47\]中的线性代数运算的反向自动微分。  
我们注意到，目前的方法也适用于有限系统或问题，具有更大的单位细胞，更复杂的哈密顿量\[97-100\]，和更复杂的压缩方案提高效率\[94\]，这有望为量子多体问题提供新的物理结果。

## V . 讨论 ##

通过自动微分来计算梯度极大地增强了现有张量网络算法的能力。研究者可以专注于核心张量网络收缩算法，而不必担心繁琐的梯度计算。自动微分的计算复杂度与张量网络的正向收缩相同。  
除了大大减少人力投入外，自动微分方法计算的梯度也与参考值略有不同。该方法通过自动微分来计算近似能量密度的精确数值梯度。\[34 35\]的方法则首先推导能量梯度的解析表达式为无穷张量网络，然后将这些网络近似收缩得到一个近似梯度。因此，这两种方法分别对近似和近似导数\[44\]求导。除了参考文献\[44\]的一般建议外，我们发现通过近似张量网络收缩的微分对于那些解析导数推导和近似都很复杂的无穷系统是有利的。  
在本文中，我们主要关注自动微分的高级应用，它通过整个压缩算法来优化张量网络和计算物理可观测值。同样的技术也适用于低级情况，如寻找最优截断基或张量网络的变分变换\[64\]。此外，该方法除了优化量子问题的期望能量外，还与张量网络的变分收缩有关\[2101\]。我们期望可微编程技术将成为标准张量网络工具箱中不可或缺的一部分。  
用深度学习框架实现张量网络程序的一个优势为可以轻松使用GPU加速。这一工作的计算是在单个GPU卡上完成的。进一步推进这一研究方向，我们认为在新兴的专业硬件上部署张量网络算法将是有益的。  
最后，对张量网络和神经网络在自动微分方面的区别进行了评述。典型的神经网络体系结构不涉及复杂的线性代数运算。然而，随着张量神经网络的发展\[102\]和各种张量网络在机器学习问题上的应用\[10- 14103\]，这两类网络之间的边界变得模糊。因此，当一个人转向更复杂的收缩方案时，本文提出的结果也将与张量网络机器学习应用相关。

[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2xpdHRsZWxhbmxhbg_size_16_color_FFFFFF_t_70]: /images/20221119/9aa31aebb9b345fa8bf0e2b68a63be73.png
[20210106145844778.png]: /images/20221119/b90617f347314a39ac35338bd33c06f6.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L2xpdHRsZWxhbmxhbg_size_16_color_FFFFFF_t_70 1]: /images/20221119/c8812dc4a3e24c8e8f88cd19113784a9.png
[20210106151530861.png]: /images/20221119/5f8b963b0ecc48da96e666aabed9ffc8.png
[20210106152023203.png]: /images/20221119/5cff44f9aece498282342e8c7e63ab67.png
[20210106152114342.png]: /images/20221119/c7bd32b2e65340998cedfddb08cf8d73.png
[20210106152130815.png]: /images/20221119/29c038c0a5044197a30ed2182462d926.png
[20210106152157932.png]: /images/20221119/326ecca2cc904050a57ca0b2e14e6a77.png
[20210106152207231.png]: /images/20221119/fc471f9b698a444290d823e99ecdbf5c.png
[20210106152343565.png]: /images/20221119/8aad736f99c84ba9b132f156304a7cad.png
[2021010615234816.png]: /images/20221119/11d9200f22b746328da8022f89c18e34.png
[2021010615235389.png]: /images/20221119/a3a8598594014b2abd1d92c8f3ea259b.png
[20210106152402929.png]: /images/20221119/4c54184314d54daa999f1d2d7d0c7c90.png
[4_a_S_1_2_33_35_QMC_95_b_QMC_95_Ref._34]: /images/20221119/0de793740be14a0bab506dc4d638e9a3.png