はじめに

べき乗法(Power Method)は、数値解析において行列の最大固有値と対応する固有ベクトルを効率的に求める手法です。特性方程式を解くことで、固有値を理論的に検証する方法も紹介します。今回はPythonを用いて具体的なコードとその解説を行います。

プログラム概要

以下の内容をPythonで実装します。

べき乗法を用いて行列の最大固有値と固有ベクトルを求める
特性方程式を解いて理論的に固有値を確認する

べき乗法の実装

Power_Method.py

import numpy as np

def power_method(A, x0, epsilon=1e-1, max_iterations=3):
    """
    べき乗法を用いて行列の最大固有値と固有ベクトルを求めます。

    Args:
        A: 行列。
        x0: 初期ベクトル。
        epsilon: 許容誤差。
        max_iterations: 最大反復回数。

    Returns:
        最大固有値と固有ベクトルのリスト。
    """
    x = x0
    lambda_old = 0
    k = 0
    eigenvalues = []  # 固有値のリスト
    eigenvectors = []  # 固有ベクトルのリスト

    while k < max_iterations:
        y = A @ x
        lambda_new = y[1]  # 固有値の推定
        x = y / y[1]  # 固有ベクトルを正規化

        eigenvalues.append(lambda_new)
        eigenvectors.append(x)

        if abs(lambda_new - lambda_old) <= epsilon:
            break
        lambda_old = lambda_new
        k += 1

    return eigenvalues, eigenvectors

# 使用例
A = np.array([[9, 4], [4, 3]])  # 行列 A
x0 = np.array([1, 1])  # 初期ベクトル

eigenvalues, eigenvectors = power_method(A, x0)

# 計算結果の表示
x1_1 = round(eigenvectors[0][0], 3)
x1_2 = round(eigenvectors[1][0], 3)
x1_3 = round(eigenvectors[2][0], 3)
lambda_1_3 = round(eigenvalues[2], 3)

print(f"(1) x1(1) = {x1_1}")
print(f"(2) x1(2) = {x1_2}")
print(f"(3) x1(3) = {x1_3}")
print(f"(4) lambda1(3) = {lambda_1_3}")

実行結果

(1) x1(1) = 1.857  
(2) x1(2) = 1.986  
(3) x1(3) = 1.999  
(4) lambda1(3) = 10.945

特性方程式で固有値を検証

べき乗法で得られた最大固有値が正しいか、特性方程式を解いて確認します。
コードはこちら↓

characteristic_equation.py

import numpy as np

# 行列 A を定義
A = np.array([[9, 4], [4, 3]])

# 特性方程式を解いて固有値を求める
eigenvalues = np.linalg.eigvals(A)

# 最大固有値を求める
max_eigenvalue = max(eigenvalues)

# 最大固有値に対応する固有ベクトルの1番目の要素を求める（𝑥2=1と仮定）
x2 = 1
x1 = (max_eigenvalue - 3) * x2 / 4

# 結果を出力
print(f"(5) |𝜆1| = {round(max_eigenvalue, 3)}")
print(f"(6) x1 = {round(x1, 3)}")

実行結果

(5) |𝜆1| = 11.0  
(6) x1 = 2.0

まとめ

べき乗法は行列の最大固有値と対応する固有ベクトルを反復計算で求める方法
特性方程式を解くことで固有値を正確に求めること

このプログラムは数値解析の基本的な学習に役立ちます。応用として、より高次元の行列や他の固有値計算法を実装することも可能です。最後まで読んでくださり、ありがとうございました。もし改善点や質問があれば、ぜひコメントしてください！