【读论文】Closed-loop Diffusion Control of Complex Physical Systems 闭环扩散控制系统

摘要

复杂物理系统的控制问题在科学和工程中具有广泛的应用。先前的研究表明，基于扩散模型的生成式控制方法在解决这些问题方面具有显著优势。然而，当将现有的生成式控制方法扩展到闭环控制场景时（这对实现有效控制至关重要），其性能和效率均面临挑战。在本文中，我们提出了一种高效的物理系统闭环扩散控制方法（CL-DiffPhyCon）。通过为不同物理时间步引入异步去噪框架，CL-DiffPhyCon 在采样过程中显著降低了计算成本，同时根据环境的实时反馈生成控制信号。此外，通过引入快速采样技术（如 DDIM），控制过程可以进一步加速。我们在两个任务上对 CL-DiffPhyCon进行了评估：一维 Burgers 方程控制和二维不可压缩流体控制。结果表明，CL-DiffPhyCon 在控制性能和采样效率方面均表现优越。

1 引言

复杂物理系统的控制问题是一个关键研究领域，涉及通过优化一系列控制动作来实现特定目标。其在科学和工程领域具有广泛的应用，例如流体控制（Verma等，2018）、等离子体控制（Degrave等，2022）和粒子动力学控制（Reyes Garza等，2023）。控制此类系统的难点在于其高维、高非线性以及随机特性。因此，为了实现良好的控制性能，闭环控制是必不可少的。具体来说，每次控制决策都应基于环境的最新状态反馈，以便能够根据任何变化持续调整控制输入。对于存在额外挑战（如部分观测或间接控制）的控制任务，这种要求尤为重要。

过去几十年中，为解决该问题开发了多种方法，包括经典控制方法、近年的强化学习方法以及最新的生成式方法。其中，扩散模型（如 DiffPhyCon，Wei 等，2024）在复杂物理系统控制中表现出竞争力，通常优于经典控制和强化学习方法。扩散模型在一般决策问题上的优越性也已被广泛证明（Ajay等，2022；Janner等，2022a）。然而，由于依赖同步去噪策略，这些扩散控制方法在处理闭环控制问题时面临显著挑战。扩散模型从纯噪声开始，为模型时域内的所有物理时间步生成去噪样本。在每个物理时间步执行完整的采样过程可以实现闭环控制，但采样成本较高。此外，这可能破坏控制信号的一致性，从而影响整体性能（Kaelbling & Lozano-Perez，2011）。另一方面，每隔几个物理时间步进行完整采样可以提高采样效率，但不再符合闭环要求，因而由于系统状态的滞后导致较差的控制决策。尽管最近提出了一种在线重新规划策略（Zhou 等，2024），用于自适应地确定重新规划的时机，但它并未建立完全的闭环框架。此外，这种方法涉及额外的似然估计或从头采样计算，并依赖于任务特定的超参数阈值，需要通过实验或大量调整确定。

在本文中，我们提出了一种新颖的 物理系统闭环扩散控制方法，命名为 CL-DiffPhyCon。其 关键思想是 将模型时域内的同步去噪解耦，允许不同物理时间步呈现不同的噪声水平。通过这种方式，自然实现了在线闭环控制序列生成：异步扩散模型按物理时间步递增的噪声水平逐步输出控制信号，使得在不等待时域内所有后续控制信号完全去噪的情况下，可以实现环境与采样控制信号的实时交互。然后，环境反馈作为后续控制信号采样的初始条件，确保生成基于可靠状态的控制信号。我们的方法也可以看作是一种 无缝的重新规划方法，以最小的采样成本利用环境的新鲜观测。因此，与现有基于扩散的控制方法相比（Wei 等，2024；Zhou 等，2024），我们的方法不仅实现了闭环控制，还显著加速了采样过程。本文方法的优势在图 1 中进行了说明。

图 1：CL-DiffPhyCon 相较以往扩散控制方法的优势。图中展示了我们的 CL-DiffPhyCon 方法（右侧）相较于以往扩散控制方法（左侧和中间）的优势。通过采用异步去噪框架，我们的方法不仅能够实现闭环控制，还能显著加速采样过程。图中的符号 DiffPhyCon-h 表示每 $h$ 个物理时间步执行一次包含 $T$ 个去噪步骤的完整采样过程。

总结而言，我们的贡献包括：（1）开发了一种新颖的复杂物理系统闭环扩散控制方法 CL-DiffPhyCon。通过对目标分布的理论分析，提出了一种闭环框架内加速扩散控制的方法。此外，通过引入快速采样技术（如 DDIM，Song 等，2020），进一步加速了控制过程。（2）我们在一维 Burgers 方程控制和二维不可压缩流体控制任务中对 CL-DiffPhyCon 进行了评估。结果表明，CL-DiffPhyCon 在控制性能和采样加速方面表现出显著优势。

2 相关工作

经典控制方法（如 PID，Li 等，2006）以其高效性、稳定性能和良好的可解释性著称，但在应用于高维复杂物理系统的控制时，在性能和效率方面面临显著挑战。近年来，模仿学习和强化学习（Pomerleau，1988；Zhuang 等，2023）在包括减阻（Rabault 等，2019；Elhawary，2020；Feng 等，2023；Wang 等，2024）、传热（Beintema 等，2020；Hachem 等，2021）以及鱼类游动（Novati 等，2017；Verma 等，2018；Feng 等，2024）等广泛的物理系统中表现出色。另一类监督学习（SL）方法（Holl 等，2020；Hwang 等，2022）通过利用神经代理模型的反向传播规划控制信号。相比之下，我们的方法并不依赖代理模型，而是同时学习物理系统的动力学和相应的控制序列。此外，物理信息神经网络（PINNs）（Willard 等，2020）最近也被用于控制（Mowlavi & Nabi，2023），但它们需要显式的偏微分方程（PDEs）。相比之下，我们的方法是数据驱动的，只要环境可访问，就能够解决范围更广的复杂物理系统控制问题。

扩散模型（Ho 等，2020）在学习高维分布方面表现优异，并在图像和视频生成领域取得了显著成功（Dhariwal & Nichol，2021；Ho 等，2022）。它们还展示了在解决科学和工程挑战中的卓越能力，例如天气预测（Price 等，2023）、机器人控制（Janner 等，2022b；Ajay 等，2022），以及偏微分方程（PDEs）中的反问题，例如设计初始或边界条件（Wu 等，2024）或规划外部时间控制信号（Wei 等，2024）。然而，由于所有物理时间步都在同一计划下去噪，它们不适合直接应用于高效的闭环控制。一些先前的工作聚焦于改进扩散生成在决策任务中的适应性（Zhou 等，2024），但这些方法并非闭环，且在重新规划的决策中需要额外的超参数和计算。相比之下，我们的方法是一种闭环方法，具有高效的采样策略。

一些最近的研究也为模型时域内的不同帧分配了不同的噪声水平（Wu 等，2023；Ruhe 等，2024）。我们的工作与这些方法的关键区别有两点：（1）我们的方法专注于复杂物理系统的闭环控制任务，而他们的方法则面向不涉及与外部环境交互的顺序内容生成。（2）我们需要学习的扩散模型来源于目标分布（请参见第 4.1 节了解详细信息），而他们的方法是启发式的。

3 背景

3.1 问题设定

我们考虑以下在线复杂物理系统控制问题：

$$\underset{{\mathbf{w}}_{1:N}}{\min }\mathcal{J}({\mathbf{u}}_0,{\mathbf{w}}_1,{\mathbf{u}}_1,\dots ,{\mathbf{w}}_N,{\mathbf{u}}_N)\text{s.t.}{\mathbf{u}}_{\tau +1}=G({\mathbf{u}}_{\tau },{\mathbf{w}}_{\tau +1}).\text{\hspace{1em}(1)}$$

这里，${\mathbf{u}}_{\tau }\in {\mathbb{R}}^{d_u}$ 和 ${\mathbf{w}}_{\tau }\in {\mathbb{R}}^{d_w}$ 分别表示物理时间步 $\tau$ 的系统状态和外部控制信号；$G$ 是由物理动力学（如偏微分方程，PDEs）决定的环境函数。$J({\mathbf{u}}_0,{\mathbf{w}}_1,{\mathbf{u}}_1,\dots ,{\mathbf{w}}_N,{\mathbf{u}}_N)$ 表示 定义在长度为 $N$ 的轨迹上的控制目标。例如，$\mathcal{J}$ 可以被设计为衡量系统状态 ${\mathbf{u}}_N$ 向目标状态 ${\mathbf{u}}^{\ast }$ 收敛的控制性能，并包含代价约束：$\mathcal{J}=\Vert {\mathbf{u}}_N-{\mathbf{u}}^{\ast }{\Vert }^2+{\sum }_{\tau =1}^N\Vert {\mathbf{w}}_{\tau }{\Vert }^2.$ 为简便起见，我们引入变量 ${\mathbf{z}}_{\tau }=[{\mathbf{w}}_{\tau },{\mathbf{u}}_{\tau }]$，表示 ${\mathbf{w}}_{\tau }$ 和 ${\mathbf{u}}_{\tau }$ 的拼接。此外，我们假设系统的转移概率具有以下马尔科夫性质：
$$p({\mathbf{z}}_{\tau +1}|{\mathbf{u}}_0,{\mathbf{z}}_{1:\tau })=p({\mathbf{z}}_{\tau +1}|{\mathbf{u}}_{\tau })=p({\mathbf{w}}_{\tau +1}|{\mathbf{u}}_{\tau })p({\mathbf{u}}_{\tau +1}|u_{\tau },{\mathbf{w}}_{\tau +1}).\text{\hspace{1em}(2)}$$

本文关注闭环控制，即每个时间步的控制信号 ${\mathbf{w}}_{\tau +1}$ 是基于当前系统状态 ${\mathbf{u}}_{\text{env},\tau }$ 规划的。闭环控制的优势显而易见：环境的实时反馈使算法能够根据环境状态的变化调整控制决策。

3.2 预备知识：扩散控制模型

DiffPhyCon（Wei 等，2024）是最近提出的一种扩散生成方法，用于以开环方式解决小 $N$ 情况下的问题（式 1）。本节简要回顾其简化版本的框架。假设我们有一个训练集 ${\mathcal{D}}_{\text{train}}$，包含来自环境的 $M$ 条轨迹 { u 0 ( i ) , z 1 : N ( i ) } i = 1 M \{\mathbf u_0^{(i)}, \mathbf z_{1:N}^{(i)}\}_{i=1}^M {u0(i)​,z1:N(i)​}i=1M​。我们在 ${\mathcal{D}}_{\text{train}}$ 上定义以下正向扩散随机微分方程（SDE，Song 等，2021）：

$$\text{d}{\mathbf{z}}_{\tau }(t)=f(t){\mathbf{z}}_{\tau }(t)dt+g(t)d{\omega }_{\tau }(t),\tau \in [1:N],t\in [0,T],\text{\hspace{1em}(3)}$$

其中 ${\mathbf{z}}_{\tau }(0)={\mathbf{z}}_{\tau }$，$f(t)$ 和 $g(t)$ 是标量函数，${\omega }_{\tau }(t)$ 是 Wiener 过程。通过式（3），我们将 ${\mathcal{D}}_{\text{train}}$ 的分布 $p_{\text{train}}({\mathbf{u}}_0,{\mathbf{z}}_{1:N})=p_{\text{train}}({\mathbf{u}}_0)p_{\text{train}}({\mathbf{z}}_{1:N}|{\mathbf{u}}_0)$ 扩展为扩散过程的分布 $p_t({\mathbf{u}}_0,{\mathbf{z}}_{1:N}(t))=p_{\text{train}}({\mathbf{u}}_0)p_t({\mathbf{z}}_{1:N}(t)|{\mathbf{u}}_0)$。终端条件为 $p_0({\mathbf{z}}_{1:N}(0)|{\mathbf{u}}_0)=p_{\text{train}}({\mathbf{z}}_{1:N}(0)|{\mathbf{u}}_0)$ 和 $p_T({\mathbf{z}}_{1:N}(T)|{\mathbf{u}}_0)\approx \mathcal{N}({\mathbf{z}}_{1:N}(T)|0,{\sigma }_T^2\mathbf{I})$。

式（3）的反向时间 SDE 为：

$$\text{d}{\mathbf{z}}_{1:N}(t)=[f(t){\mathbf{z}}_{1:N}(t)-g(t{)}^2\nabla \log p_t({\mathbf{z}}_{1:N}(t)|{\mathbf{u}}_0)]dt+g(t)d{\omega }_{1:N}(t),t\in [T,0].\text{\hspace{1em}(4)}$$

一旦获得分数函数 $\nabla \log p_t({\mathbf{z}}_{1:N}(t)|{\mathbf{u}}_0)$，便可以通过式（4）从 ${\mathbf{z}}_{1:N}(T)\sim \mathcal{N}(0,I)$ 采样得到 ${\mathbf{z}}_{1:N}(0)\sim p_{\text{train}}({\mathbf{z}}_{1:N}(0)|{\mathbf{u}}_0)$。根据 DiffPhyCon，通过解决反问题（Song 等，2021；Chung 等，2023），可实现控制目标 $\mathcal{J}({\mathbf{z}}_{1:N}(0))$：

$$\text{d}{\mathbf{z}}_{1:N}(t)=[f(t){\mathbf{z}}_{1:N}(t)-g(t{)}^2{\nabla }_{{\mathbf{z}}_{1:N}(t)}\log p_t({\mathbf{z}}_{1:N}(t)|{\mathbf{u}}_0)+g(t{)}^2\lambda \cdot {\nabla }_{{\mathbf{z}}_{1:N}(t)}\mathcal{J}({\hat{\mathbf{z}}}_{1:N}(0))]dt+g(t)d{\omega }_{1:N}(t).\text{\hspace{1em}(5)}$$

其中，${\tilde{\mathbf{z}}}_{1:N}(0)$ 是 Tweedie 估计，定义为：

$${\tilde{\mathbf{z}}}_{1:N}(0)=s(t{)}^{-1}{\mathbf{z}}_{1:N}(t)+s(t)\sigma (t{)}^2{\nabla }_{{\mathbf{z}}_{1:N}(t)}\log p_t({\mathbf{z}}_{1:N}(t)|{\mathbf{u}}_0),$$

其中 $s(t)=\exp \left({\int }_0^{t}f(\eta )d\eta \right)$ 和 $\sigma (t)=\sqrt{{\int }_0^t\frac{g(\eta {)}^2}{s(\eta {)}^2}d\eta }$。在实际实现中，这些量通过 VP SDE（Song 等，2021）专门化，即 DDPM（Ho 等，2020），扩散模型通过一个水平为 $H$ 的参数化去噪网络 ${ϵ}_{\varphi }$ 实现，其水平等于 $N$。更多细节请参考附录 A 和 B。

在 ${\mathbf{z}}_{\tau }(t)$ 中，下标 $\tau$ 表示物理时间步，而括号中的 $(t)$ 表示扩散随机微分方程（SDE）的时间步。

4 方法

在本节中，我们详细介绍了我们的方法 CL-DiffPhyCon。在第 4.1 节中，我们阐述了我们的基本思想，并推导出了需要学习的两个分布。为了从这些分布中采样，我们分别在第 4.2 节和第 4.3 节中提出了同步扩散模型和异步扩散模型。在第 4.4 节中，我们介绍了 CL-DiffPhyCon 实现的闭环控制（推理），这一过程在图 2 中也得到了说明。

图 2：用于闭环控制的 CL-DiffPhyCon。CL-DiffPhyCon 首先使用同步扩散模型进行初始化，然后迭代地利用异步扩散模型执行在线控制。每个控制信号的采样基于来自环境的最新系统状态反馈。

4.1 异步去噪框架

回顾 DiffPhyCon（Wei 等人，2024），该方法在扩散逆过程的时间范围 $H$（通常远短于 $N$）中需要对潜变量（如系统状态和控制信号）进行同步去噪。因此，为了采样控制信号，必须在整个时间范围上执行一个长度为 $T$ 的完整去噪过程，这对该时间范围内的所有潜变量引入了大量计算开销。为了解决这一问题，我们提出了 CL-DiffPhyCon 方法，采用一种异步去噪流程，使早期物理时间步的潜变量比后续时间步的潜变量提前去噪。在每个物理时间步，采样的控制信号被输入环境，其输出状态作为后续去噪过程的初始条件。因此，控制信号基于当前系统状态进行规划，实现了闭环控制。同时，与同步采样相比，连续两个时间步之间的计算成本显著降低。

形式化地，我们的目标是对马尔科夫复杂物理系统中的联合分布 $p({\mathbf{z}}_1(0),{\mathbf{z}}_2(0),\cdots ,{\mathbf{z}}_N(0)|{\mathbf{u}}_0)$ 进行建模，其中每个 ${\mathbf{z}}_{\tau }(0)=[{\mathbf{w}}_{\tau },{\mathbf{u}}_{\tau }]$ 是无噪声的控制信号和系统状态对，${\mathbf{u}}_0$ 是初始状态。定义 ${\mathbf{z}}_{\tau :\tau +H-1}(t)=[{\mathbf{z}}_{\tau }(t),{\mathbf{z}}_{\tau +1}(t),\cdots ,{\mathbf{z}}_{\tau +H-1}(t)]$ 为水平区间 $[\tau :\tau +H-1]$ 内具有相同噪声水平的潜变量序列，${\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)=[{\mathbf{z}}_{\tau }(t),{\mathbf{z}}_{\tau +1}(t+\frac{1}{H}T),\cdots ,{\mathbf{z}}_{\tau +H-1}(t+\frac{H-1}{H}T)]$ 为其对应的异步噪声水平版本。我们可以考虑以下扩展联合分布：

$$p({\mathbf{z}}_{1:N}(0),{\tilde{\mathbf{z}}}_{1:H}(\frac{1}{H}T),\cdots ,{\tilde{\mathbf{z}}}_{N+1:N+H}(\frac{1}{H}T)|{\mathbf{u}}_0)\text{\hspace{1em}(6)}$$

通过从该扩展联合分布中采样，我们可以获得期望的控制信号和状态序列 ${\mathbf{z}}_{1:N}(0)\sim p_{\text{train}}({\mathbf{z}}_{1:N}(0)|{\mathbf{u}}_0)$。通过依次对先前的变量进行条件化，我们得到以下分解定理：

定理 1. 假设联合分布 $p({\mathbf{z}}_1(0),{\mathbf{z}}_2(0),\cdots ,{\mathbf{z}}_N(0)|{\mathbf{u}}_0)$ 满足马尔科夫性质。对于任意 $\tau >0$，假设 ${\mathbf{z}}_{\tau }(T)$ 独立地服从正态分布 $N({\mathbf{z}}_{\tau }(T)|0,{\sigma }_T^2\mathbf{I})$。则扩展联合分布可以分解为：

$$p({\mathbf{z}}_{1:N}(0),{\tilde{\mathbf{z}}}_{1:H}(\frac{1}{H}T),\cdots ,{\tilde{\mathbf{z}}}_{N+1:N+H}(\frac{1}{H}T)|{\mathbf{u}}_0)=p({\tilde{\mathbf{z}}}_{1:H}(\frac{1}{H}T)|{\mathbf{u}}_0)\prod _{\tau =1}^{N}p({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(0)|{\mathbf{u}}_{\tau -1}(0),{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(\frac{1}{H}T))\mathcal{N}({\mathbf{z}}_{\tau +H}(T);0,{\sigma }_T^2\mathbf{I})\text{\hspace{1em}(7)}$$

该定理及后续命题的证明在附录 C 中给出。定理表明，为了从扩展联合分布中采样，我们只需要指定两个类型的分布用于祖先采样：$p({\tilde{\mathbf{z}}}_{1:H}(\frac{1}{H}T)|{\mathbf{u}}_0)$ ，和 $p({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(0)|{\mathbf{u}}_{\tau -1}(0),{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(\frac{1}{H}T))$ 对于 $\tau >0$。也就是说，我们需要一个扩散模型来学习初始化分布 $p({\tilde{\mathbf{z}}}_{1:H}(\frac{1}{H}T)|{\mathbf{u}}_0)$，另一个扩散模型来学习转移分布 $p({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(0)|{\mathbf{u}}_{\tau -1}(0),{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(\frac{1}{H}T))$。

4.2 从初始化分布中采样

为了从分布 $p({\hat{\mathbf{z}}}_{1:H}(\frac{1}{H}T)|{\mathbf{u}}_0)$ 中采样，我们通过以下损失函数训练一个同步扩散模型 ${ϵ}_{\varphi }$：

$${\mathcal{L}}_{\text{synch}}={\mathbb{E}}_{t,({\mathbf{u}}_0,{\mathbf{z}}_{1:H}),ϵ}\left[\Vert ϵ-{ϵ}_{\varphi }({\mathbf{z}}_{1:H}(t),{\mathbf{u}}_0,t){\Vert }_2^2\right].\text{\hspace{1em}(8)}$$

其中，${\mathbf{z}}_{1:H}(t)=s(t){\mathbf{z}}_{1:H}+s(t{)}^2\sigma (t{)}^2ϵ$，并且期望值取决于 $t\sim U(0,T)$、$({\mathbf{u}}_0,{\mathbf{z}}_{1:H})\sim {\mathcal{D}}_{\text{train}}$、$ϵ\sim N(0,I)$。由于负梯度得分函数会最小化该损失，我们得到：

$${\nabla }_{{\mathbf{z}}_{1:H}(t)}\log p_t({\mathbf{z}}_{1:H}(t)|{\mathbf{u}}_0)\approx -{ϵ}_{{\varphi }^{\ast }}({\mathbf{z}}_{1:H}(t),{\mathbf{u}}_0,t),$$

其中 ${\varphi }^{\ast }$ 是经过训练的参数。完成训练后，我们可以应用公式 (4) 来采样 { z 1 : H ( t ) } \{\mathbf z_{1:H}(t)\} {z1:H​(t)}，其中 $t$ 从 $T$ 到 $\frac{1}{H}T$ 递减，从而得到：${\tilde{\mathbf{z}}}_{1:H}(\frac{1}{H}T)=[{\mathbf{z}}_1(\frac{1}{H}T),\cdots ,{\mathbf{z}}_H(T)]$，此时，这组样本遵循所需的分布 $p({\tilde{\mathbf{z}}}_{1:H}(\frac{1}{H}T)|{\mathbf{u}}_0)$。

4.3 从过渡分布中采样

为了从分布 $p({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(0)|{\mathbf{u}}_{\tau -1}(0),{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(\frac{1}{H}T))$ 中采样，我们训练了一个异步扩散模型 ${ϵ}_{\theta }$。根据第 3.2 节中的符号，我们定义了一个异步前向扩散 SDE：
$$\begin{array}{rl}\text{d}{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)& =[\text{d}{\mathbf{z}}_{\tau }(t),\text{d}{\mathbf{z}}_{\tau +1}(t+\frac{1}{H}T),\cdots ,\text{d}{\mathbf{z}}_{\tau +H-1}(t+\frac{H-1}{H}T)]\\ & ={\tilde{f}}_{\tau :\tau +H-1}(t){\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)dt+{\tilde{g}}_{\tau :\tau +H-1}(t)d{\omega }_{\tau :\tau +H-1}(t),\end{array}\text{\hspace{1em}(9)}$$

其中 ${\tilde{f}}_{\tau :\tau +H-1}(t)=[f(t),f(t+\frac{1}{H}T),\cdots ,f(t+\frac{H-1}{H}T)]$，${\tilde{g}}_{\tau :\tau +H-1}(t)=[g(t),g(t+\frac{1}{H}T),\cdots ,g(t+\frac{H-1}{H}T)]$。方程 (9) 的反向时间 SDE 表达式为：
$$\begin{array}{rl}d{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)& =[{\tilde{f}}_{\tau :\tau +H-1}(t){\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)\\ & -{\tilde{g}}_{\tau :\tau +H-1}(t{)}^2{\nabla }_{{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)}\log p_t({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)|u_{\tau -1}(0))]dt\\ & +{\tilde{g}}_{\tau :\tau +H-1}(t)d{\omega }_{\tau :\tau +H-1}(t)\end{array}\text{\hspace{1em}(10)}$$

进一步地，包含控制目标的推断可推广为：
$$\begin{array}{rl}d{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)& =[{\tilde{f}}_{\tau :\tau +H-1}(t){\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)\\ & -{\tilde{g}}_{\tau :\tau +H-1}(t{)}^2{\nabla }_{{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)}\log p_t({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)|u_{\tau -1}(0))\\ & +{\tilde{g}}_{\tau :\tau +H-1}(t{)}^2\lambda {\nabla }_{{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)}\mathcal{J}({\hat{\mathbf{z}}}_{\tau :\tau +H-1}(0))]dt\\ & +{\tilde{g}}_{\tau :\tau +H-1}(t)d{\omega }_{\tau :\tau +H-1}(t)\end{array}\text{\hspace{1em}(11)}$$

这里，${\hat{\mathbf{z}}}_{\tau :\tau +H-1}(0)$ 是一个异步 Tweedie 估计（详见附录 A）：
$$\begin{array}{rl}{\hat{\mathbf{z}}}_{\tau +i}(0)& =s(t+\frac{i}{H}T{)}^{-1}{\mathbf{z}}_{\tau +i}(t+\frac{i}{H}T)\\ & +s(t+\frac{i}{H}T)\sigma (t+\frac{i}{H}T{)}^2\frac{\partial \log p_t({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)|{\mathbf{u}}_{\tau -1}(0))}{\partial {\mathbf{z}}_{\tau +i}(t+\frac{i}{H}T)}\end{array}$$

${\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)$ 的分布满足以下命题：

命题 1：假设联合分布 $p({\mathbf{z}}_1(0),{\mathbf{z}}_2(0),\cdots ,{\mathbf{z}}_N(0)|{\mathbf{u}}_0)$ 具有马尔科夫性质。对于任何 $t\in [0,\frac{1}{H}T]$，我们有：
$$p_t({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)|{\mathbf{u}}_{\tau -1}(0))={\mathbb{E}}_{{\mathbf{z}}_{\tau :\tau +H-1}(0)}[\prod _{i=0}^{H-1}p({\mathbf{z}}_{\tau +i}(t+\frac{i}{H}T)|{\mathbf{z}}_{\tau +i}(0)\left)\right],\text{\hspace{1em}(12)}$$

其中期望值关于 ${\mathbf{z}}_{\tau :\tau +H-1}(0)\sim p_{\text{train}}({\mathbf{z}}_{\tau :\tau +H-1}(0)|{\mathbf{u}}_{\tau -1}(0))$，并且
$$p({\mathbf{z}}_{\tau +i}(t+\frac{i}{H}T)|z_{\tau +i}(0))=\mathcal{N}({\mathbf{z}}_{\tau +i}(t+\frac{i}{H}T);s(t+\frac{i}{H}T){\mathbf{z}}_{\tau +i}(0),s(t+\frac{i}{H}T{)}^2\sigma (t+\frac{i}{H}T{)}^2\mathbf{I})$$

因此，我们可以通过以下损失函数获取 ${\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)$：
$${\mathcal{L}}_{\text{asynch}}={\mathbb{E}}_{\tau ,t,({\mathbf{u}}_{\tau -1}(0),{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)),ϵ}[\Vert ϵ-{ϵ}_{\theta }({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t),{\mathbf{u}}_{\tau -1},t\left){\Vert }_2^2\right]\text{\hspace{1em}(13)}$$

式 (13) 中的期望取决于 $\tau \sim U(1,N-H+1)$、$t\sim U(0,\frac{1}{H}T)$、$({\mathbf{u}}_{\tau -1}(0),{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t))\sim p_{\text{train}}({\mathbf{u}}_{\tau -1}(0))p_t({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)|{\mathbf{u}}_{\tau -1}(0))$，以及 $ϵ\sim \mathcal{N}(0,I)$。

完成训练后，我们可以通过以下公式计算得分函数：
$${\nabla }_{{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)}\log p_t({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t)|{\mathbf{u}}_{\tau -1}(0))\approx -{ϵ}_{{\theta }^{\ast }}({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(t),{\mathbf{u}}_{\tau -1}(0),t)，$$

其中 ${\theta }^{\ast }$ 是训练好的参数。

因此，我们可以使用公式 (10)，其中 $t$ 从 $\frac{1}{H}T$ 到 0，来采样
$${\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(0)\sim p({\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(0)|{\mathbf{u}}_{\tau -1}(0),{\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(\frac{1}{H}T))$$

4.4 闭环控制

在推理（控制）过程中，我们的目标不仅是生成符合物理动力学的控制序列和系统轨迹，还要优化控制目标。基于两个学习的模型，我们现在介绍闭环控制流程 ${ϵ}_{\varphi }$ 和 ${ϵ}_{\theta }$。首先，我们使用同步模型 ${ϵ}_{\varphi }$ 在区间 $[1,H]$ 上，通过公式 (5) 生成条件于初始状态 ${\mathbf{u}}_{\text{env},0}$ 的初始异步变量 ${\tilde{\mathbf{z}}}_{1:H}\left(\frac{1}{H}T\right)$。随后，在线控制过程开始。

在每个物理时间步 $\tau \ge 1$，我们使用公式 (11) 利用 ${\mathbf{u}}_{\text{env},\tau -1}$ 和 ${\tilde{\mathbf{z}}}_{\tau :\tau +H-1}\left(\frac{1}{H}T\right)$ 采样 ${\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(0)$。然后从采样得到的 ${\mathbf{z}}_{\tau }(0)$ 中提取控制信号 ${\mathbf{w}}_{\tau }(0)$，并将对 $({\mathbf{u}}_{\text{env},\tau -1},{\mathbf{w}}_{\tau }(0))$ 的输入提供给环境 $G$，该环境输出下一个状态 ${\mathbf{u}}_{\text{env},\tau }$。接着，我们采样一个噪声 ${\mathbf{z}}_{\tau +H}(T)\sim \mathcal{N}(0,{\sigma }_T^2\mathbf{I})$，并将其与 ${\tilde{\mathbf{z}}}_{\tau :\tau +H-1}(0)$ 的后 $H-1$ 个分量连接起来，组成 ${\tilde{\mathbf{z}}}_{\tau +1:\tau +H}\left(\frac{1}{H}T\right)$。此时，我们将 ${\mathbf{u}}_{\text{env},\tau -1}$ （而非采样得到的 ${\mathbf{u}}_{\tau -1}(0)$）和 ${\tilde{\mathbf{z}}}_{\tau +1:\tau +H}\left(\frac{1}{H}T\right)$ 传递到下一次循环。完整过程见算法 1，并在图 2 中进行了说明。因此， { w τ ( 0 ) , u env , τ } τ = 1 N \{\mathbf w_\tau(0),\mathbf u_{\text{env}, \tau}\}_{\tau=1}^N {wτ​(0),uenv,τ​}τ=1N​ 的联合分布满足以下命题。

命题 2 根据算法 1 描述的推理过程，有以下公式：

p ( w 1 ( 0 ) , u env , 1 , ⋯ , w N ( 0 ) , u env , N ∣ u env , 0 ) = ∫ p g d ( z ~ 1 : H ( 1 H T ) ∣ u env , 0 ) ∏ τ = 1 N p g d ( z ~ τ : τ + H − 1 ( 0 ) ∣ u env , τ − 1 , z ~ τ : τ + H − 1 ( 1 H T ) ) ⋅ p G ( u env , τ ∣ u env , τ − 1 , w τ ( 0 ) ) N ( z τ + H ( T ) ; 0 , σ T 2 I ) d { u 0 ( 0 ) , z ~ 1 : H ( 1 H T ) , ⋯ , u N ( 0 ) , z ~ N + 1 : N + H ( 1 H T ) } , (14) \begin{aligned}p(\mathbf w_1(0), \mathbf u_{\text{env}, 1}, \cdots, \mathbf w_N(0), \mathbf u_{\text{env}, N} | \mathbf u_{\text{env}, 0}) = \int p_{gd} ( \tilde{\mathbf z}_{1:H}(\frac{1}{H}T) | \mathbf u_{\text{env}, 0}) & \\ \prod_{\tau=1}^N p_{gd}(\tilde{\mathbf z}_{\tau:\tau+H-1}(0) | \mathbf u_{\text{env}, \tau-1}, \tilde{\mathbf z}_{\tau:\tau+H-1}(\frac{1}{H}T) ) \cdot p_G(\mathbf u_{\text{env}, \tau} | \mathbf u_{\text{env}, \tau-1}, \mathbf w_\tau(0)) &\\ \mathcal{N}(\mathbf z_{\tau+H}(T); 0, \sigma_T^2 \mathbf I) \, \text d\{\mathbf u_0(0), \tilde{\mathbf z}_{1:H}\left(\frac{1}{H}T\right), \cdots, \mathbf u_N(0), \tilde{\mathbf z}_{N+1:N+H}\left(\frac{1}{H}T\right)\}, \end{aligned} \tag{14} p(w1​(0),uenv,1​,⋯,wN​(0),uenv,N​∣uenv,0​)=∫pgd​(z~1:H​(H1​T)∣uenv,0​)τ=1∏N​pgd​(z~τ:τ+H−1​(0)∣uenv,τ−1​,z~τ:τ+H−1​(H1​T))⋅pG​(uenv,τ​∣uenv,τ−1​,wτ​(0))N(zτ+H​(T);0,σT2​I)d{u0​(0),z~1:H​(H1​T),⋯,uN​(0),z~N+1:N+H​(H1​T)},​​(14)

其中，$p_{gd}$ 表示引导采样的转移分布（公式 11），而 $p_G$ 表示环境 $G$ 的转移分布。

根据公式 (14)，对于任意 $\tau$，控制信号 ${\mathbf{w}}_{\tau }(0)$ 是条件于环境状态 ${\mathbf{u}}_{\text{env},0:\tau -1}$ 的，而不是预测的 ${\mathbf{u}}_{0:\tau -1}(0)$。因此，我们的方法实现了闭环控制。显然，CL-DiffPhyCon 比 DiffPhyCon-h 快 $H/h$ 倍，其中 DiffPhyCon-h 是一种在每 $h$ 个物理时间步重新运行一次 DiffPhyCon 采样过程（区间为 $H$）的控制方法，如图 1 所示。即便与自适应重规划扩散方法（Zhou 等人，2024）相比，CL-DiffPhyCon 仍然更高效，因为每个潜变量仅采样一次，且不涉及额外计算（如似然估计）。

此外，尽管近年来提出了一些快速采样方法（如 DDIM，Song 等人，2020）以降低扩散模型的采样成本，我们的 CL-DiffPhyCon 仍具有独立的加速效果。关键见解在于，CL-DiffPhyCon 并未减少每个物理时间步的采样总步数，因此可以结合 DDIM 方法。因此，通过在同步和异步模型的采样过程中引入 DDIM，CL-DiffPhyCon 的控制效率可以进一步提升。

5 实验

在实验中，我们旨在解答以下三个问题：（1）CL-DiffPhyCon 能否在性能上超越经典方法和最新的基线方法？（2）与扩散控制基线相比，CL-DiffPhyCon 能否如第 4.4 节分析的那样实现推理加速，并通过引入 DDIM（Song 等人，2020）获得进一步的加速？（3）CL-DiffPhyCon 能否应对部分观测和高维间接控制的挑战？为了解答这些问题，我们在两个控制任务上开展了实验：一维 Burgers 方程控制和二维不可压缩流体控制。这两个任务在科学和工程领域中具有重要的实际应用。此外，我们还提供了复现我们方法及基线方法结果的代码。

5.1 基准线

以下经典、模仿学习、强化学习和扩散控制方法被选为基准：经典控制算法，广泛应用于各个领域的比例-积分-微分（PID）控制（Li 等，2006）；一种模仿学习方法，行为克隆（BC）（Pomerleau，1988）；一种最近的强化学习方法，行为近端策略优化（BPPO）（Zhuang 等，2023）；两种扩散控制方法，包括 RDM（Zhou 等，2024），该方法自适应地决定何时进行完整的控制序列采样过程；DiffPhyCon（Wei 等，2024），其原始版本在一个采样过程中规划整个轨迹的控制序列。由于轨迹长度远大于扩散模型的时间范围 $H$（1D 和 2D 任务中分别设置为 $H=16$ 和 $H=15$），我们将 DiffPhyCon 扩展为其三个变体 DiffPhyCon-h（ h ∈ { 1 , 5 , H − 1 } h \in \{1, 5, H-1\} h∈{1,5,H−1}）作为基准，基于预定义的物理时间步长间隔 $h$ 来进行完整的采样过程。所有扩散控制基准都使用 CL-DiffPhyCon 训练的同步扩散模型进行采样，1D 和 2D 任务的步骤分别为 $T=900$ 和 $T=600$。PID 不适用于复杂的 2D 任务（Aström & Hägglund，2000）。RDM 是按照官方代码重新实现的。然而，RDM 中两个阈值超参数的默认值在我们的任务中表现不佳。因此，我们选择了一对表现最佳的值（详细信息见附录 H）。对于其他基准，我们遵循 DiffPhyCon Wei 等（2024）中的实现。

5.2 1D 伯格斯方程控制

实验设置：伯格斯方程是广泛用于描述各种物理系统的方程。我们遵循 Hwang 等（2022）；Mowlavi & Nabi（2023）的工作，考虑带有狄利克雷边界条件和外力 $\mathbf{w}(t,x)$ 的 1D 伯格斯方程，具体形式如下：

{ ∂ u ∂ t = − u ⋅ ∂ u ∂ x + ν ∂ 2 u ∂ x 2 + w ( t , x ) in  [ 0 , T ] × Ω , u ( t , x ) = 0 in  [ 0 , T ] × ∂ Ω , u ( 0 , x ) = u 0 ( x ) in  { t = 0 } × Ω , (15) \begin{cases} \frac{\partial u}{\partial t} = -u \cdot \frac{\partial u}{\partial x} + \nu \frac{\partial^2 u}{\partial x^2} + w(t, x) & \quad \text{in} \ [0, T] \times \Omega, \\ u(t, x) = 0 & \quad \text{in} \ [0, T] \times \partial \Omega, \\ u(0, x) = u_0(x) & \quad \text{in} \ \{t = 0\} \times \Omega,\\ \end{cases} \tag{15} ⎩ ⎨ ⎧​∂t∂u​=−u⋅∂x∂u​+ν∂x2∂2u​+w(t,x)u(t,x)=0u(0,x)=u0​(x)​in [0,T]×Ω,in [0,T]×∂Ω,in {t=0}×Ω,​(15)

其中，$\nu$ 是粘度，$u_0(x)$ 是初始条件。给定在 $[0,T]\times \Omega$ 范围内定义的目标状态 $u_d(t,x)$，控制目标 $J$ 是最小化 $u(t,x)$ 和 $u_d(t,x)$ 之间的误差：

$$\mathcal{J}:={\int }_T{\int }_{\Omega }|u(t,x)-u_d(t,x){|}^{2}dxdt,\text{\hspace{1em}(16)}$$

满足方程（15）。为了更好地适应实际场景的需求，我们遵循 Wei 等（2024）的做法，选择两种不同的实验设置：完全观察（FO）和部分观察（PO）。完全观察（FO）意味着所有空间位置都可以观察和控制，而部分观察（PO）意味着只有部分位置可以观察和控制，如图 3 所示。有关数据集、实验设置和实现的详细信息，请参见附录 D。

表1：1D Burgers 方程控制的比较结果。报告了在单个 NVIDIA A100 80GB GPU（16 个 GPU 核心）上的平均控制目标 $\mathcal{J}$ 和推理时间。最佳模型使用粗体表示，亚军使用下划线标出。

结果。在表 1 中，我们展示了我们提出的 CL-DiffPhyCon 和基准方法的结果。需要注意的是，在两种不同的设置 FO/PO 中报告的指标是不可直接比较的。可以看出，在 FO 和 PO 设置中，CL-DiffPhyCon 都比所有基准方法表现得更好。BC 和 BPPO 在此任务中表现较差，因为它们在很大程度上依赖于训练控制序列的质量，而这些序列远非最佳解（见附录 D.1）。与 BC 和 BPPO 相比，扩散控制基准方法表现更好，因为扩散模型在控制目标的指导下对每个时间范围进行全局优化，并且擅长从高维空间中采样（Wei 等，2024）。尽管如此，CL-DiffPhyCon 相比于 DiffPhyCon-1 和 DiffPhyCon-5（分别是 FO 和 PO 设置中的最佳基准），分别降低了 54.3% 和 51.1% 的 $\mathcal{J}$ 值。在图 3 中，我们展示了在 FO 和 PO 设置下随机选取的一个测试样本的可视化比较。结果显示，CL-DiffPhyCon 相比于最佳基准在目标状态的累计误差上更低。CL-DiffPhyCon 相较于这些扩散控制基准的优越性表明，闭环控制对于良好的控制性能是必要的，且我们的 CL-DiffPhyCon 能够有效地将来自环境的实时反馈整合到随后的控制信号采样中。更多的可视化结果见附录 E。除了控制性能的提升，表 1 还显示，CL-DiffPhyCon 显著加速了采样过程，相较于 DiffPhyCon-h 快了约 $H/h$ 倍，比 RDM 快了两倍。此外，通过与 30 步采样的 DDIM 结合，CL-DiffPhyCon 实现了额外的 5 倍加速，展示了 CL-DiffPhyCon 相比现有的快速采样方法在扩散模型中的独立加速效果。因此，CL-DiffPhyCon 的时间成本降低，已可与非扩散控制基准 BC、BPPO 和 PID 相媲美，而其性能则远超这些方法。

6 结论

在本文中，我们提出了 CL-DiffPhyCon，这是一种基于扩散的复杂物理系统闭环控制方法，基于目标分布的理论分析。通过在模型的时间范围内逐渐增加噪声水平，CL-DiffPhyCon 使得与环境的实时交互成为可能，并实现了控制信号的高效采样。在两个物理控制任务上的实验展示了其优越的控制性能和效率。我们相信，我们的方法对于物理系统之外的其他任务（如机器人和无人机控制）也具有潜在的应用价值。有关我们方法的局限性和未来工作的更多讨论，请参见附录 I。

D 1D Burgers 方程控制的详细信息

D.1 数据集

我们遵循 Wei 等（2024）中的说明生成 1D Burgers 方程数据集。具体来说，我们使用有限差分法（FDM）在空间范围 $x\in [0,1]$ 和时间范围 $t\in [0,1]$ 内生成轨迹，初始状态和控制序列遵循一定的分布。空间被离散化为 128 个网格，时间被离散化为 10000 个步骤。我们生成了 90000 个轨迹用于训练集，50 个轨迹用于测试集。对于每个训练样本，其目标状态 ${\mathbf{u}}_d(t,x)$ 是从其他训练样本中随机选择的，这意味着几乎所有训练集中的样本都是不成功的，因为它们几乎无法在随机控制下达到目标状态。因此，生成超越训练数据集性能的控制序列是一个挑战。

D.2 实验设置

与Wei 等（2024）类似，我们在第 5.2 节设计了 1D Burgers 方程控制的两种设置：

• 完全观测（Full Observation）：这是一个常见场景，系统的所有状态 $\mathbf{u}(t,x)$，其中 $x\in [0,1]$ 且 $t\in [0,1]$ 都可以被观测。
• 部分观测（Partial Observation）：在此设置下，我们隐藏了部分 $\mathbf{u}$ 的数据，控制模型的性能通过图 3 中显示的带斜线的白色区域来衡量。具体而言，训练集中的数据在区间 $x\in [\frac{1}{4},\frac{3}{4}]$ 上的状态 $\mathbf{u}(t,x)$ 被设置为零，测试时也设置 ${\mathbf{u}}_0(x)$，其中 $x\in [\frac{1}{4},\frac{3}{4}]$ 为零。仅观察、控制和评估区域 $\Omega =[0,\frac{1}{4}]\cup [\frac{3}{4},1]$ 的数据。这个设置特别具有挑战性，因为部分观测引入了不确定性。

图 3：CL-DiffPhyCon（底行）与最佳基准模型（中行）在 1D Burgers 方程控制上的视觉比较。图中，J是在可视化样本上计算的。在 PO 设置（右图）中，带斜线的白色区域为不可观察区域。

D.3 实现

我们使用 U-Net（Ronneberger 等，2015）作为模型 ${ϵ}_{\varphi }$ 和 ${ϵ}_{\theta }$ 的架构。这两个模型分别使用相同的训练数据集进行训练。需要注意的是，在部分观测设置中，未观测的数据在训练和测试过程中对模型不可见，如附录 D.2 中所述。我们在模型输入和条件的对应位置填充零，并在训练损失中排除这些位置。因此，模型仅学习观察到的状态与控制序列之间的相关性。我们使用均方误差（MSE）损失来训练模型。两个模型都有 900 个扩散步骤（T = 900）。我们使用的 DDIM（Song 等，2020）采样只有 30 个扩散步骤，超参数 $\eta =1$。模型和训练的超参数列在表 3 中。

E 1D 可视化结果

我们展示了我们的方法和基准模型在两种设置下（FO 和 PO）的可视化结果，分别在图 5、6、7 和 8 中展示。在每个设置下，我们展示了从测试集中随机选择的四个样本的结果。从这些可视化结果中可以观察到，在我们提出的 CL-DiffPhyCon 控制下，状态 $\mathbf{u}(t,x)$ 更接近目标状态 ${\mathbf{u}}_d(t,x)$，无论初始状态如何。此外，这一观察在 FO 和 PO 两种设置下都一致，表明我们的 CL-DiffPhyCon 在解决部分观测挑战方面是有效的。

图 5：1D Burgers 方程控制在 FO（完全观测）设置下的两种可视化结果。第一行是目标 ${\mathbf{u}}_d(t,x)$，误差 $\mathbf{u}(t,x)-{\mathbf{u}}_d(t,x)$ 衡量了受控状态与目标之间的差距。横轴是时间坐标，纵轴是状态。

图 7：1D Burgers 方程控制在 PO（部分观测）设置下的两种可视化结果。第一行是目标 ${\mathbf{u}}_d(t,x)$，误差 $\mathbf{u}(t,x)-{\mathbf{u}}_d(t,x)$ 衡量了受控状态与目标之间的差距。横轴是时间坐标，纵轴是状态。中间带斜线的白色区域表示不可观察区域。