#ifndef __MATRIX_H
#define __MATRIX_H

#include "common.h"
#include "complex.h"

typedef unsigned char error_code_t;

/**
 * 行列を表す構造体(struct)
 * 
 */
typedef struct {	
  double *values;       ///< 行列の成分を保存するメモリへのポインタ
  unsigned int rows;    ///< 行数   
  unsigned int columns; ///< 列数
} matrix_t;

void matrix_init(matrix_t *target);
error_code_t matrix_alloc(matrix_t *target, unsigned int rows, unsigned int columns);
error_code_t matrix_zero(matrix_t *target);
error_code_t matrix_free(matrix_t *target);
inline double *matrix_value(matrix_t *target, unsigned int row, unsigned int column);
error_code_t matrix_print(matrix_t *target);
BOOL matrix_is_square(matrix_t *target);
BOOL matrix_is_symmetric(matrix_t *target);
error_code_t matrix_copy(matrix_t *src, matrix_t *dist);
error_code_t matrix_unit(matrix_t *target, unsigned int size);
error_code_t matrix_transpose(matrix_t *target, matrix_t *result);
error_code_t matrix_swap_rows(matrix_t *target, unsigned int row1, unsigned int row2);
error_code_t matrix_swap_columns(matrix_t *target, unsigned int column1, unsigned int column2);
error_code_t matrix_inverse(matrix_t *target, matrix_t *result);
error_code_t matrix_add(matrix_t *target1, matrix_t *target2, matrix_t *result);
error_code_t matrix_sub(matrix_t *target1, matrix_t *target2, matrix_t *result);
error_code_t matrix_multi(matrix_t *target1, matrix_t *target2, matrix_t *result);
error_code_t matrix_scalar_multi(matrix_t *target, double scalar);
error_code_t matrix_scalar_div(matrix_t *target, double scalar);

/**
 * 行列を表す構造体(struct)
 * こちらは要素が複素数用。
 * 
 */
typedef struct {	
  complex_t *values;      ///< 要素を格納するためのメモリへのポインタ
  unsigned int rows;    ///< 行数
  unsigned int columns; ///< 列数
} cmatrix_t;

void cmatrix_init(cmatrix_t *target);
error_code_t cmatrix_alloc(cmatrix_t *target, unsigned int rows, unsigned int columns);
error_code_t cmatrix_zero(cmatrix_t *target);
error_code_t cmatrix_free(cmatrix_t *target);
inline complex_t *cmatrix_value(cmatrix_t *target, unsigned int row, unsigned int column);

error_code_t cmatrix_print(cmatrix_t *target);

error_code_t matrix_eigen22(
    matrix_t *target, 
    unsigned int row, unsigned int column, 
    complex_t *eigen_values);

error_code_t matrix_pivot_merge(
      matrix_t *target, 
      unsigned int row, unsigned int column, 
      matrix_t *mergee);

error_code_t matrix_hessenberg(
    matrix_t *target, 
    matrix_t *result, 
    matrix_t *transform);

#define EIGEN_VALUE_THRESHOLD 1E-10

error_code_t matrix_eigen(
    matrix_t *target, 
    complex_t *eigen_values[], 
    cmatrix_t *eigen_vectors[]);

#endif /* __MATRIX_H */