从零开始的数值分析--矩阵特征分解

数学不易猫猫叹气

于 2024-10-04 21:42:12 发布

阅读量618

点赞数 19

分类专栏：数值分析文章标签：矩阵 python 线性代数

本文链接：https://blog.csdn.net/qq_65711983/article/details/142707633

版权

# 导入numpy库，用于进行科学计算
import numpy as np

# 导入matplotlib库，用于数据可视化
import matplotlib.pyplot as plt

# 导入scipy库，用于科学计算和工程计算
import scipy as sp

# 导入sympy库，用于符号数学计算
import sympy as sym

## 约定内容

decimalNumber = 6  ##小数位数

！！！观前提醒,本章对于线性代数功底有一定要求,如果有可能还可以适当了解一下高等代数中的线性空间,若尔当标准型等内容。因本文并不涉及复矩阵与复特征值的处理,因此酉空间及其相关知识可以先放一放

对于求解方程组而言,我们的讨论已经足以解决一大部分数值解的问题,而当我们将目光放到系数矩阵的计算时,发现问题变得比较复杂.
Definition特征值与特征向量:
设 $\in \mathbb{R}^{n \times n}$ , $\lambda \in \mathbb{C}$ 和非零向量 $\in \mathbb{R}^n$ 满足 $\lambda x$ ,则 $\lambda$ 称为 $A$ 的特征值,而 $x$ 称为 $A$ 的特征向量.
因此,求解特征值就意味着求解等价问题:
$\lambda I)x = 0$
在结束了特征值和特征向量定义后,我们再来补充一些特征值的性质:
Method1线性运算性:
设 $\lambda$ 为 $\in \mathbb{R}^{n \times n}$ 的特征值,则有:

c $\lambda$ 为c $A$ 的特征值
$\lambda - \mu$ 为 $\mu I)$ 的特征值
$\lambda^k$ 为 $A^k$ 的特征值

Method2充要条件:

矩阵 $A$ 可相似对角化的充要条件是 $A$ 有n个线性无关的特征向量.
矩阵 $A$ 的不同特征值构成的特征子空间正交,即不同特征值对应的特征向量正交.

Method3Rayleigh商估计:

$\lambda_1 \leq \lambda_2 \leq ... \leq \lambda_n$ ,则 $\lambda_1 \leq \frac{(Ax,x)}{(x,x)} \leq \lambda_n$
$\lambda_1 = \min_{x \in \mathbb{R^n} , x \neq 0} \frac{(Ax,x)}{(x,x)}$ , $\lambda_n = \max_{x \in \mathbb{R^n} , x \neq 0} \frac{(Ax,x)}{(x,x)}$

def Rayleigh(A,x):
    """
    瑞利商定义

    Args:
        A (Matrix): 要求对称的系数矩阵
        x (vector): 特征向量

    Returns:
        float:瑞利商的值
    """
    
    b = np.dot(A,x)
    
    return np.dot(x,b)/np.dot(x,x)

A = np.array([
    [2,5,0],
    [5,3,7],
    [0,7,6]
])

lam,x = np.linalg.eig(A)

print(Rayleigh(A,x[2]))

-3.7754908865633916

矩阵幂法求特征值简介

矩阵幂法是一种迭代算法，主要用于求解特定类型矩阵（如幂等矩阵、对称正定矩阵等）的主特征值及其对应的特征向量。该方法尤其适用于当矩阵具有特殊性质时，能够简化计算过程并提高效率。下面简要介绍矩阵幂法的基本原理及步骤。

基本原理

矩阵幂法基于以下观察：若矩阵 $A$ 的最大特征值为 $\lambda_1$ ，对应的特征向量为 $\mathbf{v}_1$ ，那么对于任意非零初始向量 $\mathbf{x}_0$ ，通过重复乘以矩阵 $A$ （即计算 $A^n\mathbf{x}_0$ ），得到的序列往往会趋近于 $\lambda_1^n\mathbf{v}_1$ 的方向。这是因为主特征值对应的特征向量在迭代过程中会被逐渐放大，而其他特征值对应的分量会相对减小。

步骤概述

选择初始向量：选取一个非零的初始向量 $\mathbf{x}_0$ 。
迭代计算：对于迭代次数 (n=1,2,…)，执行以下操作：
计算 $\mathbf{x}_n = A\mathbf{x}_{n-1}$ 。
归一化：为了防止数值溢出，每次迭代后通常需要将 $\mathbf{x}_n$ 归一化，即 $\mathbf{x}_n := \frac{\mathbf{x}_n}{||\mathbf{x}_n||}$ 。
终止条件：设定一个收敛准则，比如当两次迭代间的向量变化小于某个阈值时，停止迭代。
估计特征值：最后得到的向量近似为最大特征值对应的特征向量，可以通过 $\mathbf{x}_n^T A \mathbf{x}_n / (\mathbf{x}_n^T \mathbf{x}_n)$ 来估算主特征值 $\lambda_1$ 。

注意事项

矩阵幂法对初值敏感，不同的初始向量可能收敛到不同的特征值。
对于非幂等或非对称矩阵，此方法可能不会直接给出所有特征值，特别是对于复数特征值或接近的特征值情况。
为了提高收敛速度和稳定性，可以采用预处理技术，如对矩阵进行对角化近似或使用加速技术（如幂迭代的改进版，如反幂迭代、共轭梯度法等）。

def PowerMethod(A,v0 = None,isRayleigh = False,tol = 1e-6,max_iter = 150,output_key = False):
    """幂法求解矩阵特征值

    Args:
        A (Matrix): 系数矩阵
        v0 (vector, optional): 迭代初始变量, Defaults to None. If None,随机初始化.
        isRayleigh (bool, optional): 采用瑞利商计算最大值. Defaults to False.
        tol (float, optional): 精度限. Defaults to 1e-6.
        max_iter (int, optional): 最大迭代上限. Defaults to 150.
        putout_key (bool, optional): 输出迭代信息. Defaults to False.
    """
    
    A = A.copy()
    
    if v0 is None:
        v0 = np.random.rand(A.shape[0])
    else:
        v0 = v0.copy()
        
        
    u0 = v0
    
    oldmu =