来源论文: https://arxiv.org/abs/2603.22710v1 生成时间: Mar 25, 2026 09:53

0. 执行摘要

量子波导电动力学（QED）中的巨腔系统因其独特的远程耦合机制和伴随的传播延迟，表现出显著的非马尔可夫动力学特性，这与传统的量子光学腔系统截然不同。为了在这种复杂的非马尔可夫背景下实现反馈控制，精确估计量子态的演化至关重要。本文深入研究了一种为巨腔系统量身定制的最优滤波器设计框架，旨在通过连续量子测量，有效地追踪系统的量子态。该滤波器通过迭代计算延迟态协方差矩阵，巧妙地解决了由延迟引起的系统算符非对易性问题，并通过数值模拟验证了其对相干态和薛定谔猫态演化的卓越跟踪性能和有效性。

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

在量子信息处理和量子计量学等前沿领域，量子过滤与控制是提升量子系统性能的核心支柱。量子过滤的根本目标是从有限的测量信息中，最优地推断量子系统的状态或动力学特性。然而，当系统与环境之间的相互作用具有记忆效应，导致系统演化不仅依赖于当前状态，还依赖于其过去状态时，经典的量子过滤理论便无法直接适用，这类系统被称为非马尔可夫量子系统。

核心科学问题：

传统的量子过滤方法主要针对马尔可夫系统，即系统未来的状态仅取决于当前状态。然而，本文聚焦的波导量子电动力学（QED）巨腔系统，因其独特的物理构造——量子腔通过多个远程耦合点与波导场相互作用——必然导致非马尔可夫动力学。这些耦合点之间的有限传播距离引入了不可忽略的时间延迟，使得腔自身的发射场在一段时间后才能再次与腔相互作用，从而产生了内在的记忆效应。因此，为这类系统设计一个能够准确估计其量子态演化的最优滤波器，是实现有效反馈控制、量子计算和量子计量学的关键前提，同时也是一个尚未完全解决的重大挑战。

理论基础：

本文的理论基础深植于量子波导电动力学（WGQED）和开放量子系统理论。

1. 哈密顿量描述： 为了精确描述巨腔系统与波导场的相互作用，首先建立了总哈密顿量 $H = H_S + H_E + H_I$。

   • $H_S = \omega_c a^\dagger a$ 表示腔模式的内部哈密顿量，其中 $\omega_c$ 是腔的角频率，$a$ 和 $a^\dagger$ 分别是腔模式的湮灭和产生算符，满足正则对易关系 $[a, a^\dagger] = 1$。
   • $H_E = \int_{-\infty}^{+\infty} d\omega b^\dagger(\omega) b(\omega)$ 表示波导场的哈密顿量，其中 $b(\omega)$ 和 $b^\dagger(\omega)$ 是频率为 $\omega$ 的波导场模式的湮灭和产生算符。
   • $H_I = \sum_{n=1}^2 \frac{V_a}{\sqrt{v_g}} \int_{-\infty}^{+\infty} d\omega (b(\omega) e^{i\omega\tau_n} a^\dagger + H.c.)$ 表示腔模式与波导场在两个空间分离点（$x_1=0$ 和 $x_2=L$）的相互作用哈密顿量。这里的 $V_a$ 是腔-波导耦合强度，$v_g$ 是波导场的群速度，$\tau_n$ 表示与第 $n$ 个耦合点相关的传播延迟（$\tau_1=0, \tau_2=T=L/v_g$）。明确的传播延迟项 $e^{i\omega\tau_n}$ 是非马尔可夫动力学产生的直接原因。

2. 海森堡绘景中的朗之万方程： 从哈密顿量出发，推导了腔模式湮灭算符 $a(t)$ 的动力学方程，即朗之万方程。在海森堡绘景中，算符随时间演化，方程形式为： $ \dot{a}(t) = -i\omega_c a(t) - \frac{\gamma}{2} [a(t) + a(t-T)] - i\sqrt{\frac{\gamma}{2}} [b_{in}(t) + b_{in}(t-T)] $ (方程6) 其中 $\gamma = 4\pi V_a^2 / v_g$ 是衰减率，$T = \tau_2 - \tau_1 = L/v_g$ 是两个耦合点之间的传播延迟。这个方程清晰地揭示了非马尔可夫的本质：腔模式的演化不仅依赖于当前时刻的腔状态 $a(t)$ 和输入场 $b_{in}(t)$，还依赖于延迟时刻的腔状态 $a(t-T)$ 和输入场 $b_{in}(t-T)$。这与传统马尔可夫系统中只依赖当前时刻状态的朗之万方程形成了鲜明对比。

3. 输入-输出关系： 为了设计量子滤波器，还需要推导系统的输入-输出关系，即输出场 $b_{out}(t)$ 与输入场 $b_{in}(t)$ 和腔状态 $a(t)$ 之间的关系。 $ b_{out}(t) = b_{in}(t-T) - i\sqrt{\frac{\gamma}{2}} [a(t-T) + a(t)] $ (方程12) 这个方程同样包含了延迟项 $b_{in}(t-T)$ 和 $a(t-T)$，进一步印证了系统的非马尔可夫特性以及由波导介导的相干自反馈机制。这表明输出场也包含了过去时刻的信息。

4. 正交分量表示： 为了将复数值的湮灭/产生算符转换为实数值的系统状态变量，从而便于应用经典的卡尔曼滤波框架，将上述方程转换为正交分量表示。 定义腔模式的正交分量 $x(t) = [q(t) \quad p(t)]^T$，其中 $q(t)$ 和 $p(t)$ 分别是腔模式振幅的实部和虚部，对应于位置和动量算符。类似地，输入场 $b_{in}(t)$ 转换为量子白噪声 $dw(t) = [\text{Re}(b_{in}(t)) \quad \text{Im}(b_{in}(t))]^T$，输出场 $b_{out}(t)$ 转换为测量输出 $dy(t) = [\text{Re}(b_{out}(t)) \quad \text{Im}(b_{out}(t))]^T$。 转换后的动力学方程变为一个线性高斯系统（方程20）： $ dx(t) = Ax(t)dt + A_d x(t-T)dt + B dw(t) + B_d dw(t-T) $ $ dy(t) = Cx(t)dt + C_d x(t-T)dt + D_d dw(t-T) $ 其中 $A, A_d, B, B_d, C, C_d, D_d$ 是由原始系数矩阵通过 $M(\cdot)$ 运算（将复数算符转换为实数矩阵）得到的实系数矩阵。这种正交分量表示的优点在于，它将非对易的量子算符问题转化为对实值矩阵的运算，同时保留了量子动力学的关键特性，为滤波器设计奠定了基础。

技术难点：

本文所提出的系统面临的主要技术难点源于其非马尔可夫特性和量子效应：

1. 非马尔可夫动力学： 系统状态的演化依赖于其历史状态，使得标准的马尔可夫量子滤波理论无法直接应用。传统的处理方法如增广希尔伯特空间法，虽然可以处理具有有理谱的非马尔可夫环境，但会导致系统维度急剧增加，带来沉重的计算负担。对于本文研究的这种直接由延迟引起的状态依赖，增广系统法会引入高维的增广延迟项，使得计算成本极高。

2. 时间延迟引起的非对易性： 这是本文工作的核心挑战。在有时间延迟的量子系统中，不同时间点（例如 $t$ 和 $t-T$）的系统算符之间通常不再满足简单的对易关系（例如 $[a(t), a^\dagger(t-T)] \neq 0$）。这种非对易性阻止了现有量子滤波方法的直接应用，因为这些方法通常依赖于特定时间点的算符代数结构。如何在滤波器的设计中显式地、系统地处理这种跨时间点的非对易相关性是关键。

3. 混合延迟： 系统方程中不仅包含状态延迟 $x(t-T)$，还包含输入延迟 $dw(t-T)$，使得问题更加复杂，需要同时考虑多种延迟效应。

4. 量子测量反作用： 连续量子测量本身会对系统状态产生扰动，这被称为测量反作用（measurement back-action）。滤波器设计必须在估计系统状态的同时，明确地将这种不可避免的量子效应纳入考虑，确保估计结果的物理一致性，特别是当涉及到共轭变量的测量时，海森堡不确定性原理限制了它们同时被精确测量的可能性。

方法细节：

本文提出的最优滤波器设计框架巧妙地解决了上述技术难点，其核心思想在于引入并迭代计算“延迟误差协方差矩阵”：

1. 最优估计量： 滤波器的目标是构建系统状态 $x(t)$ 的最优估计量 $\hat{x}(t) = E[x(t) | \mathcal{F}_t]$，即在给定直到时刻 $t$ 的测量记录 $\mathcal{F}_t$ 的条件下，系统状态的条件期望。误差定义为 $e(t) = x(t) - \hat{x}(t)$。

2. 误差协方差矩阵： 为了描述估计误差的统计特性，定义了两种误差协方差矩阵：

   • 当前误差协方差矩阵 $P(t) = E[e(t)e^T(t) | \mathcal{F}_t]$：描述当前时刻估计误差的方差。
   • 延迟误差协方差矩阵 $P_j(t) = E[e(t)e^T(t-jT) | \mathcal{F}_t]$：描述当前时刻误差 $e(t)$ 与延迟时刻 $t-jT$ 误差 $e(t-jT)$ 之间的互协方差。特别是，对于单延迟系统，$P_1(t)$ 尤为关键。当 $t < jT$ 时，设定 $P_j(t) = 0$。

3. 最优滤波器方程（定理3）： 基于上述定义，给出了巨腔系统的最优滤波器方程，其形式类似于带有延迟项的卡尔曼滤波器： $ d\hat{x}(t) = A\hat{x}(t)dt + A_d \hat{x}(t-T)dt + G(t)(D_d D_d^T)^{-1} (dy(t) - C\hat{x}(t)dt - C_d \hat{x}(t-T)dt) $ (方程22) 其中 $G(t) = P(t)C^T + P_1(t)C_d^T$ 是滤波器增益，它明确地包含了当前误差协方差 $P(t)$ 和延迟误差协方差 $P_1(t)$ 的贡献。这与传统卡尔曼滤波器中增益只依赖于当前协方差矩阵显著不同，直接体现了对延迟和非马尔可夫效应的处理。

4. 误差协方差动力学方程（定理4和定理5）： 滤波器的核心挑战在于计算 $P(t)$ 和 $P_j(t)$。本文推导了这些矩阵的微分方程：

   • $P(t)$ 的动力学方程（方程28）： 描述了当前误差协方差矩阵的演化，它依赖于 $P(t)$ 自身以及 $P_1(t)$。 $ \dot{P}(t) = (A - G(t)C)P(t) + P(t)(A - G(t)C)^T + (A_d - G(t)C_d)P_1(t) + P_1(t)(A_d - G(t)C_d)^T + BB^T + (B_d - G(t)D_d)(B_d - G(t)D_d)^T $
   • $P_j(t)$ 的动力学方程（方程36）： 这是整个框架的关键创新点，描述了延迟误差协方差矩阵的演化。 $ \dot{P}_j(t) = (A - G(t)C)P_j(t) + (A_d - G(t)C_d) P_{j-1}(t-T) + P_j(t)(A - G(t-jT)C)^T + P_{j+1}(t)(A_d - G(t-jT)C_d)^T + (B_d - G(t)D_d)B^T \mathbb{I}_{j=1} $ 其中 $\mathbb{I}_{j=1}$ 表示当 $j=1$ 时该项存在。这个方程的关键在于 $P_j(t)$ 的演化不仅依赖于 $P_j(t)$ 自身，还依赖于前一个延迟阶数 $P_{j-1}(t-T)$ 和后一个延迟阶数 $P_{j+1}(t)$。这种复杂的相互依赖性要求采用一种迭代的计算策略。

5. 区间式逆向递归算法（算法1）： 为了解决 $P_j(t)$ 方程的耦合性，本文提出了一个创新的区间式逆向递归算法。

   • 步骤1： 初始化滤波器状态 $\hat{x}(0)$ 和延迟协方差矩阵（例如 $P_j(t)=0$ 当 $t
   • 步骤2-3： 在每个时间步 $t$ 下，从最大的延迟阶数 $m$ 开始（例如 $m = \lceil T_{final}/T \rceil$），逆向递归计算 $P_j(t)$。即，先计算 $P_m(t)$，然后利用 $P_m(t)$ 计算 $P_{m-1}(t)$，以此类推，直到 $P_0(t)=P(t)$。
     • 这种逆向递归计算确保了在计算 $P_j(t)$ 时，其所需的所有未来延迟协方差（如 $P_{j+1}(t)$）都已经从上一步的递归中得到，并且过去的延迟信息（如 $P_{j-1}(t-T)$）在相应的延迟时刻也是可用的。
   • 步骤4： 一旦所有的 $P_j(t)$ (特别是 $P(t)$ 和 $P_1(t)$) 在当前时间步得到更新，便可计算滤波器增益 $G(t)$，并根据最优滤波器方程更新状态估计 $\hat{x}(t)$。

这种直接处理延迟的策略，避免了增广系统带来的维度爆炸，同时通过迭代计算延迟协方差矩阵，有效地在过滤框架中内化了非对易性和记忆效应。该方法的严谨性和系统性为非马尔可夫量子系统的实时状态估计提供了有力的工具。

2. 关键 benchmark 体系，计算所得数据，性能数据

为了全面评估所提出最优滤波器在处理量子波导QED巨腔系统时的有效性，本文采用两种具有代表性的量子态作为基准系统进行数值模拟：相干态和薛定谔猫态。这两种状态分别代表了量子力学中的“经典极限”和“非经典特性”的典范，能够充分检验滤波器在不同量子态类型下的跟踪性能。

基准系统描述：

1. 相干态 (Coherent States)： 相干态在量子力学中扮演着准经典态的角色，其在相空间中的Wigner函数呈现高斯分布。Wigner函数 $W(q,p)$ 提供了相空间中量子态的准概率分布，对于相干态，其中心位于期望值 $(\langle q \rangle, \langle p \rangle)$，表达式为： $ W(q,p) = \frac{2}{\pi} \exp(-2 [(q - \langle q \rangle)^2 + (p - \langle p \rangle)^2]) $ (方程43) 这里 $q$ 和 $p$ 分别是位置和动量正交分量。相干态的Wigner函数始终为正，表明其缺乏非经典特性。

2. 薛定谔猫态 (Schrödinger Cat States)： 薛定谔猫态是一种宏观上可区分相干态的量子叠加态（例如 $| \alpha \rangle$ 和 $| -\alpha \rangle$），是量子非经典性的重要标志。本文考虑偶数薛定谔猫态： $ |\text{cat}\rangle = \frac{1}{\sqrt{2(1+e^{-2|\alpha|^2})}} (|\alpha\rangle + |-\alpha\rangle) $ (方程44) 其中 $\alpha \in \mathbb{C}$ 是相干态的振幅。薛定谔猫态的Wigner函数更复杂，由两个相干态的Wigner函数贡献 $W_{\text{coh}}(q,p)$ 和一个干涉项 $W_{\text{int}}(q,p)$ 组成： $ W(q,p) = \frac{1}{2} (W_{\text{coh}}(q,p) + W_{\text{coh}}(q,p) + W_{\text{int}}(q,p)) $ (方程48) 其中 $W_{\text{coh}}(q,p) \propto \exp\left(-\frac{(q-q_0 \pm \beta)^2}{\sigma^2} - \sigma^2(p-p_0)^2\right)$ (方程46)，描述了两个高斯峰；$W_{\text{int}}(q,p) \propto \exp\left(-\frac{(q-q_0)^2}{\sigma^2} - \sigma^2(p-p_0)^2\right) \cos\left(\frac{2\beta}{\sigma^2}(q-q_0) - 2\beta\sigma^2(p-p_0)\right)$ (方程47)，产生相空间中的振荡条纹。干涉项的存在和Wigner函数的负值区域是薛定谔猫态非经典性的关键特征。滤波器能否准确捕捉并跟踪这些负值区域和干涉条纹，是检验其处理非经典量子态能力的重要指标。

计算参数：

数值模拟采用以下系统参数：

• 腔频率 $\omega_c = 10^9 \text{ Hz}$
• 衰减率 $\gamma = 8 \times 10^8 \text{ Hz}$
• 群速度 $v_g = 10^3 \text{ m/s}$
• 时间延迟 $T = 1.5 \times 10^{-8} \text{ s}$，对应于耦合点间距 $L = 1.5 \times 10^{-5} \text{ m}$。
• 初始系统状态和滤波器估计：[-4m, 4kg·m/s]^T 和 [4m, -4kg·m/s]^T。注意这里的单位表示可能存在排版错误，通常是 $q$ 为位置（如米），$p$ 为动量（如 $\text{kg} \cdot \text{m/s}$），但此处给出的数值单位是相同的。我们遵循原文的数值设置。
• 模拟时间间隔：$[0, 1 \times 10^{-7} \text{ s}]$，其中 $T_f = 1 \times 10^{-7} \text{ s}$ 为最终时间。

计算所得数据与性能数据：

1. 相干态的动力学演化与跟踪性能 (图2和图3)：

   • 图2：位置和动量算符的动态演化

     • 观察数据： 图2展示了系统和滤波器估计的位置 $q(t)$ 和动量 $p(t)$ 算符随时间 $t$ 的动态演化。浅蓝色虚线代表原始系统的多条随机状态轨迹，深蓝色曲线是原始系统的平均轨迹。浅红色虚线是滤波器估计的多条轨迹，红色曲线是滤波器状态的平均轨迹。
     • 性能分析：
       • 收敛性： 尽管在 $t=0$ 时刻滤波器初始状态与原始系统状态存在显著差异，但滤波器状态能够迅速收敛并准确跟踪原始系统状态的平均轨迹。这表明滤波器具有良好的初始条件鲁棒性和快速收敛能力。
       • 非马尔可夫特征捕获： 曲线包络线在大约 $1.5 \times 10^{-8} \text{ s}$（即时间延迟 $T$）附近表现出明显的振荡，这与马尔可夫系统中通常观察到的单调衰减行为截然不同。这种振荡是系统非马尔可夫动力学的直接体现，证明了所设计的滤波器能够成功捕获并跟踪这种由时间延迟引起的独特动态特征。
       • 跟踪精度： 随着时间演化，滤波器的平均轨迹与原始系统的平均轨迹高度重合，表明其在连续测量下保持了出色的跟踪精度。

   • 图3：相干态Wigner函数随时间演化

     • 观察数据： 图3(a)-(d)展示了原始系统相干态的Wigner函数在 $t=0, 0.01T_f, 0.5T_f, T_f$ 时刻的演化。图3(e)-(h)则展示了滤波器系统对应的Wigner函数演化。
     • 性能分析：
       • 初始差异与逐步匹配： 在初始时刻 $t=0$ (图3(a)与3(e))，原始系统与滤波器估计的Wigner函数中心存在明显差异。然而，随着时间推移，尤其在 $t=0.01T_f$ (图3(b)与3(f)) 和 $t=0.5T_f$ (图3(c)与3(g)) 时刻，滤波器估计的相干态Wigner函数逐渐向原始系统逼近。
       • 最终高度匹配： 到最终时刻 $t=T_f$ (图3(d)与3(h))，滤波器的Wigner函数几乎完美地匹配了原始系统的Wigner函数，无论是中心位置还是高斯分布的形状。这直观地验证了滤波器在相空间中对系统状态的准确估计能力。

2. 薛定谔猫态的动力学演化与跟踪性能 (图4)：

   • 观察数据： 图4(a)-(b)展示了原始系统薛定谔猫态的Wigner函数在 $t=0$ 和 $t=T_f$ 时刻的演化。图4(c)-(d)展示了滤波器系统对应的Wigner函数在 $t=0$ 和 $t=T_f$ 时刻的演化。图中给出了薛定谔猫态的参数 ($q_0, p_0, \beta, \sigma$)，它们随时间演化。
     • 薛定谔猫态参数在图4中的变化：
       • 原始系统在 $t=0$：$q_0 = -4, p_0 = 4, \beta = 0.8, \sigma = 0.2$。
       • 原始系统在 $t=T_f$：$q_0 = 0.29, p_0 = -0.71, \beta = 0.08, \sigma = 0.02$。
       • 滤波器在 $t=0$：$q_0 = 4, p_0 = -4, \beta = 0.8, \sigma = 0.2$。
       • 滤波器在 $t=T_f$：$q_0 = 0.29, p_0 = -0.68, \beta = 0.08, \sigma = 0.02$。
     • 这些参数的变化反映了薛定谔猫态在腔中演化（可能由于衰减或与环境的相互作用）导致其相空间分布的平移、分离度 $\beta$ 的变化以及高斯包宽 $\sigma$ 的变化。
   • 性能分析：
     • 初始差异与非经典特征的跟踪： 在初始时刻 $t=0$ (图4(a)与4(c))，原始系统和滤波器估计的猫态Wigner函数同样存在显著差异，但两者都清晰地展示了薛定谔猫态特有的干涉条纹和负值区域。这意味着滤波器能够识别和表示这种非经典特性。
     • 最终高度准确估计： 到最终时刻 $t=T_f$ (图4(b)与4(d))，滤波器估计的薛定谔猫态Wigner函数与原始系统的高度吻合。特别重要的是，滤波器不仅准确跟踪了猫态相空间中心的移动 ($q_0, p_0$)，还精确重现了其标志性的干涉条纹结构和负值区域。这表明滤波器能够成功地估计并保持了薛定谔猫态的量子相干性，这是处理非经典量子信息至关重要的能力。

结论：

通过对相干态和薛定谔猫态的数值模拟，本文所提出的最优滤波器展示了卓越的性能。它不仅能够对经典的相干态进行准确跟踪，还能有效地估计和保留非经典的薛定谔猫态的复杂量子特征，包括其特有的干涉条纹和负值区域。滤波器对非马尔可夫动力学（如振荡行为）的成功捕获进一步证实了其设计的鲁棒性和有效性。这些结果为在复杂非马尔可夫背景下实现量子态的精确控制和利用奠定了坚实的基础。

3. 代码实现细节，复现指南，所用的软件包及开源 repo link

本文作为一项理论研究，主要贡献在于提出了一个严谨的数学框架和算法，并使用数值模拟验证了其有效性。因此，原文并未提供具体的代码实现、复现指南或开源软件库链接。然而，基于论文中详尽的理论推导和算法描述（特别是算法1），我们可以清晰地构想出该滤波器的实现细节和复现步骤。

代码实现的概念性步骤：

实现该最优滤波器主要涉及数值求解一系列耦合的随机微分方程（SDEs）和常微分方程（ODEs）。

1. 系统参数定义与矩阵构建：

   • 首先，定义所有物理常数，如 $\omega_c, \gamma, v_g, T$。
   • 根据论文第2.4节的描述，从原始的复数算符方程（方程19）构建实值系数矩阵 $A, A_d, B, B_d, C, C_d, D_d$。这需要应用 $M(\cdot)$ 运算符（在第1节末定义），将复数系数转换为2x2的实数矩阵。例如，$A = M(\bar{A})$ 等。这是实现的第一步，也是确保正确性的基础。

2. 数值积分方案选择：

   • 由于系统动力学由微分方程描述，需要选择合适的数值积分方法。对于随机微分方程（如 $\hat{x}(t)$ 的方程），可以考虑Euler-Maruyama 方法。对于常微分方程（如 $P(t)$ 和 $P_j(t)$ 的方程），Runge-Kutta 方法（如4阶或5阶）通常能提供更好的精度和稳定性。时间步长 $dt$ 的选择需要足够小，以确保数值稳定性，特别是对于具有延迟的系统。

3. 最优滤波器算法 (Algorithm 1) 实现： 这是实现的核心，涉及一个嵌套的时间循环结构：

   • 外层时间循环（正向）： 迭代从 $t=0$ 到 $T_f$ 的每个时间步。

     • for t in range(0, Tf, dt):

   • 内层延迟协方差逆向递归循环（$P_j(t)$ 的计算）： 在每个外层时间步 $t$ 内部，需要计算当前时刻 $t$ 的所有延迟协方差矩阵 $P_j(t)$。论文的算法1描述了一个从最大延迟阶数 $m = \lceil T_{final}/T \rceil$ 向下到 $j=0$ 的逆向递归过程。这意味着在实际实现中，我们需要在一个小的时间窗口内（例如 $[t, t+dt]$），对所有 $P_j$ 方程进行迭代求解。

     • 初始化 $P_j$： 对于 $t < jT$ 的情况，$P_j(t)$ 设为零矩阵。这需要一个机制来存储和检索过去时刻的协方差矩阵。
     • 逆向递归求解 $dP_j(t)$：
       • for j in range(m, -1, -1):
       • 这里，$P_j(t)$ 的导数 $\dot{P}_j(t)$ (方程36) 依赖于 $P_j(t)$ 自身、$P_{j-1}(t-T)$ 和 $P_{j+1}(t)$。
       • 在实际数值求解时，这通常意味着在一个时间步长内，需要先计算所有 $P_j$ 在 $t+dt$ 的值。由于 $P_j$ 导数表达式中的相互依赖性，这可能需要一个隐式求解器，或者通过精心的离散化来处理。
       • 关键点： $P_{j-1}(t-T)$ 需要从过去的时间步中检索。这暗示着我们需要维护一个历史缓冲，存储所有 $P_j$ 在过去 $m$ 个延迟周期内的值。
       • 边界条件： $P_{m+1}(t)$ 被视为零（因为它超出了考虑的延迟范围）。$P_0(t)$ 就是 $P(t)$。

   • 计算 $P(t)$：

     • 在内层循环结束后，得到 $P_0(t)$ 和 $P_1(t)$。接着，根据方程28计算 $P(t)$ 的导数 $\dot{P}(t)$，并更新 $P(t)$ 到下一个时间步。

   • 计算滤波器增益 $G(t)$：

     • 根据 $P(t)$ 和 $P_1(t)$ (方程27) 计算 $G(t) = [P(t)C^T + P_1(t)C_d^T](D_d D_d^T)^{-1}$。

   • 更新滤波器状态 $\hat{x}(t)$：

     • 根据方程22，利用当前的 $\hat{x}(t)$、延迟的 $\hat{x}(t-T)$、计算出的 $G(t)$ 和实测的 $dy(t)$ 来更新 $\hat{x}(t)$ 到下一个时间步。
     • $dy(t)$ 必须由模拟的真实系统生成。

4. 真实系统模拟（生成 $dy(t)$）：

   • 为了验证滤波器的性能，需要一个“真实”的量子系统模型。这包括根据方程20模拟 $x(t)$ 的演化。
   • 量子白噪声生成： $dw(t)$ 需要模拟为高斯白噪声。在数值模拟中，这意味着在每个时间步生成一个方差为 $dt$ 的随机数。对于量子白噪声，它通常服从零均值高斯分布，并且其增量方差由量子涨落决定。
   • 通过真实系统方程 $dy(t) = Cx(t)dt + C_d x(t-T)dt + D_d dw(t-T)$ 生成测量输出 $dy(t)$，供滤波器使用。

5. 性能评估与可视化：

   • Wigner函数计算：
     • 对于相干态，Wigner函数（方程43）可以直接从估计的 $\langle q \rangle = \hat{x}_1(t)$ 和 $\langle p \rangle = \hat{x}_2(t)$ 计算。
     • 对于薛定谔猫态，Wigner函数（方程48, 46, 47）需要更多的参数 ($q_0, p_0, \beta, \sigma$)。论文中这些参数是给定的，但在一个真实的滤波器中，可能需要从估计的 $\hat{x}(t)$ 和 $P(t)$ 中推断或拟合出这些参数，以重构Wigner函数。最直接的方法是假定猫态的结构参数已知，滤波器仅估计其中心。如果滤波器要跟踪猫态更复杂的演化，那么需要更高级的方法来从滤波器估计中获取这些参数。
   • 绘图： 使用绘图库（如Matplotlib）生成类似论文中图2、图3和图4的跟踪曲线和Wigner函数分布图。

所用的软件包及开源 repo link：

由于没有提供实际代码，我们将列出用于实现此类复杂数值模拟和量子系统建模的常用软件包和编程语言。

• 编程语言：

  • Python： 社区庞大，拥有丰富的科学计算库。
  • MATLAB： 在控制系统和数值分析领域非常流行，内置的ODE求解器和矩阵操作非常强大。
  • C++/Fortran： 对于性能要求极高的实时系统或大规模模拟，C++或Fortran是更好的选择，但开发周期较长。

• Python 关键软件包：

  • NumPy： 基础的数值计算库，用于高效的矩阵和向量操作。
  • SciPy： 包含各种科学和工程模块，尤其是：
    • scipy.integrate.solve_ivp：用于求解ODE（如 $P(t)$ 和 $P_j(t)$ 的方程）。
    • scipy.linalg：用于线性代数运算，如矩阵求逆。
  • Matplotlib / Seaborn： 数据可视化和绘图库，用于生成论文中的图表。
  • SymPy (可选)： 符号计算库，可用于辅助推导和验证系数矩阵的转换。
  • QuTiP (可选)： 虽然该滤波器不是基于密度矩阵的，但QuTiP提供了一些量子态（如相干态、猫态）的生成和Wigner函数计算功能，可用于辅助生成真值和可视化。

• MATLAB 关键功能：

  • ode45 等ODE求解器： 用于求解微分方程。
  • 矩阵运算： MATLAB对矩阵运算原生支持良好。
  • 绘图工具： 强大的二维和三维绘图功能。

复现指南（假设使用Python）：

1. 环境准备： 安装 Python 解释器和上述列出的科学计算库（NumPy, SciPy, Matplotlib）。推荐使用 Anaconda 或 Miniconda 管理环境。

2. 文件结构： 建议创建如下文件结构：

   • main.py：主脚本，包含模拟逻辑。
   • system_model.py：定义系统常数和系数矩阵的函数。
   • filter_equations.py：定义 $\dot{x}(t), \dot{P}(t), \dot{P}_j(t)$ 等导数函数的模块。
   • utils.py：辅助函数，如Wigner函数计算。

3. 逐步实现：

   • 第一阶段：系统模型。 编写函数来计算 $A, A_d, B, B_d, C, C_d, D_d$ 矩阵。
   • 第二阶段：真实系统模拟。 编写一个模拟函数，接收初始 $x_0$ 和 $T$，使用Euler-Maruyama方法集成方程20，同时生成 $dw(t)$ 和 $dy(t)$。存储 $x(t)$ 和 $dy(t)$ 的历史轨迹。
   • 第三阶段：滤波器核心。
     • 编写函数来计算 $\dot{P}(t)$ 和 $\dot{P}_j(t)$ 的右侧（即方程28和36）。这些函数将依赖于当前 $P(t)$, $P_j(t)$ 以及历史 $P_{j-1}(t-T)$ 和 $P_{j+1}(t)$。
     • 实现Algorithm 1的嵌套循环逻辑。这需要一个数据结构（如一个列表或字典）来存储不同延迟阶数 $j$ 的 $P_j(t)$ 矩阵历史，以便在 $t-T$ 时刻检索。
     • 编写滤波器状态更新函数，实现方程22。
   • 第四阶段：Wigner函数计算。 实现Wigner函数计算（方程43, 46, 47），用于评估滤波器的输出。
   • 第五阶段：可视化。 使用Matplotlib绘制结果，复现论文中的图表。

挑战与注意事项：

• 数值稳定性： Riccati型方程的数值求解可能会遇到稳定性问题。选择合适的积分器和时间步长至关重要。
• 内存管理： 存储 $P_j(t)$ 的历史轨迹，特别是对于较长延迟或多个延迟阶数 $m$，可能会占用大量内存。需要优化数据结构。
• 计算效率： 逆向递归求解 $P_j(t)$ 是计算密集型的。对于实时应用，可能需要并行化或考虑更高效的数值方法。
• 噪声实现： 正确实现量子白噪声的统计特性对于模拟的准确性至关重要。

由于没有开源代码，复现需要对论文的数学细节有深入理解，并具备扎实的数值模拟编程能力。一个完全开源的实现将极大地促进该领域的研究和应用。

4. 关键引用文献，以及你对这项工作局限性的评论

本文在量子滤波领域做出了重要贡献，其理论框架建立在多年来量子信息和控制理论的深厚积累之上。以下是论文中一些关键引用文献的意义及其与本文工作的关联，以及我对这项工作局限性的评论。

关键引用文献及其意义：

1. [7] V. P. Belavkin (1989), “Optimal filtering of Markov signals with quantum noise.”

   • 意义： Belavkin是量子滤波理论的奠基人，他首次将经典随机过程和滤波理论扩展到非对易的量子领域，为量子卡尔曼滤波等现代量子滤波方法奠定了数学基础。
   • 关联： 本文的工作是在Belavkin开创性工作的基础上，针对更复杂的非马尔可夫量子系统进行拓展。Belavkin的原始框架主要处理马尔可夫系统，而本文则通过引入延迟协方差矩阵来处理非马尔可夫效应。

2. [5] H. M. Wiseman and G. J. Milburn (2010), “Quantum Measurement and Control.”

   • 意义： 这是一本经典的教科书，全面介绍了量子测量理论和量子反馈控制，是量子滤波和控制领域的标准参考。它深入探讨了量子测量中的反作用和量子斯托克斯方程。
   • 关联： 本文的背景理论和术语，如海森堡绘景中的朗之万方程、输入-输出关系以及量子测量反作用的考虑，都深受Wiseman和Milburn工作的影响。

3. [8] L. Bouten, R. van Handel, and M. R. James (2007), “An introduction to quantum filtering.”

   • 意义： 该综述为量子滤波理论提供了严谨的数学介绍，强调了量子力学中算符非对易性对滤波器设计的影响，这是经典滤波理论所不具备的特点。
   • 关联： 本文明确指出“耦合和传播延迟导致系统算符在不同时间的非对易性，这使得现有量子滤波技术无法直接应用”。Bouten等人的工作强调了非对易性是量子滤波的核心挑战，而本文正是通过其独特的方法来解决由延迟引起的非对易性问题。

4. [13] S. Xue et al. (2015), “Quantum filter for a class of non-Markovian quantum systems.” [14] S. Xue et al. (2019), “Modeling for non-Markovian quantum systems.”

   • 意义： 这些是本文作者团队之前关于非马尔可夫量子系统量子滤波器的研究，通常采用增广模型（即通过引入辅助系统将非马尔可夫环境转化为马尔可夫系统）。
   • 关联： 本文对比了这种“增广系统法”的局限性（高维和计算负担），并提出了新的方法，即直接在原始系统模型中处理延迟，避免了增广。这表明本文是该团队在非马尔可夫量子滤波研究上的一个重要进展。

5. [28] L. Guo et al. (2017), “Giant acoustic atom: A single quantum system with a deterministic time delay.” [21] L. Guo et al. (2020), “Oscillating bound states for a giant atom.”

   • 意义： 这些论文探讨了“巨型原子”或“巨腔”系统的独特物理特性，例如确定的时间延迟和振荡束缚态，这些特性是本文研究对象——巨腔系统——的基础。
   • 关联： 本文的巨腔系统正是这类“巨型原子”系统在波导QED中的一种实现，其非马尔可夫动力学源于这些独特的延迟相互作用。对这些物理现象的理解构成了本文系统建模的物理基础。

6. [29] C. W. Gardiner and M. J. Collett (1985), “Input and output in damped quantum systems…”

   • 意义： 这篇论文为开放量子系统的输入-输出理论和量子随机微分方程奠定了基础，引入了广为人知的量子朗之万方程。
   • 关联： 本文的朗之万方程和输入-输出关系的推导（特别是马尔可夫近似下的简化形式，如本文Remark 2所述）直接来源于Gardiner和Collett的工作。

7. [31] M. Basin et al. (2005), “Optimal filtering for linear state delay systems.”

   • 意义： 这篇论文在经典控制理论中讨论了带有状态延迟的线性系统的最优滤波问题。
   • 关联： 本文在Remark 4中将其与本文的量子滤波器进行了比较，指出本文的滤波器增益除了 $P(t)$ 外还包含 $P_1(t)$，这是经典延迟系统滤波器所不具备的，突出了量子延迟系统和非对易性的特殊处理。

对这项工作局限性的评论：

尽管本文提出了一种创新且有效的非马尔可夫量子滤波框架，并成功地通过数值模拟验证了其对相干态和薛定谔猫态的良好跟踪性能，但其仍然存在一些值得讨论的局限性：

1. 计算复杂度和可扩展性：

   • 迭代逆向递归的计算成本： 虽然本文避免了增广希尔伯特空间带来的高维问题，但其核心的区间式逆向递归算法（Algorithm 1）需要迭代计算多个延迟阶数 $j$ 的 $P_j(t)$ 矩阵。对于每个时间步长，都需要对这些耦合的Riccati型微分方程进行求解。当系统延迟周期 $T$ 很长或者模拟总时间 $T_f$ 很长，导致 $m=\lceil T_f/T \rceil$ 很大时，这种计算量将变得非常庞大。论文没有提供计算复杂度的详细分析或实际运行时间数据，这使得评估其在更大规模或实时应用中的可行性存在困难。
   • 多耦合点/多腔系统的挑战： 论文在结论中指出未来的工作将扩展到“多点和多腔耦合的复杂量子网络系统”。可以预见，随着耦合点数量或腔数量的增加，系统状态的维度会增加，更重要的是，延迟的种类和相互作用会更加复杂，可能需要更多类型的延迟协方差矩阵 $P_{j,k}(t)$，这将指数级地增加算法的计算负担。

2. 仅限于线性高斯系统：

   • 本文的滤波器是为线性系统和高斯量子噪声（量子白噪声）推导的。许多实际的量子系统可能存在显著的非线性效应，或者环境噪声并非严格高斯分布。在这种情况下，本文提出的滤波器将不再是最优的，可能需要更复杂的非线性滤波技术，如扩展卡尔曼滤波器（EKF）的量子版本或量子粒子滤波器。

3. 固定的时间延迟 $T$：

   • 模型假设存在一个固定的、确定的时间延迟 $T = L/v_g$。在现实物理系统中，传播延迟可能不是恒定的，或者存在分布式的延迟（即延迟不是单一值，而是一个范围或分布）。处理这些更复杂的时间延迟结构将需要对当前框架进行进一步的泛化。

4. 同源测量 (Homodyne Detection) 的假设：

   • 滤波器的输出测量模型是基于同源测量（homodyne detection）的，这是一种常见的连续量子测量技术。然而，其他类型的测量（如光子计数、外差探测或更广义的量子测量）可能具有不同的测量算符和噪声特性，需要对滤波器进行相应的修改。

5. 缺乏实验验证：

   • 本文的有效性完全通过数值模拟得到验证。从理论到实验的转化通常会引入新的挑战，例如实际系统参数的精确校准、环境噪声的复杂性、测量装置的效率和可靠性，以及量子测量反作用在实验中的实际影响。这些因素都可能影响理论滤波器的实际性能。

6. 薛定谔猫态Wigner函数参数处理：

   • 在猫态的Wigner函数可视化部分，论文在图注中直接给出了不同时间点的 $q_0, p_0, \beta, \sigma$ 等参数。一个“最优”的滤波器通常应该能够从测量数据中估计出这些描述猫态的关键参数，而不仅仅是中心点 $q, p$。论文似乎是模拟了猫态在特定参数下演化，然后滤波器对其进行跟踪，但如何从滤波器估计的 $\hat{x}(t)$ 和 $P(t)$ 中推断或拟合出完整的猫态Wigner函数参数，或者滤波器是否能直接跟踪这些参数的动态演化，没有详细说明。这可能意味着滤波器在处理猫态时，其估计的完整性仍有提升空间。

总的来说，本文为处理非马尔可夫量子系统中的延迟效应提供了一个强有力的新范式，但未来的研究需要解决其计算效率、对非线性的适应性以及实验验证等方面的挑战，以使其在更广泛的量子技术应用中发挥潜力。

5. 其他你认为必要的补充

本文在量子滤波领域取得的突破性进展，不仅为波导QED巨腔系统提供了精确的量子态估计工具，更对整个量子信息科学和工程领域具有深远的影响。以下将从更广泛的视角对本文的重要性、潜在应用、未来展望以及其理论深度进行补充阐述。

A. 广泛影响和重要性

1. 超越马尔可夫近似的量子物理理解： 大多数量子开放系统理论和控制方法都基于马尔可夫近似，即假设环境无记忆或记忆时间极短。然而，本文研究的巨腔系统明确展示了由远程耦合和有限传播速度引起的内在记忆效应。通过直接在动力学方程中引入延迟并设计相应的滤波器，本文提供了一种处理非马尔可夫效应的通用范式，这对于理解和控制与结构化或强耦合环境相互作用的真实量子系统至关重要。它推动了量子物理学对开放系统动力学边界的探索。

2. 量子技术发展的基石：

   • 量子计算与容错： 精确的量子态估计是量子计算中进行量子态层析、错误检测和校正的基础。在具有长程相互作用和延迟的量子计算架构（如基于波导连接的超导量子比特）中，本文的滤波技术可以提供实时的反馈信息，对于实现容错量子计算具有不可估量的价值。
   • 量子计量学与传感： 量子滤波器能够从噪声测量中提取最优状态信息，从而提高量子传感器的精度和灵敏度。在巨腔系统中，这种增强的滤波能力可以用于开发更精确的量子计时、重力梯度测量或其他物理量传感。
   • 量子网络： 巨腔系统及其远程耦合是构建分布式量子网络和量子中继器的潜在构件。在这些网络中，节点间的信号传输必然存在延迟。本文的滤波框架为在这种延迟环境下维护和估计网络中量子态的纠缠和相干性提供了关键工具。

3. 量子控制理论的新范式： 本文的滤波器设计方法为延迟量子反馈控制系统奠定了基础。在反馈控制中，控制器接收系统输出测量，然后生成控制输入。如果测量或控制信号的传播存在延迟，那么标准的马尔可夫控制理论会失效。本文的滤波框架，通过准确估计延迟系统的状态，可以与后续的量子控制器结合，实现对延迟量子系统的稳定化、态制备或门操作。这种将经典控制理论（如卡尔曼滤波）扩展到复杂的延迟量子系统，是控制理论领域的一项重要突破。

B. 理论深度与创新点

1. 直接建模延迟而非增广： 这是本文最核心的理论创新之一。传统的非马尔可夫问题通常通过将环境的记忆效应“增广”到系统 Hilbert 空间中来解决，这导致系统维度急剧膨胀，计算成本高昂。本文则选择直接在朗之万方程中引入延迟项，并以此为基础推导出包括延迟协方差矩阵 $P_j(t)$ 在内的全套滤波方程。这种方法在保留系统物理特性（非局域耦合、多延迟动力学）的同时，避免了不必要的维度增加，是解决复杂非马尔可夫问题的优雅方案。

2. 延迟误差协方差矩阵的系统性处理： 引入并系统性地计算 $P_j(t) = E[e(t)e^T(t-jT) | \mathcal{F}_t]$ 是解决延迟引起非对易性的关键。这些矩阵刻画了当前时刻的估计误差与过去延迟时刻的误差之间的统计相关性。其动力学方程（定理5）清晰地展示了不同延迟阶数 $P_j(t)$ 之间的耦合关系。通过设计精巧的“区间式逆向递归算法”，本文提供了一个在计算上可行且理论上一致的策略来更新这些协方差矩阵，从而使滤波器能够准确地追踪系统状态，即使在存在显著记忆效应的情况下。

3. 对非对易性的内化处理： 量子系统中的核心特征是算符的非对易性。在有延迟的系统中，不同时间点的算符之间的对易关系会变得复杂。本文通过在海森堡绘景中进行推导，然后转换到实值正交分量表示，并推导包含所有交叉协方差项的Riccati型方程，巧妙地内化了这种非对易性。这意味着在实值矩阵运算层面，滤波器已经“感知”并处理了底层量子算符的非对易特性，无需在每个时间步进行显式的量子对易子计算，从而使得算法在形式上保持了卡尔曼滤波的简洁性和可操作性。

C. 潜在的未来研究方向和应用拓展

1. 实验验证与平台实现： 将本文的理论框架应用于实际的波导QED实验平台是至关重要的下一步。这可能涉及超导电路QED、光子晶体波导或声学表面波器件。实验验证将需要解决诸如精确控制延迟、高保真测量、以及在真实噪声背景下维持系统相干性等挑战。这将不仅验证理论，也将推动实验技术的发展。

2. 与量子控制的融合： 本文的滤波器是量子控制回路中的“观测器”。未来的工作可以将此滤波器与动态规划、强化学习或最优控制理论相结合，设计出能够主动利用延迟效应、甚至通过反馈来塑造非马尔可夫动力学的控制器。例如，利用延迟实现稳定化、纠缠生成或量子门操作。

3. 非线性与非高斯扩展： 尽管本文聚焦于线性高斯系统，但现实世界的量子器件可能存在非线性响应，或者环境噪声可能呈现非高斯特性。未来的研究可以探索将本文的延迟滤波框架扩展到非线性系统，例如通过量子扩展卡尔曼滤波器或量子粒子滤波器的概念进行推广。

4. 多体与多延迟网络： 正如论文结论所述，将框架扩展到多腔、多耦合点的复杂量子网络是重要方向。这需要开发更通用的延迟模型，例如分布式延迟或时间依赖的延迟，并处理由此带来的更高维度的延迟协方差矩阵系统。

5. 实时计算优化： 为了实现真正的实时反馈控制，需要对滤波器算法的计算效率进行优化。这可能包括开发近似滤波方案、利用并行计算架构（如GPU）或量子计算机辅助计算等，以加速延迟协方差矩阵的迭代求解。

6. 噪声谱与环境工程： 本工作为理解和利用波导QED中结构化环境的记忆效应提供了工具。未来的研究可以探索如何通过工程化波导环境的色散关系或耦合点布局，来主动设计有利的非马尔可夫动力学，并利用滤波器对其进行精确跟踪和控制，从而实现新颖的量子态制备或量子信息处理协议。

总而言之，本文不仅解决了量子波导QED巨腔系统中一个复杂的非马尔可夫滤波问题，更提供了一个处理广泛延迟量子系统的新思路和强大的工具集。其深远的理论影响和广阔的应用前景，预示着它将在未来的量子技术发展中扮演关键角色。