べき乗法

Power_Method.cpp

//
//  Power_Method.cpp
//  Power Iteration
//	べき乗法
//  Created by knuu on 2014/06/23.
//
//

#include <iostream>
#include <iomanip> // std::setpricision
#include <Eigen/Dense> // 行列計算のライブラリ
#include <cmath> // pow

#define M_PRE 10 // 小数点以下の桁数
#define SIZE_A 2 // 行列のサイズ
#define EPS 0.000000001 // 誤差

// べき乗法
void power_method(Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> A, int count)
{
  std::cout<<std::fixed<<std::setprecision(M_PRE); // 小数点以下の桁数設定
	
  Eigen::Matrix<double, Eigen::Dynamic, 1> X_n; // X^(k+1)
  X_n = Eigen::MatrixXd::Ones(A.cols(), 1);	// 初期値はX^(0) = [1 1　... 1]^T
  Eigen::Matrix<double, Eigen::Dynamic, 1> X = X_n; // X^(k)
	
  // count回繰り返しかける
  for (int i = 0; i < count; i++) {
    X = X_n;
    X_n = A * X;	// X^(k+1) = A * X^(k)
  }
	
  double lambda = X_n(0,0) / X(0,0);
  for (int i = 0; i < A.cols(); i++) {
    X_n(i, 0) = X_n(i, 0) / pow(lambda, count);
  }
  std::cout<<"最大固有値: "<<lambda<<std::endl;
  std::cout<<"最大固有値の固有ベクトル"<<std::endl;
  std::cout<<X_n<<std::endl;
}

int main()
{
  Eigen::Matrix<double, SIZE_A, SIZE_A> A;
  A << 1,2, 1,1;
  int count = 500; // 繰り返し回数
  power_method(A, count);
}