lu.cpp

// ////////////////////////////////////////////////////////////////
// LU分解を用いて連立一次方程式を解く
// ////////////////////////////////////////////////////////////////
#include <iostream>

#include <tnt_array1d.h>
#include <tnt_array1d_utils.h>
#include <tnt_array2d.h>
#include <tnt_array2d_utils.h>
#include <jama_lu.h>

int
main(int, char *[])
{
  TNT::Array2D<double> A(5, 5);
  A[0][0] = 7.0; A[0][1] = 4.0; A[0][2] = 1.0; A[0][3] = 8.0; A[0][4] = 3.0;
  A[1][0] = 1.0; A[1][1] = 8.0; A[1][2] = 6.0; A[1][3] = 5.0; A[1][4] = 8.0;
  A[2][0] = 4.0; A[2][1] = 2.0; A[2][2] = 7.0; A[2][3] = 0.0; A[2][4] = 1.0;
  A[3][0] = 3.0; A[3][1] = 7.0; A[3][2] = 5.0; A[3][3] = 1.0; A[3][4] = 1.0;
  A[4][0] = 2.0; A[4][1] = 5.0; A[4][2] = 0.0; A[4][3] = 9.0; A[4][4] = 4.0;
  std::cout << "A: " << std::endl;
  std::cout << A << std::endl;

  TNT::Array1D<double> b(5);
  b[0] = 3.0;
  b[1] = 1.0;
  b[2] = 8.0;
  b[3] = 4.0;
  b[4] = 9.0;
  std::cout << "b: " << std::endl;
  std::cout << b << std::endl;

  JAMA::LU<double> lu(A);
  std::cout << "solve: " << std::endl;
  std::cout << lu.solve(b) << std::endl;

  return 0;
}