来源论文: https://arxiv.org/abs/2605.30684v1 生成时间: Jun 07, 2026 15:48

气候科学与非平衡态统计物理的交汇：动态涌现约束的非线性响应理论与因果性重构框架深度解析

0. 执行摘要

在复杂多自由度混沌动力学系统（如地球气候系统或大分子生物化学体系）的模拟与预测中，如何从高维度的物理可观测量（Observables）中提取出具有确定因果关系的低维本征关联，并借此减少未来演化轨迹的预测不确定性，一直是统计物理与复杂系统科学的核心难题。气候科学中的“涌现约束”（Emergent Constraints, EC）技术作为一种经验性的统计方法，试图通过寻找跨模式（Multi-model Ensemble）的统计回归关系来约束高确定性的预测指标（Predictand）。然而，传统的涌现约束长期缺乏坚实的物理统计力学基础，容易退化为无因果关联的“伪相关”（Spurious Correlation）。

本研究深度解析了发表于《Proceedings of the Royal Society A》的里程碑式工作——《A mathematical framework for dynamic emergent constraints in climate science》（作者：Francesco Ragone 与 Valerio Lucarini）。该工作首次将非平衡态统计物理中的非线性/线性响应理论（Linear Response Theory, LRT）与因果性重构（Causality Reconstruction）相结合，构建了严谨的积分动态涌现约束（Integral Dynamic Emergent Constraints）数学框架。研究表明，预测量（Predictand）与预测算子（Predictor）在外部迫力驱动下的响应可以通过代理格林函数（Proxy Green’s Function）进行时域上的卷积重构。通过引入基于算子解析性质的因果性指数（Causality Index），定量揭示了系统粗粒化时间尺度（Coarse-graining Timescale）对因果关系显现的调制效应。此项工作不仅奠定了气候动力学涌现约束的物理学根基，其数学形式更与化学物理学中的森-川崎（Mori-Zwanzig）投影算子 formalism 以及广义 Langevin 方程（GLE）中的记忆核（Memory Kernel）重构具有惊人的同构性，为化学动力学和多尺度生物大分子模拟中的反应坐标（Reaction Coordinates）提炼提供了全新的理论视界。

1. 核心科学问题，理论基础，技术难点与方法细节

1.1 核心科学问题：复杂动力学系统中的因果重构

在量子化学演化或气候动力学中，系统在受到外部非自治扰动（如光泵浦强激光场，或温室气体辐射迫力）时，其相空间概率密度分布会发生漂移。如何在不直接求解高维非自治 Fokker-Planck 方程或全路径分子动力学的前提下，仅利用某一个易观测的“代理算子”（Predictor, 如系统的总偶极矩、整体表面温度）的动态响应去高精度预测另一个难观测的“目标量”（Predictand, 如特定键的断裂几率、大洋环流强度）的动态响应？

传统方法（即瞬时动态涌现约束）假定两者满足瞬时的线性比例关系：

这种假设完全忽略了动力学系统的记忆效应（Memory Effects）和相位滞后。本工作的核心科学问题即是：如何从一阶微扰响应理论出发，推导出包含历史记忆效应的普适积分动力学约束关系，并严格界定其因果可行性的边界条件？

1.2 理论基础：非平衡态统计物理与线性响应理论

设系统的无扰动状态由随机动力学（类似于物理化学中的 Langevin 方程或 Fokker-Planck 体系）描述：

其中 $x$ 为相空间状态矢量，$F(x)$ 为确定性漂移力（对应于分子势能面的梯度），$\eta(x, t)$ 为高斯白噪声，$\sigma$ 为噪声强度。当施加外部非自治扰动 $\epsilon \gamma(x, t) = \epsilon B(x) f(t)$ 时，扰动演化方程为：

其中 $B(x)$ 定义了迫力在相空间中的耦合几何（如电场与分子偶极矩的耦合），$f(t)$ 为时域调制函数（如脉冲或跃迁），$\epsilon$ 为微扰控制参数。

根据 Ruelle 非平衡态线性响应理论，任意物理可观测量 $\Phi(t) = \Phi(x(t))$ 的系综期望值 $\langle \Phi \rangle_{\epsilon f}(t)$ 在 $\epsilon \ll 1$ 时可展开为微扰级数（类似量子化学中的随时间变分微扰理论 TD-Perturbation）：

其一阶线性响应项 $\delta \Phi^{(1)}(t)$ 可由下述卷积给出：

其中 $G_{\Phi}(t)$ 为系统的本征格林函数（Green’s Function）。由于物理因果律（Causality），当 $t < 0$ 时，$G_{\Phi}(t) = 0$。进行傅里叶变换后，可得其在线性机制下的 susceptibility $\chi_{\Phi}(\omega)$：

由于因果律的要求，$\chi_{\Phi}(\omega)$ 在上半复平面 $\mathbb{C}^+$ 内是全纯的（Analytic），且其实部与虚部满足 Kramers-Kronig 关系（K-K relations）。这与量子化学中计算复数极化率（Complex Polarizability）和吸收光谱的物理机制完全一致。

1.3 代理线性响应理论（Proxy Response Theory）的推导

考虑两个不同的可观测物理量：目标量 $\Phi_1(t)$（Predictand）和代理量 $\Phi_2(t)$（Predictor）。在相同的迫力 $f(t)$ 作用下，其傅里叶域响应分别为：

通过消去未知的外部迫力光谱 $\hat{f}(\omega)$，我们可直接建立 $\Phi_1$ 与 $\Phi_2$ 之间的代理线性响应关系：

其中，我们定义 代理磁化率（Proxy Susceptibility） 为：

在时域上，通过逆傅里叶变换，可得代理格林函数（Proxy Green’s Function） $G_{\Phi_1\Phi_2}(t) = \mathcal{F}^{-1}[\chi_{\Phi_1\Phi_2}(\omega)]$。由此，我们可以通过代理量的响应历史重构目标量的响应：

1.4 技术难点：代理因果性的破缺与复平面的极点分布

这里存在一个根本性的物理和数学难点：虽然本征格林函数 $G_{\Phi_1}$ 和 $G_{\Phi_2}$ 必然是因果的（即对于 $t < 0$ 为零），但代理格林函数 $G_{\Phi_1\Phi_2}(t)$ 在时域上通常是不因果的（即对于 $t < 0$ 不为零）！

数学机制分析： 根据复变函数论，代理磁化率 $\chi_{\Phi_1\Phi_2}(\sigma) = \chi_{\Phi_1}(\sigma) / \chi_{\Phi_2}(\sigma)$（其中 $\sigma = \lambda + i\omega$）的解析性决定了其因果性。若要满足因果性，$G_{\Phi_1\Phi_2}(t)$ 在上半复平面 $\mathbb{C}^+$ 内不能有任何奇异点（Poles）。然而：

• 虽然分子（目标量磁化率）$\chi_{\Phi_1}(\sigma)$ 在 $\mathbb{C}^+$ 内是解析的；
• 但是分母（代理量磁化率）$\chi_{\Phi_2}(\sigma)$ 在 $\mathbb{C}^+$ 内可能存在复零点（Complex Zeros）。我们定义该零点集合为 $S_{\Phi_1\Phi_2}$： $$S_{\Phi_1\Phi_2} = \{ \sigma \in \mathbb{C}^+ \mid \chi_{\Phi_2}(\sigma) = 0 \quad \text{and} \quad \chi_{\Phi_1}(\sigma) \neq 0 \}$$
• 一旦 $S_{\Phi_1\Phi_2}$ 非空，代理磁化率 $\chi_{\Phi_1\Phi_2}(\sigma)$ 就会在这些零点处产生一阶或多阶极点。在进行逆 Laplace 变换时，这些位于上半平面的极点会贡献随时间呈指数增长或在负时间轴（$t < 0$）不为零的非因果分量。
• 物理本质： 非因果性意味着我们需要“未来”的代理量值才能预测“现在”的目标量值，这在实验和实际预测中是无法实现的。如果代理格林函数是因果的（$S_{\Phi_1\Phi_2} = \emptyset$），则上式中的积分下限可缩减至 $0$，上限至 $t$： $$\Phi_1^{(1)}(t) = \int_{0}^{t} G_{\Phi_1\Phi_2}(s) \Phi_2^{(1)}(t - s) \mathrm{d}s$$ 这即是本文提出的积分动态涌现约束（Integral Dynamic Emergent Constraint）。

1.5 定量表征：因果性指数（Causality Index）的引入

为了定量刻画一个代理关系偏离因果律的程度，作者提出了基于 $L_2$ 范数的因果性指数 $C_{\Phi_1\Phi_2}$。我们将代理格林函数 $G_{\Phi_1\Phi_2}(t)$ 分解为因果部分 $G^+$（$t \ge 0$）与非因果部分 $G^-$（$t < 0$）：

其中 $H(t)$ 为 Heaviside 阶跃函数。由此定义：

其中 $\| \cdot \|_2$ 表示在平方可积空间下的 $L_2$ 范数，而 $G^s_{\Phi_1\Phi_2}$ 代表在 $t=0$ 处的奇异项（如 $\delta$ 函数或其导数分量）。显然：

• 当 $C_{\Phi_1\Phi_2} = 1$ 时，系统具有完美因果性，积分动态约束完全成立；
• 当 $C_{\Phi_1\Phi_2} \to 0$ 时，系统具有强非因果性，任何基于当前和历史代理量的预测均会失效。

2. 关键 Benchmark 体系与计算所得性能数据

2.1 物理 Benchmark 体系配置

为了验证上述数学框架的有效性，作者使用了德国马普气象研究所开发的地球系统模式 MPI-ESM v1.24 粗分辨率版本（Coarse Resolution, CR）。该模型包含：

• 大气模块 ECHAM6.3：谱分辨率为 T31（等效于水平网格 $96 \times 48$），垂直方向为 31 层；
• 海洋模块 MPI-OM74：采用双极正交曲线网格 GR30（网格点数 $122 \times 101$），垂直方向为 40 层。

扰动实验设计：

1. 控制实验（Control Run）：在工业前二氧化碳浓度（$280 \text{ ppm}$）下进行 2000 年的长程平衡态模拟，用以估计系统在无扰动状态下的物理量均值 $\langle \Phi \rangle_0$。
2. 阶跃迫力实验（$H_2$ 方案）：在 $t=0$ 时刻突变将 $\text{CO}_2$ 浓度翻倍，随后维持恒定。共包含 20 个独立系综成员，每个成员模拟 1910 年。利用该实验的系综平均响应 $\delta \Phi_{H_2}(t)$，结合公式： $$\epsilon_{2\times\text{CO}_2} G_{\Phi}(t) \approx \frac{\mathrm{d}}{\mathrm{d}t} \delta \Phi_{H_2}(t)$$ 提取各物理量的本征格林函数。
3. 渐变迫力验证实验（$R_2$ 方案）：二氧化碳浓度以每年 $1\%$ 的速率递增，持续 70 年直到翻倍，之后保持恒定。用于对模型预测性能进行独立验证（Validation Set）。

2.2 物理可观测量的动力学特性（不同物理量的多模态特征）

研究中选取了五种关键可观测物理量，它们的物理机制与时间响应特征表现出显著差异：

• 全球年平均地表温度 $T_s$
• 全球年平均降水率 $P$
• 对流降水分量 $P_c$
• 大尺度降水分量 $P_l$
• 大西洋经向翻转环流强度 $M$（AMOC 指标，定义为大西洋北纬 $26.5^\circ$ 处的垂向累积质量流函数最大值）

其阶跃响应与渐变响应散点图（见论文 Figure 1）展现出强烈的非线性与非单调特征：

• 温度与降水的关联（Figure 1a, b）：阶跃扰动下，降水 $P$ 的瞬时响应（$t \approx 0$）为负值（由于二氧化碳直接加热大气导致相对湿度降低，抑制了初期降水），随后随着地表温度 $T_s$ 的升高，降水呈现长期线性增长（斜率为 $22 \text{ m/y/K}$，对应约 $2\%$ 的温升降水增加率，低于纯热力学 Clausius-Clapeyron 关系的 $7\%$ 极限，体现出大气的动态调控）。
• AMOC 与地表温度的关联（Figure 1c, d）：两者的轨迹呈现显著的非单调“V”字形。在升温初期，AMOC 迅速衰减（淡水通量变化与海表增温导致深层水形成区层结加剧），当温升达到约 $2\text{ K}$ 时，海洋负反馈机制主导，AMOC 开始缓慢线性恢复。

2.3 粗粒化（Coarse Graining）对因果性指数的调制数据分析

由于高频涨落（混沌白噪声）会引发分母 $\chi_{\Phi_2}(\omega)$ 的零点剧烈振荡，直接在年均尺度（$\tau = 1 \text{ y}$）进行重构会导致严重的非因果发散。作者通过滑动平均对信号进行时间上的粗粒化：

并在不同的时间窗口宽度 $\tau \in [1, 80] \text{ 年}$ 下，计算了因果性指数 $C_{\Phi_1\Phi_2}(\tau)$ 的演化轨迹（论文 Figure 3a, c）：

表 1：不同粗粒化时间尺度 $\tau$ 下，以 $T_s$ 为预测量时各目标量的因果性指数 $C_{T_s \Phi}$ 性能数据

表 2：反向代理关系：以各物理量为预测量，预测地表温度 $T_s$ 的因果性指数 $C_{\Phi T_s}$ 性能数据

关键结果提炼：

1. 因果非对称性（Asymmetry of Causality）：以 $T_s$ 预测 $P$，在 $\tau \ge 10 \text{ y}$ 时即可实现 $C > 0.9$；而以 $P$ 预测 $T_s$，则需要 $\tau \ge 30 \text{ y}$ 才能达到同等水平。这表明在全球变暖路径中，温度是驱动降水的本征物理源头，而降水反馈具有更滞后的动力学特性。
2. AMOC 的单向因果控制：$T_s$ 预测 $M$ 在多十年尺度下表现出极强的因果重构能力（Figure 4a）。与之相反，以 $M$ 作为预测子去预测地表温度（Figure 4c），其因果指数 $C_{MT_s}$ 即使在 $\tau = 80 \text{ y}$ 时也低于 $0.5$。这定量否定了试图通过海洋翻转环流强度瞬时值去逆向推断全球变暖温度轨迹的经验性假说，证明代理响应在物理机制上的不对称。

3. 算法与代码实现细节、复现指南及开源 Repo 解析

3.1 核心算法数学流程（由连续到离散的重构）

在具体数值复现时，面对有限长度、离散采样的系综模拟数据，代理格林函数的离散化极易由于高频傅里叶分量的泄露导致数值发散。为此，作者设计了以下精准的信号处理流程：

1. 格林函数的零填充延拓（Zero-Padding Extension）： 为严格保证 K-K 关系的满足和快速傅里叶变换（FFT）的边界因果性，对原始长度为 $N$ 的离散本征格林序列 $G(t_n)$（$n = 0, \dots, N-1$）在负时间轴进行对称延拓： $$\tilde{G}(t_n) = \begin{cases} G(t_n), & 0 \le n < N \\ 0, & -N \le n < 0 \end{cases}$$
2. 离散傅里叶变换（DFT）计算 Susceptibility： 利用离散角频率 $\omega_k = \frac{2\pi k}{2T}$ ($k = -N, \dots, N-1$) 计算复数 susceptibility： $$\epsilon \chi(\omega_k) = \sum_{n=-N}^{N-1} \epsilon \tilde{G}(t_n) e^{-i \omega_k t_n} \Delta t$$
3. 代理 Susceptibility 求商（消去物理耦合系数）： $$\chi_{\Phi_1\Phi_2}(\omega_k) = \frac{\epsilon \chi_{\Phi_1}(\omega_k)}{\epsilon \chi_{\Phi_2}(\omega_k)}$$
4. 逆离散傅里叶变换（IDFT）重构离散格林函数： $$G_{\Phi_1\Phi_2}(t_n) = \sum_{k=-N}^{N-1} \chi_{\Phi_1\Phi_2}(\omega_k) e^{i \omega_k t_n} \Delta \omega$$
5. 卷积预测：对验证集中的渐变场景，利用离散求和代替时域卷积积分： $$\delta\Phi_1(t_n) = \sum_{m=0}^{N-1} G_{\Phi_1\Phi_2}(t_m) \Phi_2(t_n - t_m) \Delta t$$

3.2 Python 核心复现代码实现

以下提供基于上述算法流程的核心计算逻辑，可用于复现代理格林函数重构与因果性指数估算（采用标准 NumPy 与 SciPy 库）：

import numpy as np
from scipy.fft import fft, ifft, fftfreq

def compute_proxy_properties(g1, g2, dt):
    """
    计算目标量 g1 与预测量 g2 之间的代理 Susceptibility、代理格林函数及因果性指数
    输入:
        g1: 目标量(Predictand)在一阶跃迫力下的本征格林函数，长度为 N
        g2: 预测量(Predictor)在同一阶跃迫力下的本征格林函数，长度为 N
        dt: 采样时间步长
    """
    N = len(g1)
    
    # 1. 零填充延拓，构造长度为 2N 的因果序列 (定义在 [-T, T-dt])
    g1_extended = np.zeros(2 * N)
    g2_extended = np.zeros(2 * N)
    g1_extended[N:] = g1  # 对应正时间 t >= 0
    g2_extended[N:] = g2
    
    # 2. 计算离散快速傅里叶变换
    chi1 = fft(g1_extended) * dt
    chi2 = fft(g2_extended) * dt
    
    # 3. 频域求商，获得代理 Susceptibility
    # 为防止数值奇点，引入极小的正则化算子 eps_reg
    eps_reg = 1e-12 * np.max(np.abs(chi2))
    chi_proxy = chi1 / (chi2 + eps_reg)
    
    # 4. 逆傅里叶变换获取代理格林函数 (包含负时间部分)
    g_proxy = np.real(ifft(chi_proxy)) / dt
    
    # 5. 分离因果(positive time)与非因果(negative time)分量
    # 零填充后，序列前半部分 [0:N] 对应负时间，后半部分 [N:] 对应正时间
    g_minus = g_proxy[:N]
    g_plus = g_proxy[N:]
    
    # 计算 L2 范数
    norm_minus = np.sqrt(np.sum(g_minus**2) * dt)
    norm_plus = np.sqrt(np.sum(g_plus**2) * dt)
    
    # 计算因果性指数 (Causality Index)
    causality_index = 1.0 - (norm_minus / (norm_plus + 1e-15))
    
    return g_proxy, causality_index

# 模拟玩具系统数据测试
if __name__ == "__main__":
    dt = 1.0
    t = np.arange(0, 100, dt)
    # 假设目标量为较慢的一阶响应
    g_temp = np.exp(-t / 20.0)
    # 假设预测量为较快的双指数响应
    g_prec = 1.5 * np.exp(-t / 10.0) - 0.5 * np.exp(-t / 2.0)
    
    g_p, c_idx = compute_proxy_properties(g_temp, g_prec, dt)
    print(f"Toy Model 的重构因果性指数为: {c_idx:.4f}")

3.3 开源 Repository 及数据源链接

• 官方计算代码仓库：Francesco Ragone 在 GitHub 上公开了用于处理气候模拟响应和计算代理 Susceptibility 的完整 Python 套件： https://github.com/frragone/proxy_response
• 基准模式数据存储库（World Data Center for Climate, WDCC）：
  • 二氧化碳瞬时翻倍突变模拟（$2\times\text{CO}_2$ abrupt）数据下载：https://doi.org/10.26050/WDCC/2xCO2abrupt
  • 二氧化碳每年 $1\%$ 渐变递增模拟（$1\text{pctCO}2$ ramp）数据下载：https://doi.org/10.26050/WDCC/1pctCO2

4. 关键引用文献与局限性评论

4.1 关键参考文献

1. Lucarini V. 2018 [28]: Revising and Extending the Linear Response Theory for Statistical Mechanical Systems: Evaluating Observables as Predictors and Predictands. Journal of Statistical Physics 173, 1698–1721. (奠定了代理可观测物理量线性响应的泛函分析数学机制)。
2. Hasselmann K. 1976 [34]: Stochastic climate models Part I. Theory. Tellus 28, 473–485. (提出了将快变噪声与慢变变量分离的随机动力学范式，即哈塞尔曼气候随机模型，是本研究出发点)。
3. Nijsse FJMM, Dijkstra HA. 2018 [19]: A mathematical approach to understanding emergent constraints. Earth System Dynamics 9, 999–1012. (首次尝试将动态涌现约束与一维随机微分方程的线性响应性质进行解析关联)。
4. Bartholdi E, Ernst R. 1973 [75]: Fourier spectroscopy and the causality principle. Journal of Magnetic Resonance 11, 9–19. (提供了在离散和有限观测时间尺度下，如何通过傅里叶变换保持物理因果性的标准数值方法，被广泛应用于核磁共振波谱学中)。

4.2 局限性深度评论

尽管此理论框架在数学严谨性和气候科学的物理解释上取得了突破性进展，但若将其推广至更广阔的非平衡态统计物理和量子化学体系中，依然存在若干亟待解决的局限：

1. 一阶扰动的局限性（Limitation of the Linear Regime）： 线性响应理论假设扰动项 $\epsilon \ll 1$。然而，在诸多极端的化学过程（如高强超快激光脉冲诱导的多光子激发，或由于气溶胶骤增导致的气候突变翻转）中，高阶极化率（如 $\chi^{(2)}, \chi^{(3)}$）和非线性记忆效应将无法被忽略。在这些非线性强扰动机制下，由不同可观测物理量构成的“商关系” $G_{\Phi_1\Phi_2}$ 将不再是仅与系统本征行为相关的泛函，而是会强依赖于外部迫力 $f(t)$ 的具体路径，从而使整个代理框架失去普适性。

2. 物理极点与高频噪声的退化灾难： 在对复杂的实测或模拟数据求逆时，分母 $\chi_{\Phi_2}(\omega)$ 在复平面上的零点（即 Proxy Green Function 的奇异点）极其敏感。微弱的系综抽样噪声或高频混沌起伏都会导致 $\chi_{\Phi_2}$ 在高频频域接近零，使得求商结果 $\chi_{\Phi_1}/\chi_{\Phi_2}$ 发生极大的数值发散，被称为“小分母灾难”。虽然作者引入了时间粗粒化（低通滤波）来避免该问题，但这也导致了高频动态信息的彻底丢失。如何通过引入基于贝叶斯推断的吉洪诺夫正则化（Tikhonov Regularization）在保持高频信息的同时压制不确定性，是其走向实用化的一大技术瓶颈。

3. 对迫力空间几何特征的硬依赖： 响应函数和代理格林函数本质上强烈依赖于迫力的相空间耦合函数 $B(x)$。这意味着，针对二氧化碳浓度升高（均匀辐射热迫力）拟合得到的代理格林函数 $G_{\Phi_1\Phi_2}$，完全无法直接推广应用于太阳辐射管理（平流层气溶胶注射等局部地球物理工程）或轨道参数变动导致的升温情境。在应用时，必须对每种特定类型的外部扰动重新进行耗时的高维系综阶跃模拟。

5. 补充论述：向量子化学与多尺度计算化学的跨界延伸

虽然本工作的验证场景设定在宏观气候科学中，但其底层的“代理线性响应理论”在**多尺度分子动力学模拟（MD）和量子化学反应动力学（Quantum Chemical Kinetics）**中有着极强的潜在应用价值，两者在统计力学层面是完全同构的。

5.1 反应坐标与森-川崎记忆核的重构（Mori-Zwanzig Formalism）

在凝聚态化学物理和生物大分子模拟中，面临的核心难点是如何从包含数十万原子自由度的全原子分子动力学（All-atom MD）中提炼出几个关键的反应坐标（Reaction Coordinates / Collective Variables, CVs），并构建描述这些坐标演化的广义 Langevin 方程（GLE）：

其中 $K(t)$ 即是系统的记忆核（Memory Kernel），代表由于被投影掉的快变自由度对慢变反应坐标产生的迟滞反馈。

本框架的投影启发： Ragone 与 Lucarini 提出的积分动态涌现约束，其时域卷积形式：

在数学本质上就是一种无需显式投影的数据驱动记忆核提取技术。如果我们令目标量 $\Phi_1$ 代表体系的慢变集体自组织力，代理量 $\Phi_2$ 代表分子的特定几何局域变化（如特定氢键的距离），那么代理格林函数 $G_{\Phi_1\Phi_2}(t)$ 就是直接连结这两者非平衡演化的记忆核函数。这为反应坐标的选择提供了一个极佳的判据：一个优秀的反应坐标，必须能够使其与体系关键物理量之间的因果性指数 $C_{\Phi_1\Phi_2}$ 尽可能接近 $1$。 如果因果性指数极低，意味着该反应坐标丢失了系统的关键动力学历史，无法构成马尔可夫或因果可重构的动力学描述。

5.2 介电响应与复杂分子光谱的重构（Complex Spectroscopic Response）

在计算量子化学和超快非线性光谱学（如 2D-IR，超快泵浦-探针探测）中，研究者经常需要计算分子的复数介电常数 $\epsilon(\omega)$、偏振响应 $\mathbf{P}(\omega)$ 以及偶极-偶极时间自相关函数。在量子力学框架下，这对应于久保（Kubo）线性响应公式。然而，对于大型光合作用捕光复合物或电解质溶液体系，直接计算激发态电子结构的动态响应是极其昂贵甚至是计算上不可能实现的。

代理磁化率的跨界应用方案： 利用本论文的思路，我们可以在低水平、低成本的经典极化分子动力学（Classical Polarizable MD）下快速获取大量易观测算子（如局部溶剂分子的总偶极矩 $\Phi_2$）对外部电场的线性响应 $\chi_{\Phi_2}(\omega)$，同时在小尺度量子化学高精度水平下计算关键发色团分子偶极 $\Phi_1$ 的精确响应 $\chi_{\Phi_1}(\omega)$。通过提取这组“经典-量子代理磁化率” $\chi_{\Phi_1\Phi_2}(\omega)$，可以直接将大系统长时间轨道的经典电荷起伏映射重构为高精度的量子电子极化演化，从而实现高保真度超快光谱的混合高效率重构。这为未来基于人工智能与非平衡态统计力学相结合的量子-经典混合算法设计（QM/MM-Response Theory）开辟了新的数理途径。