src/cs/Matrix/Matrix.cs
using System.Collections.Generic;
using System.Linq; using static System.Math; namespace Prelude { public class Matrix { public int[] Size { get; private set; } private double[][] _Rows; public double[][] Rows { get { return _Rows; } set { int rows = Size[0], cols = Size[1]; if (value.Length > rows) { var limit = Min(value.Length, rows * cols); for (var i = 0; i < limit; ++i) { int row = (int)Floor((double)(i / cols)); int col = i % cols; _Rows[row][col] = value[i][0]; } } else { double[][] temp = Create(rows, cols); for (var row = 0; row < rows; ++row) temp[row] = value[row].Take(cols).ToArray(); _Rows = temp; } } } public Matrix(int n) { Size = new int[] { n, n }; Rows = Create(n, n); } public Matrix(int rows, int cols) { Size = new int[] { rows, cols }; Rows = Create(rows, cols); } public static double[][] Create(int rows, int cols) { double[][] result = new double[rows][]; for (int i = 0; i < rows; ++i) result[i] = new double[cols]; return result; } public static Matrix Unit(int n) { var temp = new Matrix(n); foreach (var index in temp.Indexes()) { int i = index[0], j = index[1]; temp.Rows[i][j] = 1; } return temp; } public static Matrix Identity(int n) { var temp = new Matrix(n); for (int i = 0; i < n; ++i) temp.Rows[i][i] = 1; return temp; } public static Matrix Transpose(Matrix a) { var temp = new Matrix(a.Size[1], a.Size[0]); foreach (var index in a.Indexes()) { int i = index[0], j = index[1]; temp.Rows[j][i] = a.Rows[i][j]; } return temp; } public static Matrix Add(params Matrix[] addends) { var size = addends[0].Size; var sum = new Matrix(size[0], size[1]); foreach (Matrix matrix in addends) foreach (var index in matrix.Indexes()) { int i = index[0], j = index[1]; sum.Rows[i][j] += matrix.Rows[i][j]; } return sum; } public static Matrix Adj(Matrix a) { var temp = a.Clone(); foreach (var index in temp.Indexes()) { int i = index[0], j = index[1]; temp.Rows[i][j] = a.Cofactor(i, j); } return Transpose(temp); } public static double Det(Matrix a) { int rows = a.Size[0]; switch (rows) { case 1: return a.Rows[0][0]; case 2: return a.Rows[0][0] * a.Rows[1][1] - a.Rows[0][1] * a.Rows[1][0]; default: double sum = 0; for (int i = 0; i < rows; ++i) sum += a.Rows[0][i] * a.Cofactor(0, i); return sum; } } public static Matrix Dot(Matrix a, Matrix b) { int m = a.Size[0], p = a.Size[1], n = b.Size[1]; var product = new Matrix(m, n); foreach (var index in product.Indexes()) { int i = index[0], j = index[1]; double sum = 0; for (int k = 0; k < p; ++k) { sum += a.Rows[i][k] * b.Rows[k][j]; } product.Rows[i][j] = sum; } return product; } public static Matrix Invert(Matrix a) { var adjugate = Adj(a); double det = Det(a); return Multiply(adjugate, 1 / det); } public static Matrix Multiply(Matrix a, double k) { var clone = a.Clone(); foreach (var index in clone.Indexes()) { int i = index[0], j = index[1]; clone.Rows[i][j] *= k; } return clone; } public static double Trace(Matrix a) { double trace = 0; foreach (var index in a.Indexes()) { int i = index[0], j = index[1]; if (i == j) { trace += a.Rows[i][j]; } } return trace; } public Matrix Clone() { var original = this; int rows = original.Size[0], cols = original.Size[1]; var clone = new Matrix(rows, cols); foreach (var index in clone.Indexes()) { int i = index[0], j = index[1]; clone.Rows[i][j] = original.Rows[i][j]; } return clone; } public double Cofactor(int i = 0, int j = 0) => Pow(-1, i + j) * Det(RemoveRow(i).RemoveColumn(j)); public List<int[]> Indexes(int offset = 0) { int rows = Size[0], cols = Size[1]; var pairs = new List<int[]>(); for (var i = 0; i < rows; ++i) for (var j = 0; j < cols; ++j) { int[] pair = { i + offset, j + offset }; pairs.Add(pair); } return pairs; } public Matrix RemoveColumn(int index) { var original = Clone(); int rows = original.Size[0], cols = original.Size[1]; if (index < 0 || index >= cols) { return original; } else { var temp = new Matrix(rows, cols - 1); for (var i = 0; i < rows; ++i) for (var j = 0; j < index; ++j) temp.Rows[i][j] = original.Rows[i][j]; for (var i = 0; i < rows; ++i) for (var j = index; j < cols - 1; ++j) temp.Rows[i][j] = original.Rows[i][j + 1]; return temp; } } public Matrix RemoveRow(int index) { var original = Clone(); int rows = original.Size[0], cols = original.Size[1]; if (index < 0 || index >= rows) { return original; } else { var temp = new Matrix(rows - 1, cols); for (var i = 0; i < index; ++i) for (var j = 0; j < cols; ++j) temp.Rows[i][j] = original.Rows[i][j]; for (var i = index; i < rows - 1; ++i) for (var j = 0; j < cols; ++j) temp.Rows[i][j] = original.Rows[i + 1][j]; return temp; } } public override string ToString() { var matrix = this; int rank = matrix.Size[0]; var rows = new string[rank]; for (var i = 0; i < rank; ++i) rows[i] = string.Join(",", matrix.Rows[i]); return string.Join("\r\n", rows); } } } |