using System;
using System.Collections.Generic;
namespace NDimensionalProjection
{
public class Matrix
{
private readonly double[,] _data;
public int Rows { get; }
public int Cols { get; }
public Matrix(int rows, int cols)
{
Rows = rows;
Cols = cols;
_data = new double[rows, cols];
}
public Matrix(double[,] data)
{
Rows = data.GetLength(0);
Cols = data.GetLength(1);
_data = (double[,])data.Clone();
}
public double this[int r, int c]
{
get { return _data[r, c]; }
set { _data[r, c] = value; }
}
// --- 演算子オーバーロード ---
public static Matrix operator +(Matrix a, Matrix b)
{
if (a.Rows != b.Rows || a.Cols != b.Cols)
throw new ArgumentException("行列のサイズが一致しません。");
var result = new Matrix(a.Rows, a.Cols);
for (int i = 0; i < a.Rows; i++)
for (int j = 0; j < a.Cols; j++)
result[i, j] = a[i, j] + b[i, j];
return result;
}
public static Matrix operator -(Matrix a, Matrix b)
{
if (a.Rows != b.Rows || a.Cols != b.Cols)
throw new ArgumentException("行列のサイズが一致しません。");
var result = new Matrix(a.Rows, a.Cols);
for (int i = 0; i < a.Rows; i++)
for (int j = 0; j < a.Cols; j++)
result[i, j] = a[i, j] - b[i, j];
return result;
}
public static Matrix operator *(Matrix a, Matrix b)
{
if (a.Cols != b.Rows)
throw new ArgumentException("左の行列の列数と右の行列の行数が一致しません。");
var result = new Matrix(a.Rows, b.Cols);
for (int i = 0; i < result.Rows; i++)
{
for (int j = 0; j < result.Cols; j++)
{
double sum = 0;
for (int k = 0; k < a.Cols; k++)
sum += a[i, k] * b[k, j];
result[i, j] = sum;
}
}
return result;
}
public static Matrix operator *(Matrix a, double scalar)
{
var result = new Matrix(a.Rows, a.Cols);
for (int i = 0; i < a.Rows; i++)
for (int j = 0; j < a.Cols; j++)
result[i, j] = a[i, j] * scalar;
return result;
}
public static Matrix operator *(double scalar, Matrix a)
{
return a * scalar;
}
// --- 転置行列 (Transpose) ---
public Matrix Transpose()
{
var result = new Matrix(Cols, Rows); // 行と列のサイズを入れ替え
for (int i = 0; i < Rows; i++)
{
for (int j = 0; j < Cols; j++)
{
result[j, i] = _data[i, j]; // 成分を入れ替え
}
}
return result;
}
// --- 逆行列 (Inverse) ---
public Matrix Inverse()
{
if (Rows != Cols)
throw new InvalidOperationException("正方行列ではありません。");
int n = Rows;
double[,] augmented = new double[n, 2 * n];
for (int i = 0; i < n; i++)
{
for (int j = 0; j < n; j++)
{
augmented[i, j] = _data[i, j];
augmented[i, j + n] = (i == j) ? 1.0 : 0.0;
}
}
for (int i = 0; i < n; i++)
{
double pivot = augmented[i, i];
int pivotRow = i;
for (int k = i + 1; k < n; k++)
{
if (Math.Abs(augmented[k, i]) > Math.Abs(pivot))
{
pivot = augmented[k, i];
pivotRow = k;
}
}
if (Math.Abs(pivot) < 1e-10)
throw new InvalidOperationException("特異行列のため逆行列が存在しません。");
if (pivotRow != i)
{
for (int j = 0; j < 2 * n; j++)
{
double temp = augmented[i, j];
augmented[i, j] = augmented[pivotRow, j];
augmented[pivotRow, j] = temp;
}
}
for (int j = 0; j < 2 * n; j++) augmented[i, j] /= pivot;
for (int k = 0; k < n; k++)
{
if (k != i)
{
double factor = augmented[k, i];
for (int j = 0; j < 2 * n; j++) augmented[k, j] -= factor * augmented[i, j];
}
}
}
var inverse = new Matrix(n, n);
for (int i = 0; i < n; i++)
for (int j = 0; j < n; j++)
inverse[i, j] = augmented[i, j + n];
return inverse;
}
public override string ToString()
{
string s = "";
for (int i = 0; i < Rows; i++)
{
s += "| ";
for (int j = 0; j < Cols; j++) s += $"{_data[i, j],6:F2} ";
s += "|\n";
}
return s;
}
// 1. 各要素に関数を適用する (Map)
// 例: 行列の全要素にシグモイド関数を掛ける時などに使用
public Matrix Map(Func<double, double> func)
{
var result = new Matrix(Rows, Cols);
for (int i = 0; i < Rows; i++)
for (int j = 0; j < Cols; j++)
result[i, j] = func(_data[i, j]);
return result;
}
// 2. 要素ごとの掛け算 (アダマール積)
// 通常の行列積(*)とは異なり、同じ位置の成分同士を掛け合わせる
public static Matrix ElementWiseMultiply(Matrix a, Matrix b)
{
if (a.Rows != b.Rows || a.Cols != b.Cols)
throw new ArgumentException("行列のサイズが一致しません。");
var result = new Matrix(a.Rows, a.Cols);
for (int i = 0; i < a.Rows; i++)
for (int j = 0; j < a.Cols; j++)
result[i, j] = a[i, j] * b[i, j];
return result;
}
}
/// <summary>
/// n次元空間からm次元空間への投影を行うクラス
/// </summary>
public class NDimensionalProjector
{
public int SourceDimension { get; } // 投影元の次元数 (n)
public int TargetDimension { get; } // 投影先の次元数 (m)
// 投影時に残す軸のインデックス
private int[] _projectionAxes;
// 回転させる平面（軸iと軸j）と、その回転角度(ラジアン)のリスト
// 回転は順番によって結果が変わるため、リストで保持します
private readonly List<(int axis1, int axis2, double angle)> _rotations;
public NDimensionalProjector(int sourceDimension, int targetDimension)
{
if (sourceDimension <= 0 || targetDimension <= 0)
throw new ArgumentException("次元数は1以上である必要があります。");
if (targetDimension > sourceDimension)
throw new ArgumentException("投影先の次元(m)は、投影元の次元(n)以下である必要があります。");
SourceDimension = sourceDimension;
TargetDimension = targetDimension;
_rotations = new List<(int, int, double)>();
// デフォルトでは先頭から m 個の軸を残す
_projectionAxes = new int[targetDimension];
for (int i = 0; i < targetDimension; i++) _projectionAxes[i] = i;
}
/// <summary>
/// 投影時に残す（スクリーンとする）軸を指定します。
/// </summary>
/// <param name="axes">残す軸のインデックス（例: XとZを残すなら 0, 2）</param>
public void SetProjectionAxes(params int[] axes)
{
if (axes.Length != TargetDimension)
throw new ArgumentException($"指定する軸の数は投影先の次元数({TargetDimension})と一致する必要があります。");
foreach (int axis in axes)
{
if (axis < 0 || axis >= SourceDimension)
throw new ArgumentOutOfRangeException("指定された軸のインデックスが投影元の次元数の範囲外です。");
}
_projectionAxes = (int[])axes.Clone();
}
/// <summary>
/// 指定した2つの軸で構成される平面での回転を追加します。
/// </summary>
/// <param name="axis1">軸のインデックス (0からn-1)</param>
/// <param name="axis2">軸のインデックス (0からn-1)</param>
/// <param name="angleRadians">回転角度（ラジアン）</param>
public void AddRotation(int axis1, int axis2, double angleRadians)
{
if (axis1 < 0 || axis1 >= SourceDimension || axis2 < 0 || axis2 >= SourceDimension)
throw new ArgumentOutOfRangeException("軸のインデックスが範囲外です。");
if (axis1 == axis2)
throw new ArgumentException("異なる2つの軸を指定してください。");
_rotations.Add((axis1, axis2, angleRadians));
}
/// <summary>
/// 登録されたすべての回転をクリアします。
/// </summary>
public void ClearRotations()
{
_rotations.Clear();
}
/// <summary>
/// n x n の合成回転行列を生成します。
/// </summary>
public Matrix GetCombinedRotationMatrix()
{
// 初期状態は単位行列
Matrix combined = CreateIdentityMatrix(SourceDimension);
foreach (var rot in _rotations)
{
Matrix r = CreateGivensRotationMatrix(SourceDimension, rot.axis1, rot.axis2, rot.angle);
// 回転行列を掛け合わせる (後から追加した回転が左から掛かるようにする)
combined = r * combined;
}
return combined;
}
/// <summary>
/// m x n の射影行列を生成します。(正投影)
/// 指定されたm個の成分を残し、残りを捨てる行列です。
/// </summary>
public Matrix GetProjectionMatrix()
{
var proj = new Matrix(TargetDimension, SourceDimension);
for (int i = 0; i < TargetDimension; i++)
{
proj[i, _projectionAxes[i]] = 1.0;
}
return proj;
}
/// <summary>
/// n次元の座標を回転させ、m次元に投影します。
/// </summary>
/// <param name="point">n行1列の行列（列ベクトル）</param>
/// <returns>m行1列の行列（列ベクトル）</returns>
public Matrix Project(Matrix point)
{
if (point.Rows != SourceDimension || point.Cols != 1)
throw new ArgumentException($"入力点は {SourceDimension}行 1列 の列ベクトルである必要があります。");
Matrix rotationMatrix = GetCombinedRotationMatrix();
Matrix projectionMatrix = GetProjectionMatrix();
// P * R * x の計算
Matrix rotatedPoint = rotationMatrix * point;
Matrix projectedPoint = projectionMatrix * rotatedPoint;
return projectedPoint;
}
// --- ヘルパーメソッド ---
private Matrix CreateIdentityMatrix(int size)
{
var m = new Matrix(size, size);
for (int i = 0; i < size; i++) m[i, i] = 1.0;
return m;
}
// 2つの軸で構成される平面での回転行列 (Givens回転行列)
private Matrix CreateGivensRotationMatrix(int size, int i, int j, double angle)
{
Matrix m = CreateIdentityMatrix(size);
m[i, i] = Math.Cos(angle);
m[j, j] = Math.Cos(angle);
m[i, j] = -Math.Sin(angle);
m[j, i] = Math.Sin(angle);
return m;
}
}
}