[论文解读] On Computation of Matrix Mittag-Leffler Function

Ivan Matychyn|arXiv (Cornell University)|Jun 5, 2017

Matrix Theory and Algorithms参考文献 10被引用 4

一句话总结

本文提出了一种基于若尔当代数标准型分解的数值方法，用于计算矩阵梅特莱夫函数，该方法已实现为 MATLAB 例程（mlfm.m）。该方法能够精确且高效地计算求解分数阶微分方程所必需的矩阵梅特莱夫函数，验证结果表明与解析解相比绝对误差低于 10⁻¹⁵。

ABSTRACT

A method for computation of the matrix Mittag-Leffler function is presented. The method is based on Jordan canonical form and implemented as a Matlab routine.

研究动机与目标

• 开发一种可靠且高效的矩阵梅特莱夫函数计算方法，该函数是求解分数阶微分方程（FDEs）的关键组成部分。
• 克服现有方法的局限性，例如缺乏半群性质，导致无法直接应用适用于矩阵指数的缩放与平方法。
• 提供一种基于谱分解的数值稳定算法，特别适用于具有分数阶动力学特性的系统。
• 通过基准问题对方法进行实现与验证，包括具有解析参考解的巴格利-托尔维克方程。
• 提供一个公开可用的 MATLAB 例程（mlfm.m），以供科学计算与分数阶动力学研究中的实际应用。

提出的方法

• 该方法将矩阵 A 分解为其若尔当代数标准型 J = P⁻¹AP，从而通过 A 的谱集上的插值计算 f(A)。
• 对于与特征值 λₖ 关联的每个若尔当块，利用标量梅特莱夫函数 Eα,β(λₖtα) 的导数（最高至 mk−1 阶）来计算矩阵函数。
• 标量广义梅特莱夫函数 Eα,β+(m−1)αᵐ(λₖ) 通过最优抛物线路径（OPC）算法计算，该算法由 R. 加拉帕实现为 MATLAB 例程。
• 该算法同时处理求解卡普托型 FDEs 所需的 Eα,₁ 和 Eα,α 函数，例如 z(t) = Eα,₁(Atα)z⁰ + ∫₀ᵗ(t−τ)α⁻¹Eα,α(A(t−τ)α)f(τ)dτ。
• 该方法应用于 c=0 时的巴格利-托尔维克方程，利用误差函数与互补误差函数的性质，推导出 E½,₁(B) 和 E½,½(B) 的解析表达式。
• 通过将计算得到的矩阵梅特莱夫函数与涉及 erfc(−p) 和指数积分的闭式表达式进行比较，对实现进行验证，误差分析证实了其高精度。

实验结果

研究问题

• RQ1 当标准矩阵指数技术（如缩放与平方法）不适用时，如何准确且高效地计算矩阵梅特莱夫函数？
• RQ2 是否可以有效利用若尔当代数标准型来高精度计算涉及广义梅特莱夫函数的矩阵函数？
• RQ3 基于 MATLAB 的矩阵梅特莱夫函数实现，在具有已知解析解的基准问题上可达到何种数值精度？
• RQ4 该方法在不同矩阵结构（如高阶分数阶微分方程中出现的矩阵）下，其稳定性和精度在多大程度上得以保持？
• RQ5 在矩阵函数计算的背景下，与替代数值方法相比，基于 OPC 的标量函数评估性能如何？

主要发现

• mlfm.m MATLAB 例程在与 c=0 时巴格利-托尔维克方程的解析解比较时，计算矩阵梅特莱夫函数的绝对误差低于 10⁻¹⁵。
• 利用互补误差函数的性质和恒等式 E½,₁(z) = eᶻ² erfc(−z)，推导出 E½,₁(B) 和 E½,½(B) 的解析表达式，从而实现严格验证。
• 该方法成功计算了具有非对角化若尔当块的矩阵梅特莱夫函数，证明了其对病态矩阵的鲁棒性。
• 使用最优抛物线路径（OPC）算法计算标量梅特莱夫函数，确保了所需广义函数的高精度评估。
• 该实现正确处理了矩阵 α-指数函数 eₐᴬᵗ = tᵃ⁻¹Eα,α(Atᵃ)，这是求解黎曼-刘维尔与卡普托形式 FDEs 所必需的。
• 验证结果表明，该方法实现了机器精度的准确性，表现为计算矩阵与参考矩阵之间的逐元素绝对差低于 10⁻¹⁵。