(19)国家知识产权局 (12)发明 专利申请 (10)申请公布号 (43)申请公布日 (21)申请 号 202210550149.1 (22)申请日 2022.05.20 (71)申请人 大连理工大 学 地址 116024 辽宁省大连市甘井 子区凌工 路２ 号 (72)发明人 郑勇刚　吴哲同　刘振海　叶宏飞　 (74)专利代理 机构 辽宁鸿文知识产权代理有限 公司 21102 专利代理师 许明章　王海波 (51)Int.Cl. G06F 30/27(2020.01) G06F 30/17(2020.01) G06F 17/16(2006.01) G06F 17/13(2006.01) G06N 3/04(2006.01)G06N 3/08(2006.01) G06F 111/10(2020.01) G06F 119/14(2020.01) G06F 119/08(2020.01) (54)发明名称 数据驱动识别偏微分方程的序列奇异值过 滤方法 (57)摘要 本发明属于计算力学领域， 提供一种数据驱 动识别偏微 分方程的序列奇异值过滤方法， 该方 法使用奇异值分解方法和强秩揭示正交三角分 解方法提取隐藏在数据集中的本征线性结构， 提 高了从观测数据提取非线性偏微分方程算法的 计算效率、 准确度和鲁棒性。 本发明通过设计奇 异值过滤的算法， 对数据矩阵进行反复缩减并计 算矩阵的奇异值， 可以有效识别控制方程的个数 并提取有关项。 采用强秩揭示正交三角分解方法 为每个控制方程设置左端项， 并通过奇异值分解 和最小二乘法有效识别控制方程组的最简形式。 本发明所提出的方法作为一种新的控制方程数 据驱动识别方法， 借助高效的矩阵分解技术， 快 速准确， 并且可以通过编写简单循环程序实现， 降低了数值实施复杂度。 权利要求书4页 说明书10页 附图3页 CN 114970339 A 2022.08.30 CN 114970339 A 1.一种数据驱动识别偏微分方程的序列奇异值过 滤方法， 其特 征在于， 步骤如下： 首先， 给出序列奇异值过滤方法识别偏微分方程组的基本格式； 对于一个在n维的时空 坐标(t x1 x2…xn‑1)下， 由m个 变量(u1 u2…um)组成的系统， 设该系统由一组微分方程控制： 将微分方程中所有可能存在的项罗列出来， 如式(2)所示： 其中， 为待定系数； 设采样点的数量 为M， 式(2)应在所有采样点上处 处成立， 即： 其中， X1～XN代表所有假设存在的项， 包括高阶导数和非线性项， 因此根据在各个采样 点处测得的变量值， 通过数值计算得到数据 矩阵A所有 元素的大小； 序列奇异 值过滤方法的 目标是求解式(3)中的系数矩阵A： 识别控制方程的个 数P， 并计算待定系数 的值， 同时使 每一列系数向量尽可能稀疏； 接着， 使用 奇异值分解方法来确定控制方程的个数， 并过滤掉函数库中冗余的项； 式 (3)表明， 由于控制方程的存在， 使得数据矩阵A分布在低维的子空间上， 而不是整个N维空 间： 对于P个控制方程控制的系统， 控制方程会对数据集产生P维的约束， 使得数据矩阵A分 布在N‑P维的子空间上； 使用奇异值分解方法提取出P的值； 奇异值分解是将M ×N的数据矩 阵A分解为正交矩阵U， V和对角矩阵∑的乘积： A＝U∑VT                        (4) 由于采样点个数M远大于函数库的规模N， 将上式展开有： 其中， σ1， σ2， ...， σN称为A的奇异值， 满足σ1≥σ2≥…≥σN≥0； 当A不是满秩矩 阵时， 即 rank(A)＜N， 此时A的奇异 值只有前rank(A)个非零， 即： σrank(A)+1＝σrank(A)+2＝…＝σN＝0； 此 时得到约化的奇异值分解： 测量数据会受到噪声的影响， 并且数值微分也会带来难以避免的数值误差， 导致实 际 的数据矩阵A并不是完全亏秩的， 是一个近似奇异的病态矩阵； 此时后N ‑rank(A)个奇异值 的大小不严格等于0， 是与前面奇异值相比的小量， 它们的大小关系取决于噪声的量级大权　利　要　求　书 1/4 页 2 CN 114970339 A 2小； 当噪声控制在一定范围时， A的后N ‑rank(A)个奇异值会远小于前面的奇异值， 即 σrank(A)+1， σrank(A)+2， ...， σN＜＜σrank(A)； 基于上述原理， 通过判断矩阵的小奇异值的数量来确 定控制方程的个数， N ‑rank(A)即为P的大小； 此外， 为了防止各项的数量级相差过大对计算 结果造成影响， 在预处 理阶段对数据矩阵A的每一列实施L2范 数的归一 化； 为了提高后续流程的计算效率和精度， 通过设计迭代流程， 反复对数据矩阵进行缩减 并实施奇异值分解， 来筛选所有与控制方程相关的项； 具体过程如下： 依 次删掉数据矩阵A 中的每一列进 行奇异值分解， 删掉某列后如果仍然识别出了P 个控制方程， 说明这一列对应 的项不包含在 控制方程组里面； 如果识别出了P ‑1个控制方程， 说明这项存在于控制方程组 里面； 计算N次后可筛选出N ′个与控制方程有关的项， 将这些项对应的列保留下来得到重构 的数据矩阵A ′， 以供后续 步骤使用； 之后， 使用强秩揭示正交三角分解方法sRRQR提取每个方程的左端项； 给出sRRQR的计 算格式， 对于任一数据矩阵A， sR RQR将其分解成式(7)形式： 其中， Q为正交矩阵， R为上三角矩阵， Q1， Q2， B， C， D分别为 Q， R的子块； 并满足： σ(D)≤σQ+1(A)·q1(Q， N)               (9) |(B‑1C)1～Q， 1～N ‑Q|≤q2(Q， N)                   (10) 其中σ(B)是矩阵B的奇异值， Q是需要预先设定的值， q1(Q， N)， q2(Q， N)具有多项式级的 上界； 式(10)指矩阵B‑1C每个元素的大小都被限制不大于q2(Q， N)； 式(7)中的Π是一个置换 矩阵； 将前面通过奇异值分解确定的矩阵的秩rank(A)代入为Q的值， 此时Q＝rank(A)＝N ′ ‑ P， 所以有N ′ ‑Q＝P， 对缩减后的数据矩阵A ′使用sRRQR方法； 根据计算得到的Π来提取A ′Π 的后N′ ‑Q列在原矩阵A′中的位置信息， 即找出了N ′ ‑Q个即P个左端 项在A′中的位置； 最后， 使用奇异值分解方法和最小二乘法确定稀疏形式的控制方程组； 确定控制方程 稀疏形式的方法与用奇异值分解方法过滤函数库中冗余项的方法相同： 将P个左端项分别 逐个与N′ ‑P个候选右端项组成矩阵， 此时的矩阵受一个方程控制， 即这个左端项对应的方 程； 依次删掉每个候选右端项进行奇异值分解， 删掉某项后如果仍然识别出了1个控制方 程， 说明这一列对应的项不包含在此控制方程里面； 如果识别出了0个控制方程， 说明这项 存在于此控制方程中； 进行P(N ′ ‑P)次奇异值分解后， 所有方程的左端项和右端项都已确 定； 使用最小二乘法计算对应系数， 结果如式(1 1)所示： 根据上述理论推导得出的奇异值分解和强秩揭示正交三角分解方法的基本格式， 嵌入权　利　要　求　书 2/4 页 3 CN 114970339 A 3