玻尔百科在科学计算领域，许多物理现象通过庞大的线性方程组进行建模，这些方程组通常由巨大的矩阵表示。这些矩阵的庞大规模带来了一个根本性挑战：显式地组装、存储和操作它们可能会耗尽最强大超级计算机的内存和计算资源。这一瓶颈限制了科学家和工程师能够解决问题的规模和复杂性。无矩阵方法为这一问题提供了一种范式转换的解决方案。它们完全绕开了构建矩阵的需要，转而专注于其根本目的：对向量的作用。本文将深入探讨这一优雅而强大的技术。第一章 原理与机制 将揭示无矩阵方法背后的核心思想，解释它们如何与迭代求解器协同工作，以及为何它们在现代硬件上异常高效。随后的 应用与跨学科联系 一章将展示其在从计算力学到量子化学等不同领域中的变革性影响。
想象一下你想描述一个人。你可以创建一个详尽无遗的清单，列出关于他的一切可以想到的事实——身高、体重、发色、每颗雀斑的确切位置。这将是一个庞大到令人不知所措的数据集合。或者，你可以描述他做什么：他是一位杰出的音乐家，一个富有同情心的朋友，一个会讲精彩故事的人。第二种描述，即基于行动的描述，往往比静态的属性列表更能有效地抓住一个人的本质。
在计算科学的世界里，我们处理源于物理现象建模的巨型矩阵时也面临类似的选择。矩阵，我们通常将其想象成一个巨大的数字网格，其本质上是对一个线性变换的描述。它真正的目的，它的本质，在于它对一个向量的作用——它如何拉伸、旋转和反射向量以形成一个新的向量。如果我们能直接捕捉这种作用，而无需写下那庞大的数字列表，会怎么样？这就是无矩阵方法背后核心的、极具解放性的思想。
让我们考虑一个线性方程组，这是科学计算的基础，写作 。我们寻找的是向量 ，当矩阵 作用于它时，会产生向量 。我们的传统直觉是把 看作一个有形的物体，一个存储在计算机内存中的矩阵。
无矩阵方法挑战了这一点。它认为：只要我们有一个函数，一个“黑箱”，能够为我们给定的任何输入向量 计算出乘积 ，我们就拥有了所需的一切。矩阵 本身可以仍然是“机器中的幽灵”——它的作用被明确定义，但其形式从未被显式组装。
例如，想象一个简单的物理系统，比如一串由弹簧连接的质量块。力和位移可能导出一个矩阵，其对向量 的作用由一组简单的规则定义：
为了计算向量 ，我们不需要看到矩阵。我们只需遵循这个“配方”。这个“配方”就是算子。这种从将矩阵视为静态对象到将其视为动态作用的视角转变，是第一个关键原理。
如果我们无法访问 的元素，我们怎么可能解出 中的 呢？像高斯消元法这类涉及系统性地修改矩阵元素的方法，就不在考虑之列了。然而，这扇门为一类庞大而强大的技术——迭代方法——敞开了。
可以把迭代方法想象成一个复杂的“猜温度”游戏。我们从一个解的初始猜测开始，称之为 。这个猜测几乎肯定是错的。为了找出它错得有多离谱，我们计算残差向量：。残差告诉我们我们得到的结果（）与我们想要的结果（）之间的差异。如果残差是一个零向量，我们就大功告成了！如果不是，迭代算法会利用 中的信息，智能地选择一个移动方向，从而产生一个更好的猜测 。这个过程不断重复，生成一系列猜测，并有望越来越接近真实解。
关键的洞见在于，在这场“舞蹈”的每一步，算法从矩阵 那里唯一需要的就是它的作用。为了计算残差，它需要执行一次矩阵向量乘积 。像共轭梯度（CG）法（用于对称矩阵）或广义最小残差（GMRES）法（用于一般矩阵）这样的算法，完全是围绕这个操作构建的。它们仅使用矩阵向量乘积和向量-向量运算（如内积和加法）来构造一系列搜索方向和最优步长。
例如，著名的Arnoldi 迭代是 GMRES 的核心，它为Krylov 子空间——一个由 张成的空间——构建一个基。为了生成这个基，它唯一需要涉及 的操作就是计算矩阵向量乘积的能力。矩阵本身可以保持为一个幻影，仅由其行为定义。
这种“仅需作用”的算子可能看起来像一个聪明的数学技巧，但在实践中它从何而来？它自然地产生于我们建模物理世界的方式。许多由偏微分方程（PDE）控制的复杂物理系统，都使用有限元方法（FEM）或有限体积法等方法进行分析。
这些方法的核心思想是“分而治之”。一个复杂的物体，如飞机机翼或汽车发动机缸体，在计算上被分解成一个由数百万个称为单元的简单小块（如微小的砖块或金字塔）组成的网格。整个系统的行为是通过理解这些简单单元如何与其邻居相互作用来确定的。
当我们计算全局矩阵 对向量 的作用时，我们不需要先通过考虑所有相互作用来构建 。相反，我们可以对网格中的每一个单元执行一个循环：
当我们遍历完所有单元后，全局向量 将持有 的正确结果。我们完美地再现了全局矩阵的作用，而从未构建过它。物体到单元的物理分解直接映射到了矩阵向量乘积的计算分解。物理模型与计算算法之间的这种美妙统一是现代模拟的基石。
这个原理是如此强大，以至于可以无缝地扩展到非线性问题。在像 Newton-Raphson 方法这样的方法中，必须求解一个涉及切线矩阵（雅可比矩阵）的线性系统。即使是这个复杂的、依赖于状态的矩阵，也可以以无矩阵的方式施加其作用，要么通过解析地推导其作用，要么通过巧妙地使用有限差分来近似其作用。
这一切似乎是一种优雅的计算方式，但这仅仅是品味问题吗？完全不是。无矩阵方法不仅仅是一种奇技淫巧；它们是一种必需品，由计算中两个最基本的约束所驱动：内存和速度。
内存方面的论点是直截了当的。对于大规模三维模拟或在每个单元内使用高阶多项式的方法，矩阵 中非零元素的数量可能会变得天文数字般巨大。存储这个矩阵可能需要太字节（TB）的内存，超出了即使是最大型超级计算机的容量。无矩阵方法通过完全不存储矩阵，优雅地绕过了这道“内存墙”。
速度方面的论点则更为微妙，它揭示了现代计算机架构的一个深刻事实。把一个高性能计算机处理器想象成一个工厂。工厂有一条极其快速的装配线（其计算单元，以 FLOPs，即每秒浮点运算次数来衡量），但用于运送零件的传送带系统相对较慢（其内存带宽）。工厂的效率取决于保持装配线的繁忙。
组装矩阵 (SpMV)： 一个标准的稀疏矩阵向量乘积（SpMV）就像一个装配线工人，每次操作都需要一个独特的、特定的零件。矩阵的值和它们的位置存储在内存的各个角落。为了执行每次乘法，处理器必须从内存中请求一个值，等待缓慢的传送带将其送达，执行一两次快速操作，然后请求下一个。工人大部分时间都在等待零件。这个过程是内存受限的。计算与数据移动的比率，即算术强度，低得可怜。对于一个典型的 SpMV，它可能在每字节数据移动 FLOPs 左右。
无矩阵方法： 一个无矩阵方法，特别是对于高阶单元，就像一个一次性接收一整盘零件的工人（单个单元的数据，这些数据很小且在内存中是连续的）。然后，他们可以在需要下一盘零件之前，对这些本地数据进行大量的计算。这让装配线一直高效运转。这个过程可以变成计算受限的。它的算术强度要高得多，并且至关重要地是，它通常随着方法的复杂性（例如，多项式阶数 ）的增加而增加。达到 FLOPs/字节或更高的算术强度并不少见。
结果是惊人的。在现代处理器上，计算速度远快于内存访问速度，一个计算受限的无矩阵算子可以比其内存受限的组装对应物快 10 倍、100 倍，甚至更多。我们不仅节省了内存；我们还释放了硬件真正的计算能力。
一个快速的矩阵向量乘积很棒，但对于困难的问题，我们需要一个预处理器来确保迭代求解器在合理的步数内收敛。一个预处理器 将系统转换为一个更容易求解的系统，如 。应用预处理器意味着我们需要一种方法来求解与矩阵 相关的系统。
这就提出了一个引人入胜的新难题：如果我们费了这么大劲去避免构建 ，我们又怎么可能构建一个本应近似 的预处理器 呢？然而，这个约束催生了创造力，并导向了更为优雅的解决方案。
策略1：构建一个更简单的幽灵。 我们不能使用像不完全 LU 分解（ILU）或标准代数多重网格（AMG）这样的“黑箱”代数预处理器，因为它们需要访问矩阵元素。这代表了在灵活性上的一个权衡。但我们可以回到源头：物理。我们可以构建一个简化的物理模型——例如，通过忽略次要的物理效应或使用平均的材料属性——并为这个更简单的问题组装矩阵 。这个 就成了我们的预处理器 。这就是基于物理的预处理的精髓，我们近似的是物理，而不仅仅是代数。
策略2：以火攻火。 一个更优美的方法是使用本身就是无矩阵的预处理器。一个典型的例子是多项式预处理器。在这里， 的作用被一个 的多项式所近似，例如 。应用这个预处理器仅仅意味着多次应用 的作用！这完美地融入了我们的计算策略。我们用一系列我们已经完善的高效、计算受限的矩阵向量乘积，来替代一个可能缓慢、内存受限的预处理步骤（如三角求解）。
这导致了一种对算法设计的整体性观点。选择一个无矩阵算子不是一个孤立的决定。它影响预处理器的选择和求解器的整个结构。为了达到峰值性能，应该将一个计算受限的无矩阵算子与一个同样高性能且避免重新引入内存带宽瓶颈的预处理器配对。这就创造了一个良性循环，最终得到一个不仅在数学上合理，而且完美契合现代计算机架构优势的求解器。
当我们从将矩阵视为静态的数字网格转变为将其视为动态的算子——一个接收向量并将其转换为另一个向量的机器时，会发生一种深刻的视角转变。这不仅仅是一个语义游戏；它正是无矩阵哲学的核心。要领略其威力，我们必须踏上一场穿越现代科学与工程领域的旅程，去看看这个单一的思想如何解锁了那些曾被认为庞大到不可能或复杂到无法解决的问题。这不是一个关于蛮力计算的故事，而是一个关于数学优雅和对真正重要事物深刻理解的故事。
想象一下试图理解一台复杂的机器，比如一架大型喷气式客机。一种方法是为每一个零件——每一根电线、每一颗螺丝、每一个铆钉——创建一份完整的蓝图，并将其存放在一个巨大的图书馆里。这就是“显式矩阵”方法。但如果零件数量达到数十亿呢？图书馆将大到不可能。另一种无矩阵的方法是直接与机器互动。我们推动一个杠杆（我们的输入向量），观察控制面如何移动（输出向量）。通过执行一系列巧妙的“戳探”，我们可以在无需完整蓝图的情况下，学到解决问题所需的一切。
这正是计算力学许多领域中的情况。当工程师模拟一座桥梁在负载下的行为或一架飞机机翼上的气流时，他们将物体离散化为数百万甚至数十亿个微小的有限元。这些单元之间的相互作用由一个“切线刚度矩阵”来描述。对于一个几何复杂、非线性的材料，这个矩阵不仅巨大，而且在模拟的每一步都会改变。显式地组装和存储它会消耗太字节（TB）的内存，即使是超级计算机也会因此而停滞。
在这里，Newton-Krylov 方法作为一个无矩阵范式的优美范例前来救场。为了求解非线性平衡方程，该方法用一个小的扰动向量“戳探”系统，并通过简单地用轻微扰动后的状态重新评估系统的控制方程来计算其响应——即雅可比-向量积。这完全避免了构建矩阵，将内存需求从与相互作用数量成正比降低到与单元数量本身成正比。它将一个棘手的内存问题转变为一个可管理的计算问题。
当我们模拟多物理场现象时，挑战变得更加尖锐，例如模拟流体饱和的多孔材料（如土壤或骨骼）的行为——这个领域被称为多孔弹性力学。在这里，固体骨架变形，同时内部的流体流动，它们的行为错综复杂地耦合在一起。系统矩阵呈现出块状结构，不仅描述了固体-固体和流体-流体的相互作用，还描述了关键的固体-流体耦合。对于一个大型三维模拟，存储这个完整矩阵所需的内存可以轻易达到数百吉字节（GB），远超单个计算节点的能力。在这里，无矩阵方法不仅仅是一种便利；它们是一种赋能技术。没有它们，这种规模的模拟根本不可能实现。此外，这种哲学甚至延伸到了预处理——加速迭代求解的艺术——在这里，人们可以设计出强大的、基于物理的无矩阵预处理器，如多重网格方法，它们像“智能戳探”一样引导解快速收敛。
除了求解方程组，科学中最基本的任务之一是寻找算子的特征值和特征向量。这些代表了系统的特殊“模式”或“状态”——小提琴弦的共振频率、分子的振动模式或原子的能级。在许多情况下，我们不需要数百万个所有可能的模式；我们只关心处于极端位置的少数几个——最低的能态、最高的频率，或者最容易失稳的模式。
考虑大规模优化的任务：寻找一个具有数百万变量的函数的最小值，就像调整一个复杂机器学习模型的参数一样。一个驻点只有在景观向各个方向都向上弯曲时，才是一个真正的局部最小值。这个性质被编码在 Hessian 矩阵（二阶导数矩阵）中，它必须是正定的。我们如何能在不构建一个百万乘百万的 Hessian 矩阵的情况下检查这一点呢？我们可以使用一个无矩阵的迭代特征求解器，比如 Lanczos 算法，来寻找它的最小特征值。这个寻找过程所需要的只是一个计算 Hessian-向量积的方法，这通常可以高效完成，这种技术有时被称为“无 Hessian”优化。如果找到的最小特征值为正，我们就可以确信我们找到了一个最小值。
这种对极端特征值的追求在量子化学中得到了最壮观的体现。薛定谔方程，在一个离散的基组中，变成了一个巨大的特征值问题。
在所有这些特征值问题中，核心思想是相同的：我们在一个小的、可管理的子空间中构建一个近似解，我们不是通过了解整个矩阵来完善这个子空间，而是通过迭代地将其作用应用于我们当前的最佳猜测来完善它。
到目前为止，我们将无矩阵方法视为解决过大问题的一种方式。但出人意料的是，它们也成为了使计算更快的关键，即使在矩阵可以被存储的情况下也是如此。原因在于现代计算机的架构，特别是图形处理单元（GPU）。
现代 GPU 是一个计算巨擘，每秒能够执行数万亿次计算。然而，它与内存的连接是一个相对的瓶颈。这就像一个工厂，拥有一条快得惊人的装配线，却总是在等待由慢得多的传送带运送零件。这种关系被“屋顶线模型”所描述，它告诉我们性能要么受限于计算速度，要么受限于内存带宽。关键指标是算术强度——执行的计算量与从内存移动的数据量之比。
使用显式的、存储的稀疏矩阵（SpMV）的传统方法具有非常低的算术强度。对于从矩阵中读取的每一个数字，计算机只执行一两个操作。它长期处于数据“饥饿”状态，即内存受限，因此无法达到其峰值性能。
现在，考虑高阶有限元方法，如等几何分析（IGA），它使用光滑、复杂的基函数来达到高精度。在这里，无矩阵方法不存储预先计算的矩阵元素。相反，它利用被称为求和-分解的巧妙张量代数，在每个单元内部即时重新计算算子的作用。这种方法对它从内存中读取的每块数据执行更多的计算。它的算术强度很高，并且至关重要地是，它随着方法的复杂性（多项式次数 ）的增加而增长。对于一个足够高阶的方法，无矩阵核心变得计算受限。它有足够的工作来让 GPU 强大的处理器保持满负荷运转。令人震惊的结果是，即使无矩阵方法可能执行更多的总浮点运算，但它可能快上几个数量级，因为它充分利用了硬件的潜力。这是数学算法与计算机架构的美妙共同进化。
我们的旅程展示了矩阵向量乘积的力量，但留下了一个关键问题未解：我们从哪里得到这个乘积？对于复杂的非线性问题，像我们在高分子物理例子中看到的那样，手工推导雅可比-向量积的表达式可能是一项艰巨且易于出错的任务。
这就是与计算机科学的深刻联系提供了一个近乎神奇的解决方案的地方：自动微分（AD）。想象一下，你可以教你的计算机微积分的规则。当它执行计算你的物理模型方程的程序时，它可以同时、精确地计算出导数。这就是 AD 所做的。
通过使用一种特殊的数字类型来增强代码，这种数字类型同时携带一个值和一个导数（一个“对偶数”），我们可以使用前向模式 AD。当我们运行我们的模拟代码时，我们用我们的扰动向量 来“播种”输入的导数部分。然后，代码根据链式法则将导数传播到每一个数学运算中，最终输出的导数部分恰好就是雅可比-向量积 。
更强大的是，我们可以使用反向模式 AD，这与驱动深度神经网络训练的技术相同（在那里它被称为反向传播）。这使我们能够计算向量-雅可比积 ，其计算成本与输入变量的数量惊人地无关。这个操作是*伴随方法*的核心，对于大规模灵敏度分析、数据同化和优化至关重要。
其意义是惊人的。我们可以为我们复杂的物理模拟编写一次代码。然后，AD 框架提供了快速 Newton-Krylov 求解器（需要 ）和强大的基于伴随的优化算法（需要 ）所需的“无矩阵”机制，所有这一切都无需我们编写一行解析导数代码。它代表了物理模拟、数值分析和编译器技术的宏大统一，让科学家能够向他们的模型提出越来越复杂的问题。
从处理大到不可能的系统，到寻找量子态，再到在超级计算机上最大化性能，无矩阵哲学揭示了它不是一个单一的技巧，而是一种强大而统一的思维方式。它教导我们关注算子的动态作用而非其静态表示，这种视角的转变为计算科学不断拓展边界提供了持续的动力。