企业官网建设流程全解析

第一次用 DMD 在 Matlab 里算出库普曼算子谱的时候，我盯着特征值图看了很久——一个混沌的 Lorenz 系统，居然也能有这么漂亮的谱结构。原来非线性系统在“函数空间”这一层处理之后，可以变成一个线性算子的问题。这篇文章想把遍历论、动态模态分解和库普曼算子谱特性这三件事串在一起讲清楚，并给出相应的 Matlab 代码实现。目标是让你在现有数据的基础上，算出系统的 Koopman 特征值、特征函数与模态，并用来判断系统是周期的、拟周期的、耗散的还是混合/混沌的。适合正在研究非线性动力学、流场降阶、系统辨识，或者单纯想知道 DMD 背后那套数学究竟在算什么的人。

1. 线性化的另一条路：Koopman 算子为什么值得做

做动力系统研究的人，早就习惯了“状态空间”视角。拿到一个非线性系统，要么画相图，要么算雅可比矩阵看局部稳定性，要么求李雅普诺夫指数判断混沌。这些方法都很有用，但它们本质上都围绕“状态点”和“切空间”打转。雅可比线性化只能告诉你某个不动点附近的局部行为；李雅普诺夫指数只能告诉你相邻轨迹的发散速率。问题在于，真实系统里我们往往想要的是一个“全局”解释：这个系统的长期吸引子长什么样，系统里有哪些稳定的振荡频率，瞬态过程会怎么衰减，冷不丁一个扰动会把状态带到哪里去。状态空间视角很难给出这种全局图景。

Koopman 算子的思路是把视角从状态空间切换到函数空间。

考虑离散时间动力系统：

x_{k+1} = F(x_k)

对任意可观测量 g(x)，定义库普曼算子 U 为：

(Ug)(x) = g(F(x))

也就是说，U 把“读取当前状态 x 的函数 g”变成“读取下一时刻状态 F(x) 的函数”。这个定义里没有任何线性化操作，但 U 本身却天然是线性算子，因为函数空间里加法是逐点定义的，标量乘也保持不变。非线性系统 F 哪怕再复杂，它在函数空间诱导出来的 U 都是线性的。代价是函数空间通常是无穷维的，一个有限维非线性系统对应一个无穷维线性算子。这正是整个领域最迷人的地方：用无穷维线性算子的谱理论，去解释有限维非线性系统的全局动力学。

在 Koopman 文献里经常看到“升维”这个说法，就是指这事。我们不是把系统变成更高维的线性系统，而是直接认为系统在函数空间里本来就是线性的，只是我们过去一直站在状态空间里看问题。实际建模中，我们只保留有限个可观测量，比如坐标函数、多项式、傅里叶基，于是无穷维线性算子被投影到有限维子空间上。动态模态分解（DMD）以及扩展 DMD（EDMD），做的就是这件事——用数据去估计这个投影后的有限维矩阵，再求它的谱。

一个极简单的例子可以说明这个思想：单位圆上的旋转：

z_{k+1} = e^{iω} z_k

状态 z 本身是复数，但如果我们定义可观测量 g(z) = z^{-1}，那么 g(z_{k+1}) = e^{-iω} g(z_k)。这个可观测量沿轨迹只乘一个固定相位因子，正好是 Koopman 特征函数的样子，对应的特征值是 λ = e^{-iω}。系统在状态空间里旋转一圈很简单，但它告诉我们的是一种普遍规律：任何系统只要存在可辨认的周期性，Koopman 谱里就一定有落在单位圆上的特征值与之对应。

2. 遍历论给计算带来的约束：不变测度与谱分类

2.1 为什么一定要讲遍历论

很多人学 DMD 时会把谱理论当纯线性代数，跳过遍历论直接进入代码。这个做法在简单周期系统上通常没问题，但一碰到混沌系统就会出乱。原因在于：Koopman 算子的谱特性并不仅仅由映射 F 决定，还取决于我们讨论哪个函数空间，而函数空间又和一个不变测度绑在一起。

不变测度是指满足：

μ(F^{-1}(A)) = μ(A)

的概率测度。它描述了系统长期访问不同区域的时间占比。对任意一条初始轨迹，我们可以用经验平均构造不变测度：

μ̂_N = (1/N) Σ_{k=0}^{N-1} δ_{x_k}

N 足够大时，经验测度会收敛到某个不变测度。这就是为什么我们在做谱估计时，本质上是在对一条长时间轨迹做时间平均。

Birkhoff 遍历定理告诉你，如果系统是遍历的，那么对几乎处处初始点，时间平均等于空间平均：

lim_{N→∞} (1/N) Σ_{k=0}^{N-1} g(x_k) = ∫ g(x) dμ(x)

这个结论在数值上非常实用。它意味着我们可以用一条轨迹的时间平均，去估计系统在全局测度空间里的统计量，而不需要同时模拟成千上万个初始点。

2.2 不变测度的数值感受

拿 Logistic 映射来感觉一下。系统：

x_{k+1} = 4x(1-x)

理论上有著名的绝对连续不变测度：

ρ(x) = 1 / (π√(x(1-x)))

用一条轨迹去逼近它，只需要短短几万步。下面这段简单的 Matlab 代码可以直接验证：

N = 2e5; x = zeros(1, N); x(1) = 0.2; for k = 1:N-1 x(k+1) = 4*x(k)*(1-x(k)); end % 时间平均 S = mean(x); fprintf('Logistic 时间平均 = %.4f, 理论空间平均 = %.4f\n', S, 0.5); % 不变密度直方图与理论密度对比 edges = linspace(0, 1, 61); counts = histcounts(x, edges); pdf = counts / N / (edges(2)-edges(1)); figure; bar(edges(1:end-1), pdf, 'hist', 'FaceAlpha', 0.5); hold on; xx = linspace(0.01, 0.99, 200); plot(xx, 1./(pi*sqrt(xx.*(1-xx))), 'r-', 'LineWidth', 1.6); xlabel('x'); ylabel('密度'); legend('实测', '理论');

运行下来你会发现，时间平均非常接近 1/2，而且直方图在两端凸起，和理论密度的边缘奇异性吻合。这个例子说明，混沌系统不是没有规律，而是它的规律不在单条轨迹上，而在测度意义上。

2.3 谱分类与系统行为对偶

在保测变换的前提下，库普曼算子 U 在 L²(μ) 上是一个等距算子，如果 F 可逆，它还是酉算子。酉算子的谱落在单位圆周 |λ|=1 上。谱的类型和动力学行为之间有一个很粗略但极其有用的对应：

遍历论里有一个重要结论：如果一个系统是混合的，那么除了常数特征函数，L² 里不存在可积的特征函数。所以混合系统对应的 Koopman 谱必然是连续谱。Logistic 映射在参数为 4 时就是典型的混合系统，它的谱不应该用一两个离散特征值去理解，而应该看作一个连续谱的数值近似。

这个判断对实际项目非常关键。很多初学者拿到 DMD 特征值图，看到混沌系统的特征值不像旋转系统那样落在几个孤立点上，就会怀疑代码写错了。实际上不是，这个“云”本身就是混沌系统谱特性的真实表现。后面第五节我会用 Matlab 复现这个现象。

3. 谱特性解读：特征函数如何解码长期动力学

3.1 特征函数是“非线性坐标”

假设我们找到了库普曼算子的特征函数 φ 和特征值 λ，满足：

φ(x_{k+1}) = λ φ(x_k)

如果 |λ| = 1，令 λ = e^{iω}，那么这个特征函数沿轨迹只做相位旋转，幅值不变。这等于在说：存在一个特殊的非线性坐标，在这个坐标下系统的演化竟然是线性的。系统在原始状态空间里绕来绕去，但在 φ 这个坐标下就是匀速旋转。这个 φ 把复杂的非线性运动“拉直”了。

如果 |λ| < 1，特征函数对应的模态随时间指数衰减，它描述了系统的瞬态弛豫过程；如果 |λ| > 1，对应发散模态，这种模态通常意味着不稳定，或者你的数据采集时长不够，把部分趋势误当成模态。

需要注意，特征函数并不是随便找个可观测量就能得到的。只有落在 Koopman 算子谱分解中的那些可观测分量才具有干净的单频演化性质。一个任意可观测量 g 通常同时含有多个特征函数分量。

3.2 特征展开就是“谱域预测公式”

Koopman 谱最实用的地方在于：任意可观测量 g 可以按特征函数展开：

g(x) ≈ Σ_{j=1}^{r} b_j φ_j(x)

沿轨迹演化后：

g(x_k) ≈ Σ_{j=1}^{r} b_j λ_j^k φ_j(x_0)

这本质上是一个指数和模型。给定初始时刻的 φ_j(x_0) 和系数 b_j，就可以预测任意未来时刻的可观测量值。DMD 做的就是把这个展开限制在“坐标函数”张成的子空间里，并用快照数据去拟合 λ_j、φ_j 和 b_j。你不需要知道 φ 的显式解析形式，数据会告诉你一切。

这个公式也是流场降阶模型的理论基础。圆柱绕流中涡脱落对应的 Koopman 特征值，落在单位圆上某个 Strouhal 频率附近，对应两个共轭模态；这两个模态就是流场里周期性上下交替的涡结构。用谱域预测比直接在状态空间做非线性仿真要便宜得多。

3.3 连续谱的坑：不要把每个特征值都当成物理模态

一个很容易犯的错，是把 DMD 画出来的每个特征值都解读为一个真实存在的物理模态。对离散谱系统，这个解读没太大问题；但连续谱系统里，DMD 算出的特征值密集分布，与其说每个点是一个“模式”，不如说这些点的集合是在逼近一个连续频带。有限数据、有限维截断决定了我们只能得到连续谱的离散采样。物理上真正重要的可能是这个谱带位于哪些频率区间，能量集中在哪，而不是单看某一两个特征值。判断一个 DMD 结果靠不靠谱，一定要结合谱的类型：先分清是点谱还是连续谱近似，再去讨论模态的物理含义。

4. 动态模态分解：从数据快照估计 Koopman 谱

4.1 矩阵构造与线性模型假设

动态模态分解的输入是两个快照矩阵：

X = [x_0, x_1, ..., x_{m-1}]，Y = [x_1, x_2, ..., x_m]

其中每个 x_k 是 n 维状态向量。DMD 假设存在一个近似线性算子 A 使得：

x_{k+1} ≈ A x_k

这个 A 在 Koopman 视角下，可以理解为库普曼算子在坐标函数子空间上的一个有限维近似。目标是最小化：

min_A ||Y - AX||_F^2

解是伪逆形式：A = YX^+，其中 X^+ 是 X 的 Moore-Penrose 伪逆。直接对 X 求伪逆在数值上是不稳定的，尤其当系统维数高、数据存在噪声或者状态矩阵接近奇异时。标准做法是用 SVD 把 X 截断到一个低秩子空间再求解。

4.2 SVD 截断背后的动机

对 X 做奇异值分解：

X = UΣV^H

截断保留前 r 个主奇异值：

U_r ∈ C^{n×r}, Σ_r ∈ R^{r×r}, V_r ∈ C^{m×r}

然后投影矩阵：

à = U_r^H Y V_r Σ_r^{-1}

这个 Ã 就是 A 在低秩子空间上的表示。有人会问：为什么不直接用完整 X 的伪逆 A = YX^+？原因主要有两个。第一，真实数据几乎都含噪声，SVD 截断相当于把信噪比低的奇异方向丢弃，避免噪声被逆运算放大；第二，n 通常很大，完整伪逆矩阵既占内存又不稳定，低秩投影能显著降低计算量。秩 r 的选取常用能量准则：保留的奇异值平方和占全部奇异值平方和的 99%，但后面我会说，这个准则在混沌系统中需要额外小心。

4.3 谱分解与模态回映

求 Ã 的特征值分解：

à w_j = λ_j w_j

λ_j 就是 DMD/Koopman 离散时间特征值的估计。对应的高维模态是：

Φ_j = Y V_r Σ_r^{-1} w_j

这里用 Y 而不是 X，是 DMD 算法里一个容易被忽略但很重要的细节：模态必须在“下一步状态空间”里正确定标。初始系数 b 由最小二乘给出：

b = Φ^+ x_0

于是重建任意时刻状态：

x_k = Σ_{j=1}^{r} b_j Φ_j λ_j^k = Φ(b ∘ λ^k)

其中 ∘ 表示逐元素乘法。这就是谱域时间传播公式。

4.4 EDMD：把可观测量选择权还给使用者

DMD 默认只把状态向量的坐标函数作为可观测量，这对强非线性系统来说限制很大。EDMD 的思想是引入一个字典可观测量：

θ(x) = [θ_1(x), θ_2(x), ..., θ_p(x)]^T

例如包含 x、x²、x³、x₁x₂ 等多项式。然后构造两个矩阵：

Θ(X) = [θ(x_0), ..., θ(x_{m-1})]，Θ(Y) = [θ(x_1), ..., θ(x_m)]

Koopman 矩阵的估计是：

K = G^+ A

其中 G = Θ(X)Θ(X)^T，A = Θ(Y)Θ(X)^T。理论上，字典越完备，K 的特征值越接近真实 Koopman 谱。但字典也不是越大越好：G 矩阵会越来越病态，噪声会被放大，最终过拟合训练数据。所以 EDMD 实际用起来要加正则化，或者用正交化做降维。后面 Matlab 部分我给一个带 Tikhonov 正则化的最小实现。

5. Matlab 实现：从系统仿真到谱绘制的完整流程

5.1 DMD 核心函数

先给出一个可直接复用的核心函数。它实现了上一节描述的流程，并且支持中心化处理。

function [Phi, lam, b, Xrec] = koopman_dmd(X, Y, r, center) % KOOPMAN_DMD 动态模态分解/库普曼算子谱估计 % 输入: % X, Y : n x m 快照对，Y(:,k) 是 X(:,k) 的下一时刻状态 % r : 截断秩 % center: 逻辑值，是否先去掉时间均值 % 输出: % Phi : n x (r+1) 模态矩阵，center=true 时会补一列均值模态 % lam : (r+1) x 1 Koopman 特征值 % b : 模态系数 % Xrec : 用 DMD 重建的状态序列 [n, m] = size(X); if nargin < 3 || isempty(r) r = min(n, m-1); end if nargin < 4 || isempty(center) center = true; end if center mu = mean(X, 2); Xc = X - mu; Yc = Y - mu; else mu = zeros(n, 1); Xc = X; Yc = Y; end % 截断 SVD [U, S, V] = svd(Xc, 'econ'); r = min(r, length(S)); U = U(:, 1:r); S = S(1:r, 1:r); V = V(:, 1:r); % 投影矩阵 Atilde = U' * Yc * V / S; % 特征分解 [W, D] = eig(Atilde); lam = diag(D); % 模态回映 Phi = Yc * V / S * W; % 中心化时补上均值模态，对应 λ=1 if center Phi = [mu, Phi]; lam = [1; lam]; end % 初始系数 b = Phi \ X(:, 1); % 重建 k = 0:m; Xrec = Phi * (b .* (lam .^ k)); end

几个实现细节需要解释清楚。第一，Phi = Yc * V / S * W中用的是 Yc 而不是 Xc，这是 DMD 的标准做法，能让模态在时间传播中保持正确的相位信息。第二，中心化之后补一个均值模态 λ=1，这一步在流场 DMD 里几乎是必须的，因为它把“平均流”和“脉动模态”分开。否则你会发现一个特征值 λ≈1 的大能量模态吸收掉所有均值分量，其他物理模态反而被压缩到噪声级别。第三，求初始系数直接用了左除\，Matlab 会自动选最小二乘解，比手写pinv(Phi) * X(:,1)更稳定也更快。

5.2 主程序：Lorenz 系统的谱估计

下面这段代码用 Lorenz 系统生成数据，调用上面的 DMD 函数，并绘制特征值分布。

%% koopman_lorenz.m clear; close all; clc; % Lorenz 系统参数 sigma = 10; rho = 28; beta = 8/3; f = @(t, x) [sigma*(x(2)-x(1)); x(1)*(rho-x(3)) - x(2); x(1)*x(2) - beta*x(3)]; dt = 0.02; t = (0:dt:30)'; [~, X] = ode45(f, t, [8; 8; 27]); X = X.'; % 3 x N X1 = X(:, 1:end-1); X2 = X(:, 2:end); r = 20; [Phi, lam, b, Xrec] = koopman_dmd(X1, X2, r, true); % 绘制离散特征值 thetaL = linspace(0, 2*pi, 200); figure; plot(cos(thetaL), sin(thetaL), 'k--', 'LineWidth', 1); hold on; plot(real(lam), imag(lam), 'ro', 'MarkerFaceColor', 'r'); axis equal; grid on; xlim([-1.2, 1.2]); ylim([-1.2, 1.2]); xlabel('Re(\lambda)'); ylabel('Im(\lambda)'); title('Lorenz 系统 DMD 特征值'); % 连续时间频率-增长率图 lam_c = log(lam) / dt; figure; plot(imag(lam_c)/(2*pi), real(lam_c), 'b*'); xlabel('频率 (Hz)'); ylabel('增长率'); grid on; title('Lorenz 系统连续时间 Koopman 谱');

Lorenz 系统在 rho=28 时处于混沌状态。我们预期特征值不是孤立的几个点，而是形成一团不规则的分布，这对应连续谱的数值近似。你会看到大部分特征值靠近单位圆但并非精确落在圆上，体现了混沌系统对初始条件的敏感性和有限时间数据的影响。

5.3 旋转系统：验证代码是否正确

验证代码是否正确的最好方法，是先跑一个理论谱已知的系统。单位圆旋转系统就是这样：

%% koopman_rotation.m alpha = 2*pi * 0.17; k = 0:999; Z = exp(1i * alpha * k); Xrot = [real(Z); imag(Z)]; X1 = Xrot(:, 1:end-1); X2 = Xrot(:, 2:end); [PhiR, lamR, bR, ~] = koopman_dmd(X1, X2, 4, true); figure; plot(cos(thetaL), sin(thetaL), 'k--', 'LineWidth', 1); hold on; plot(real(lamR), imag(lamR), 'bs', 'MarkerFaceColor', 'b'); axis equal; grid on; xlabel('Re(\lambda)'); ylabel('Im(\lambda)'); title('旋转系统 DMD 特征值（应落在单位圆上）');

注意这里把复信号 Z 拆成实部和虚部两个分量，构成二维状态。系统本质上是单频旋转，因此 Koopman 谱应该出现一对共轭特征值 e^{±iα}，再加上均值模态 λ=1。如果代码正确，这三个特征值会清清楚楚落在单位圆上。

5.4 混沌 Logistic 映射的延迟嵌入谱

Logistic 映射本身是标量系统，只有一个状态变量，直接用 DMD 会显得单薄。一个更合理的做法是把它嵌入到延迟坐标空间里，用 Hankel 矩阵作为快照矩阵。延迟嵌入在遍历论里可以理解为用历史信息扩充状态空间，让高维 Koopman 近似变得可能。

%% koopman_logistic_delay.m N = 20000; x = zeros(1, N); x(1) = 0.2; for n = 1:N-1 x(n+1) = 4*x(n)*(1-x(n)); end d = 40; % 延迟深度 H = zeros(d, N-d+1); for i = 1:d H(i, :) = x(i:N-d+i); end X1 = H(:, 1:end-1); X2 = H(:, 2:end); [PhiH, lamH, bH, ~] = koopman_dmd(X1, X2, 20, true); figure; plot(cos(thetaL), sin(thetaL), 'k--', 'LineWidth', 1); hold on; plot(real(lamH), imag(lamH), 'g.', 'MarkerSize', 10); axis equal; grid on; xlabel('Re(\lambda)'); ylabel('Im(\lambda)'); title('Logistic 映射延迟嵌入后的 DMD 特征值');

这个例子最能说明“连续谱”是什么意思。Logistic 映射是混合系统，理论上除了常值特征函数，L² 里没有孤立特征函数。延迟嵌入 DMD 得到的特征值会散布在单位圆附近，形成云状结构。这不是数值错误，而是用有限维矩阵逼近连续谱的必然现象。

5.5 EDMD 的简化实现

如果要对强非线性系统做更精细的谱估计，可以用 EDMD。这里给一个带正则化的简化版本：

function [lam, PhiVal] = koopman_edmd(X, Y, dictFun, reg) % KOOPMAN_EDMD 扩展 DMD % dictFun 把状态向量映射为字典可观测向量 % reg 为正则化系数 ThX = dictFun(X); % p x m ThY = dictFun(Y); % p x m G = ThX * ThX'; A = ThY * ThX'; if nargin < 4 || isempty(reg) reg = 1e-8; end K = (A / (G + reg * eye(size(G)))).'; [W, D] = eig(K); lam = diag(D); % 特征函数在采样点上的取值 PhiVal = W.' * ThX; end

一个简单的字典可以用多项式：

dict = @(X) [X; X(1,:).^2; X(2,:).^2; X(1,:).*X(2,:); X(1,:).^3];

EDMD 和 DMD 的核心差别在于：DMD 的字典就是坐标函数本身，EDMD 可以自由选择任意非线性函数作为字典项。字典选择需要针对系统特点来定，没有万能答案。太小的字典近似能力差，太大的字典引入数值病态。实际项目里我会先用延迟嵌入 DMD 把谱的基本结构看清楚，再根据谱的形状决定要不要上 EDMD 和上什么字典。

5.6 运行环境与工具箱说明

上面所有代码只用到了 Matlab 最基础的矩阵计算和绘图命令，不需要任何额外工具箱。版本方面，R2021a 以后都可以直接运行；如果用的是 Octave，把eig(A)的结果手动用diag(D)提出来即可，其余语法完全兼容。装环境不需要折腾，核心就一句话：只要是能跑矩阵分解的 Matlab 环境就够用了。

6. 数值实验：周期、耗散与混沌三大典型系统的谱对比

6.1 实验一：单位圆旋转——纯点谱

用 5.3 节代码跑完，特征值图里应该出现三个点：λ=1 是均值模态，一对共轭特征值位于单位圆上，角度对应 α=0.17×2π。这说明系统只有单频振荡，谱特性极其干净。用 DMD 重建未来 100 步，你会发现相位继续保持，没有振幅衰减。这类系统最简单的谱判据是：有效特征值全部落在单位圆上，且数量级很小。

6.2 实验二：Lorenz 混沌——连续谱近似

Lorenz 实验的特征值云通常呈现两个特点。第一，大部分特征值模长小于 1，但靠近单位圆的数据点明显多于远离单位圆的点；第二，特征值分布存在若干“聚集区”，对应 Lorenz 系统在混沌吸引子上翻转时的优势频率。需要注意的是，这些优势频率并不是调和周期，而是连续谱中能量比较集中的区间。

如果用重建轨迹和真实轨迹对比，前几个主模态就能在短时间内（通常几个 Lyapunov 时间）跟踪得很准，再往后误差迅速放大。这个现象本身就是混沌系统谱特性的一部分：连续谱意味着没有有限个离散模式能长期支撑预测。想靠增加 r 来延长预测时长是不现实的，加再多模态也只是在拟合噪声和轨迹细节。

6.3 实验三：Duffing 耗散——单位圆内的衰减谱

再补一个强耗散系统。无强迫的 Duffing 振子：

%% koopman_duffing.m f = @(t, x) [x(2); -0.3*x(2) - x(1) + x(1)^3]; dt = 0.05; t = (0:dt:25)'; [~, X] = ode45(f, t, [2.5; 0]); X = X.'; X1 = X(:, 1:end-1); X2 = X(:, 2:end); [PhiD, lamD, bD, ~] = koopman_dmd(X1, X2, 8, true); figure; plot(cos(thetaL), sin(thetaL), 'k--', 'LineWidth', 1); hold on; plot(real(lamD), imag(lamD), 'm^', 'MarkerFaceColor', 'm'); axis equal; grid on; xlabel('Re(\lambda)'); ylabel('Im(\lambda)'); title('Duffing 振子 DMD 特征值');

Duffing 振子在初值 2.5 时最终收敛到一个稳恒状态，瞬态过程很快。DMD 特征值会明显落在单位圆内，模长越大表示衰减越慢。这种系统用谱去拟合瞬态特别合适，几个衰减模态就能把整个收敛过程描述得相当精确。

6.4 三类系统谱特征总结

这三类对比放在一起，就能很清楚看出 Koopman 谱特性到底在区分什么：它区分的是系统长期行为中“可分解为有限个稳定模式”的程度。周期系统最可分解，混沌系统最不可分解，耗散系统夹在中间。

7. 实战中的坑与参数选择建议

7.1 先判断谱类型，再谈论模态数

我见过不少项目，拿到数据上来就自动保留 99% 能量，然后开始数 DMD 模态。这顺序反了。正确做法是先看谱的整体形态：如果特征值集中在少量孤立点，说明系统谱接近点谱，模态数取前几个大的就行；如果特征值分布呈云状，说明系统有连续谱成分，模态数更多是控制近似误差，而不是在找“真实物理模式”。连续谱系统里，刻意保留 99% 能量可能会把大量噪声模态也保留进来。

7.2 中心化不是细枝末节

前面提过均值模态的问题，这里再强调一次。对于非零均值的数据，如果不中心化，DMD 的第一模态几乎必然是一个 λ≈1 的“准静态”模态，它把时间平均分量全吃掉了。这样做本身不算错，但它会掩盖真正的脉动结构。我的习惯是标配 center=true，除非你明确知道均值项就是你要研究的对象。中心化之后的均值模态单独放出来，既能清晰看到平均流，又不干扰其他模态。

7.3 连续谱的结果不要过度解读

混沌系统的特征值云，容易让人产生两种极端反应。一种是把每个小圆点都当成一个真实的物理振荡，试图逐个解释；另一种是看到云就认为 DMD 失效了。两种都不对。连续谱系统的 DMD 谱，本质上是有限维投影对无限维连续谱的离散采样。它只能告诉你谱带的位置和大致形状，不能告诉你每一个“伪点”的物理意义。如果项目要求给出准确的连续谱密度，应该配合更强的统计手段，比如矩方法或加长 Hankel 窗口，而不是单纯增加 DMD 的秩。

7.4 时间步长与特征值混叠

连续时间系统的 DMD 特征值 λ 是按照离散时间步长采样的。λ 的角度对应频率的方式是：

ω = imag(log(λ)) / Δt

因此 Δt 过大时，高频成分会折叠到低频区间，产生混叠。经典奈奎斯特准则在这里仍然适用：Δt 必须小于最小关注周期的一半。做流场数据时尤其要注意，高速摄影下采样率如果偏慢，高斯特罗哈尔数的涡脱落模态会被脸部映射到其他频率上。

7.5 噪声与奇异值截断的折中

噪声对 DMD 特征值的影响很微妙。轻微噪声会把单位圆上的特征值往圆内拉，表现为额外的稳定阻尼；强噪声则可能让特征值向外散开，产生虚假不稳定模态。SVD 截断能压住一部分噪声，但截断太狠会把弱模态也丢掉。我的经验是：先画奇异值谱，找明显的“肘部”；如果信号接近混沌，奇异值衰减很慢，此时宁可保留多些模态，再用“特征值落在圆盘内且成对出现”的规律去剔除明显伪迹。

7.6 EDMD 字典不是越大越好

EDMD 最诱人的地方是可以用任意非线性字典，但越丰富的字典越容易把数据中的噪声也建模进来。G 矩阵的条件数会迅速增大，最终导致特征值发生虚假分裂。我在代码里加了一个最小的 Tikhonov 正则化项 reg=1e-8，这只能应付轻微病态；更复杂的场景建议先对字典做正交化，或者用稀疏回归的方式自动筛字典项。EDMD 的正确打开方式是“小步快跑”：从低阶字典开始，逐步加项，观察谱是否稳定收敛。如果加一项字典谱就剧烈变化，说明当前数据支撑不起这么复杂的字典。

收尾：我实际做谱分析时的一个习惯

最后分享一个我自己的操作习惯。每次拿到新的动力系统数据，我不是直接跑完整谱分析，而是先做三件事：画一条轨迹判断系统是收敛、周期还是混沌；用快速傅里叶变换看频谱的能量集中区间；再用 DMD 算一次低秩谱看特征值的整体形态。这三步加起来不超过半小时，但能避免后面 90% 的误判。尤其是当你分不清系统是“带噪声的周期”还是“真的混沌”时，Koopman 特征值在单位圆上的分布形态比单纯看轨迹图要可靠得多。这个方法我用了挺久，希望对你也有用。