使用同伦-切比雪夫方法求解参数特征值问题

《Mathematics and Computers in Simulation》：Solving parametric eigenvalue problems using a Homotopy-Chebyshev method

　　吴云云 | 李雅云
南京财经大学数学系，中国南京文苑路8号，210046

**摘要**
本文提出了一种用于求解参数化特征值问题的鲁棒数值方法，其中系统矩阵随参数平滑变化，其特征值和特征向量也随之变化。虽然基于泰勒或切比雪夫展开的现有方法能够提供精确的局部近似，但它们往往难以处理大参数变化、特征值交叉和接近退化的情况。为了解决这些问题，我们提出了一种同伦-切比雪夫方法，该方法沿着连续的参数路径追踪特征值对解，从而在广泛的参数范围内提高鲁棒性。该方法结合了切比雪夫多项式的全局逼近能力和同伦延续的稳定性，克服了纯展开方法的固有缺陷。数值实验表明，我们的同伦方法在广泛的参数区间内优于纯展开方法，并且对于中等规模的系统仍然具有较高的计算可行性。在结构动力学等实际应用中，该方法能够准确预测不同材料属性下的振动频率，显著提高了可靠性。该算法的鲁棒性使其特别适用于不确定性量化、结构动力学等参数依赖性起关键作用的应用领域。

**引言**
参数化特征值问题（PEPs）的特点是方程 A(μ)v(μ) = λ(μ)v(μ)，其中 A(μ)∈?^n×n 取决于参数 μ。这类问题在计算物理、结构工程和数据驱动建模中无处不在。在本研究中，我们假设矩阵元素 A(μ)ij 是 μ 的足够平滑的函数。如果所有元素 A(μ)ij 关于 μ 是连续的，那么 n 个特征值和特征向量函数 λ(μ) 和 v(μ) 也是连续的。此外，如果 A(μ) 由 μ 的全纯函数组成，则特征值 λ(μ) 是解析函数，除了在孤立的代数奇点处外。
在从材料变化引起的振动分析到参数依赖图谱聚类的多种应用中，都需要高效地追踪特征值对轨迹 (λ(μ), v(μ))。传统的数值方法面临三个基本挑战：（1）在大的参数区间内保持精度；（2）解决特征值交叉和分岔；（3）扩展到高维系统。
最近在微分方程数值方法方面的进展促使人们开发出用于参数化特征值问题的鲁棒求解器。例如，分数阶微积分和基于同伦的方法已成功应用于复杂动态系统 [1]、[2]、[3]。这些研究强调了对于能够处理参数依赖奇点和非线性变化的可靠数值技术的需求。在此背景下，我们将同伦分析框架扩展到参数化特征值问题，提供了一种结合局部精度和全局收敛性的统一方法。

解决 PEPs 的传统方法大体上可以分为局部展开方法和全局逼近技术。局部方法（如泰勒级数展开 [4]）在固定参数值附近提供高效的近似，但在较大的参数变化下精度会下降。这些方法通常采用张量化Arnoldi方法，通过形式为 A(μ) = ∑_?=0^M A_? φ_?(μ) 的多项式混沌展开来计算极端特征值，从而在参数空间中对问题进行离散化。
使用多项式或有理逼近的全局方法 [5]、[6] 在参数范围内提供更好的覆盖，但在处理奇点和大型系统时计算成本较高。这些方法采用类似的离散化策略，但通常使用牛顿-拉弗森迭代来解决由此产生的张量方程。Sousedik 和 Elman [7] 的最新工作展示了基于张量方法的改进性能，而 Williams [8] 则开发了专门用于随机特征值问题的技术。Rahman [9] 计算了广义随机特征值问题 A(μ)v(μ) = λ(μ)B(μ)v(μ) 的统计矩。Verhoosel、Gutiérrez 和 Hulshoff [10] 假设 A(μ) = A_0 + ∑_i=1^M A_i μ_i，其中 μ 是参数向量，采用有限元离散化方法处理 λ(μ) 和 v(μ)。这种离散化方法通过逆迭代求解张量结构方程。Williams [8]、[11] 采取了不同的方法，引入了人为的时间依赖性将问题转化为积分方程。Williams 还使用了多项式混沌展开，但适用于特征向量而非 A(μ)。
最近，Mach 和 Freita [12] 提出了两种求解参数化特征值问题的数值方法。第一种方法使用关于参考点 μ0 的泰勒级数展开，通过迭代计算泰勒系数来逼近特征值对。这种方法在 μ0 附近提供局部精确的近似。第二种技术利用切比雪夫展开在更宽的参数区间内构建特征值对的多项式近似。与局部化的泰勒级数方法不同，切比雪夫展开在扩展范围内提供了更高的精度，同时保持计算效率。然而，每次迭代都需要求解线性系统，计算负担仍然较大。为了解决这些挑战，我们开发了一种将延续技术与切比雪夫逼近相结合的同伦-切比雪夫方法。该方法自适应地推进参数，并沿着延续路径构建局部谱表示，从而在更宽的参数范围内稳定追踪特征值对，而无需依赖纯粹的局部展开。

所提出的同伦-切比雪夫方法通过几个关键创新点与现有方法区别开来。虽然纯切比雪夫方法能够提供全局逼近，但它们难以处理非解析点和特征值交叉。我们的方法将切比雪夫逼近整合到同伦延续框架中，实现全局多项式逼近与局部路径跟踪校正之间的自动切换。这种自适应策略克服了基于泰勒的方法（收敛半径有限）和纯切比雪夫方法（难以处理奇点）的局限。与传统遵循预定路径的延续方法不同，我们的算法根据局部误差估计动态地调整逼近和校正阶段。
从计算角度来看，我们的同伦方法通过几种创新机制实现了显著的效率提升。问题自然地分离为预测器和校正器步骤，使得通过智能路径规划可以最小化昂贵的矩阵分解。此外，该方法的结构允许与现有技术无缝集成：局部展开在规则区域内提供高精度，而全局多项式逼近在较大区间内保持分辨率，并通过严格的误差估计器实现自动切换。

**数值实现**
数值实验证明了该方法在多个关键方面的优越性：它能够在传统方法失效的非解析点追踪特征值轨迹，有效处理多个特征值对，以及在大规模问题上的计算可扩展性。即使在涉及非正规矩阵或参数依赖谱结构的挑战性场景下，该方法也能保持稳定的性能，使其特别适用于传统扰动方法不可靠的问题。
该方法主要用于离线预计算参数范围内的特征值对轨迹，从而可以通过多项式评估实现快速实时查询。然而，该算法也可以通过增量更新支持对缓慢变化参数的实时分析。这种灵活性使其适用于从设计优化（离线）到结构系统在线监控的各种应用。

**问题陈述**
给定一个平滑参数化的矩阵 A(μ)∈?^n×n，我们希望计算满足 A(μ)v(μ) = λ(μ)v(μ) 且 ‖v(μ)‖^2 = 1 的相应特征值对 (λ(μ), v(μ))。假设 A(μ) 关于 μ 是可微的，并且已知在 μ=μ0 处的初始特征值对 (λ_0, v_0)，我们的目标是追踪 μ 变化时 λ(μ) 和 v(μ) 的演变。我们引入了零阶变形方程 H(v, λ, t) = (1?t)A_μ_0v ? λv + t(A(μ)v ? λv) = 0，其中 t∈[0,1] 是同伦参数。

**数值策略**
预测器-校正器方法为追踪解路径提供了一个有效的框架，结合了欧拉预测和牛顿校正。预测器步骤通过求解包含雅可比矩阵的边界系统来计算切线方向，并加入归一化约束，得到导数 dv/dt 和 dλ/dt：
dv/dt = (1?tk)A_μ_0 + tkA(μ)?λ_kI?v_kvk
H_0^?1?A(μ)?A(μ_0)vk_0。
然后通过显式欧拉积分更新预测值：
vk+1 = vk + Δt·dv/dt, λ_k+1 = λ_k + Δt·dλ/dt, tk+1 = tk

**收敛性分析**
Bernstein 椭圆 Eρ 是复平面上的一个区域，定义为圆 |z| = ρ > 1 在映射 x = (z+z?1)/2 下的像。这种构造在逼近理论中广泛用于分析切比雪夫展开的收敛性。参数 ρ 控制椭圆的大小，并与展开系数的指数衰减速率直接相关。基于这种几何设置，我们建立了关于解析特征值切比雪夫展开系数衰减的以下引理。

**数值实验和结果**
本节验证了所提出的同伦-切比雪夫算法在求解参数化特征值问题方面的有效性和鲁棒性。实验通过与参考解的比较来验证精度。

**结论**
本文提出了一种同伦-切比雪夫方法，用于求解参数化特征值问题，解决了大参数变化、特征值交叉和接近退化情况等挑战。该方法将切比雪夫全局逼近与同伦延续相结合，确保在广泛的参数范围内获得鲁棒和准确的解。与纯展开方法不同，我们的方法根据误差估计动态地在全局逼近和局部校正之间切换。

**作者贡献声明**
吴云云：撰写 – 审稿与编辑、原始稿撰写、方法论、形式分析。
李雅云：撰写 – 审稿与编辑、可视化、软件开发、调查。