Commit 496d88fe authored by Шарафутдинов Ратмир's avatar Шарафутдинов Ратмир
Browse files

Matrix class was improved

parent 32f62fda
......@@ -12,7 +12,8 @@ namespace SVD
public static bool output3 = false; //пошаговый вывод для метода HouseCol
public static bool output4 = true; //пошаговый вывод для метода Diag
//Вычисление вектора Хаусхолдера (работает корректно, проверено на бумаге)
//Вычисление вектора Хаусхолдера
public static Matrix House(Matrix vector1) //стр. 183
{
Matrix vector = new Matrix(vector1);
......@@ -23,7 +24,7 @@ namespace SVD
}
int n = vector.M;
float q = Matrix.VectorNorm(vector);
double q = vector.CalculateVectorNorm();
//>>>
if (output1)
......@@ -37,8 +38,8 @@ namespace SVD
if (q != 0 )
{
float el = vector[0, 0];
float beta = el + Math.Sign(el) * q;
double el = vector[0, 0];
double beta = el + Math.Sign(el) * q;
//>>
if (output1)
......@@ -86,8 +87,8 @@ namespace SVD
public static Matrix HouseRow(Matrix A1, Matrix v) //стр. 184
{
Matrix A = new Matrix(A1);
Matrix AT = Matrix.Transp(A);
Matrix vT = Matrix.Transp(v);
Matrix AT = A.CreateTransposeMatrix();
Matrix vT = v.CreateTransposeMatrix();
//>>
if (output2)
......@@ -102,9 +103,9 @@ namespace SVD
}
//<<
float temp = Matrix.Multiply(vT, v)[0, 0]; //берем значение из полученной матрицы 1 на 1
double temp = Matrix.Multiply(vT, v)[0, 0]; //берем значение из полученной матрицы 1 на 1
float beta = -2 / temp;
double beta = -2 / temp;
//>>
if (output2)
......@@ -131,14 +132,14 @@ namespace SVD
}
//<<
Matrix tempA = Matrix.Multiply(v, Matrix.Transp(w));
A = Matrix.Addition(A, tempA);
Matrix tempA = Matrix.Multiply(v, w.CreateTransposeMatrix());
A = Matrix.CreateSumMatrix(A, tempA);
//>>
if (output2)
{
Console.WriteLine(" 6. Транспонируем w:");
Matrix.Transp(w).Output();
w.CreateTransposeMatrix().Output();
Console.WriteLine(" 7. Вычисляем v * wT:");
tempA.Output();
Console.WriteLine(" 8. Вычисляем A:\n");
......@@ -157,7 +158,7 @@ namespace SVD
public static Matrix HouseCol(Matrix A1, Matrix v) //стр. 184
{
Matrix A = new Matrix(A1);
Matrix vT = Matrix.Transp(v);
Matrix vT = v.CreateTransposeMatrix();
//>>
if (output3)
......@@ -172,9 +173,9 @@ namespace SVD
}
//<<
float temp = Matrix.Multiply(vT, v)[0, 0];
double temp = Matrix.Multiply(vT, v)[0, 0];
float beta = -2 / temp;
double beta = -2 / temp;
//>>
if (output3)
......@@ -200,7 +201,7 @@ namespace SVD
//<<
Matrix tempA = Matrix.Multiply(w, vT);
A = Matrix.Addition(A, tempA);
A = Matrix.CreateSumMatrix(A, tempA);
//>>
if (output3)
......@@ -226,8 +227,6 @@ namespace SVD
int n = A.N;
Matrix v = new Matrix(m, 1);
//>>
if (output4)
{
......@@ -248,7 +247,7 @@ namespace SVD
v.Output();
//<<
Matrix subMatrix1 = Matrix.SubMatrix(A, j, m - 1, j, j);
Matrix subMatrix1 = A.CreateSubMatrix(j, m - 1, j, j);
Matrix house1 = HouseHolder.House(subMatrix1);
int t = j;
......@@ -273,8 +272,18 @@ namespace SVD
}
//<<
Matrix subMatrix2 = Matrix.SubMatrix(A, j, m - 1, j, n - 1);
Matrix subMatrix3 = Matrix.SubMatrix(v, j, m - 1, 0, 0);
//TEST>>
//У нас есть A и A1. Нужно найти U из U = A*A1^-1 //не работает не для квадр.
//PA = P*A => P = PA*A^(-1)
//TEST<<
Matrix subMatrix2 =A.CreateSubMatrix(j, m - 1, j, n - 1);
Matrix subMatrix3 = v.CreateSubMatrix(j, m - 1, 0, 0);
Matrix rowHouse = HouseHolder.HouseRow(subMatrix2, subMatrix3);
for (int i = j; i < m; i++)
......@@ -317,8 +326,8 @@ namespace SVD
if (j <= n - 3)
{
Matrix subMatrix4 = Matrix.SubMatrix(A, j, j, j + 1, n - 1);
Matrix house2 = HouseHolder.House(Matrix.Transp(subMatrix4));
Matrix subMatrix4 = A.CreateSubMatrix(j, j, j + 1, n - 1);
Matrix house2 = HouseHolder.House(subMatrix4.CreateTransposeMatrix());
for (int i = j + 1; i < n; i++)
{
......@@ -332,7 +341,7 @@ namespace SVD
Console.WriteLine(" 10. Выделяем из матрицы А вектор в строке {0} в диапозоне столбцов {1}-{2}:", j, j + 1, n - 1);
subMatrix4.Output();
Console.WriteLine(" 11. Транспонируем:");
Matrix.Transp(subMatrix4).Output();
subMatrix4.CreateTransposeMatrix().Output();
Console.WriteLine(" 12. Вычисляем от этого вектора вектор Хаусхолдера:");
house2.Output();
Console.WriteLine(" 13. Присваиваем элементы вектора Хаусахолдера элементам вектора v в диапозоне {0}-{1}:", j + 1, n - 1);
......@@ -341,8 +350,8 @@ namespace SVD
}
//<<
Matrix subMatrix5 = Matrix.SubMatrix(A, j, m - 1, j + 1, n - 1);
Matrix subMatrix6 = Matrix.SubMatrix(v, j + 1, n - 1, 0, 0);
Matrix subMatrix5 = A.CreateSubMatrix(j, m - 1, j + 1, n - 1);
Matrix subMatrix6 = v.CreateSubMatrix(j + 1, n - 1, 0, 0);
Matrix colHouse = HouseHolder.HouseCol(subMatrix5, subMatrix6);
for (int i = j; i < m; i++)
......@@ -368,8 +377,8 @@ namespace SVD
}
//<<
Matrix subMatrix7 = Matrix.SubMatrix(v, j + 2, n - 1, 0, 0);
Matrix subMatrix8 = Matrix.Transp(subMatrix7);
Matrix subMatrix7 = v.CreateSubMatrix(j + 2, n - 1, 0, 0);
Matrix subMatrix8 = subMatrix7.CreateTransposeMatrix();
for (int i = j + 2; i < n; i++)
{
......
......@@ -7,28 +7,32 @@ namespace SVD
{
class Matrix
{
//Fields
protected float[,] data;
//==================== Fields ====================
protected double[,] data;
protected int m;
protected int n;
//Properties
//==================== Properties ====================
public int M { get => this.m; }
public int N { get => this.n; }
//Constructor
public bool IsSquare { get => this.M == this.N; }
//==================== Constructors ====================
public Matrix(int m, int n)
public Matrix(int m, int n = 1)
{
this.m = m;
this.n = n;
this.Initialize();
//this.Initialize();
this.FillMatrix(0);
}
public Matrix(int m, int n, float value)
public Matrix(int m, int n, double value)
{
this.m = m;
this.n = n;
......@@ -40,7 +44,7 @@ namespace SVD
{
m = matrix.m;
n = matrix.n;
data = new float[m, n];
data = new double[m, n];
for (int i = 0; i < m; i++)
{
......@@ -51,9 +55,9 @@ namespace SVD
}
}
//Overloading
//==================== Overloading ====================
public float this[int index1, int index2]
public double this[int index1, int index2]
{
get
{
......@@ -65,12 +69,25 @@ namespace SVD
}
}
//Methods
public double this[int index]
{
get
{
return n == 1 ? data[index, 0] : throw new Exception("This matrix is not a vector. Error.");
}
set
{
data[index, 0] = (n == 1) ? value : throw new Exception("This matrix is not a vector. Error.");
}
}
protected void Initialize()
//==================== Non-Static Methods ====================
public void RandomInitialize() //заполнение элементами от 1 до 10
{
var rand = new Random();
data = new float[m, n];
data = new double[m, n];
for (int i = 0; i < m; i++)
{
......@@ -81,9 +98,9 @@ namespace SVD
}
}
protected void FillMatrix(float value)
protected void FillMatrix(double value)
{
data = new float[m, n];
data = new double[m, n];
for (int i = 0; i < m; i++)
{
......@@ -97,98 +114,130 @@ namespace SVD
public void Output()
{
Console.WriteLine();
for (int i = 0; i < m; i++)
{
Console.Write("{0, 15}", "║");
Console.Write("{0, 18}", "║");
for (int j = 0; j < n; j++)
{
Console.Write("{0, 15}", data[i, j]);
Console.Write("{0, 18}", Math.Round(data[i, j], 11));
}
Console.WriteLine("{0, 15}", "║");
Console.WriteLine("{0, 18}", "║");
}
Console.WriteLine();
}
public static Matrix Multiply(Matrix A, Matrix B)
private void ProcessFunctionOverData(Action<int, int> func)
{
if (A.N == B.M)
for (var i = 0; i < this.M; i++)
{
Matrix C = new Matrix(A.M, B.N);
for (int i = 0; i < A.M; i++)
for (var j = 0; j < this.N; j++)
{
for (int j = 0; j < B.N; j++)
{
C[i, j] = 0;
for (int k = 0; k < A.N; k++)
{
C[i, j] = C[i, j] + A[i, k] * B[k, j];
}
}
func(i, j);
}
return C;
}
else
}
private Matrix CreateMatrixWithoutColumn(int column)
{
if (column < 0 || column >= this.N)
{
Console.WriteLine("================ ERROR ==================");
Console.WriteLine("=========== Incorrect input =============");
Console.WriteLine("=========================================");
Console.WriteLine();
throw new ArgumentException("Invalid column index.");
}
throw new Exception();
Matrix result = new Matrix(this.M, this.N - 1);
result.ProcessFunctionOverData((i, j) => result[i, j] = j < column ? this[i, j] : this[i, j + 1]);
return result;
}
private Matrix CreateMatrixWithoutRow(int row)
{
if (row < 0 || row >= this.M)
{
throw new ArgumentException("Invalid row index.");
}
Matrix result = new Matrix(this.M - 1, this.N);
result.ProcessFunctionOverData((i, j) => result[i, j] = i < row ? this[i, j] : this[i + 1, j]);
return result;
}
public static Matrix MultiplyByNumber(Matrix A1, float a)
public double CalculateDeterminant()
{
Matrix A = new Matrix(A1);
if (!this.IsSquare)
{
throw new InvalidOperationException("Determinant can be calculated only for square matrix.");
}
for (int i = 0; i < A.M; i++)
if (this.N == 2)
{
for (int j = 0; j < A.N; j++)
{
float value = A[i, j];
A[i, j] = value * a;
return this[0, 0] * this[1, 1] - this[0, 1] * this[1, 0];
}
}
double result = 0;
for (var j = 0; j < this.N; j++)
{
result += (j % 2 == 1 ? 1 : -1) * this[1, j] * this.CreateMatrixWithoutColumn(j).CreateMatrixWithoutRow(1).CalculateDeterminant();
}
return (double)result;
}
return A;
private double CalculateMinor(int i, int j)
{
return CreateMatrixWithoutColumn(j).CreateMatrixWithoutRow(i).CalculateDeterminant();
}
public static Matrix Addition(Matrix A, Matrix B)
public Matrix CreateInvertibleMatrix()
{
if (A.M == B.M && A.N == B.N)
if (this.M != this.N)
{
int m = A.M;
int n = A.N;
return null;
}
Matrix C = new Matrix(m, n);
double determinant = CalculateDeterminant();
for (int i = 0; i < m; i++)
{
for (int j = 0; j < n; j++)
{
float value = A[i, j] + B[i, j];
C[i, j] = value;
if (determinant == 0)
{
return null;
}
}
}
Matrix result = new Matrix(M, M);
return C;
}
else
ProcessFunctionOverData((i, j) => result[i, j] = ((i + j) % 2 == 1 ? -1 : 1) * CalculateMinor(i, j) / determinant);
result = result.CreateTransposeMatrix();
return result;
}
public Matrix Transp()
{
int m = this.M;
int n = this.n;
Matrix result = new Matrix(n, m);
for (int i = 0; i < result.M; i++)
{
throw new Exception();
for (int j = 0; j < result.N; j++)
{
result[i, j] = this[j, i];
}
}
return result;
}
public static Matrix SubMatrix(Matrix A, int a1, int a2, int b1, int b2)
public Matrix CreateSubMatrix(int a1, int a2, int b1, int b2)
{
Matrix B = new Matrix(a2 - a1 + 1, b2 - b1 + 1);
Matrix result = new Matrix(a2 - a1 + 1, b2 - b1 + 1);
int a = a1;
int b = b1;
......@@ -197,25 +246,70 @@ namespace SVD
{
for (int j = b1; j <= b2; j++)
{
B[i - a, j - b] = A[i, j];
result[i - a, j - b] = this[i, j];
}
}
return B;
return result;
}
public static float VectorNorm(Matrix vector) //стр. 16
public Matrix CreateTransposeMatrix() //возвращает матрицу, которую можно присвоить
{
int m = vector.M;
int n = vector.N;
int m = this.M;
int n = this.n;
Matrix result = new Matrix(n, m);
for (int i = 0; i < result.M; i++)
{
for (int j = 0; j < result.N; j++)
{
result[i, j] = this[j, i];
}
}
return result;
}
public double[] VectorToArray()
{
double[] array = new double[this.M];
for (int i = 0; i < this.M; i++)
{
array[i] = this[i];
}
return array;
}
public double[,] MatrixToArray()
{
double[,] array = new double[this.M, this.N];
for (int i = 0; i < this.M; i++)
{
for (int j = 0; j < this.N; j++)
{
array[i, j] = this[i, j];
}
}
return array;
}
public double CalculateVectorNorm() //стр. 16
{
int m = this.M;
int n = this.N;
if (n == 1)
{
Matrix vectorT = Matrix.Transp(vector);
Matrix vectorT = this.Transp();
float element = Matrix.Multiply(vectorT, vector)[0, 0];
double element = Matrix.Multiply(vectorT, this)[0, 0];
element = (float)Math.Sqrt(element);
element = Math.Sqrt(element);
return element;
}
......@@ -230,22 +324,86 @@ namespace SVD
}
}
public static Matrix Transp(Matrix A)
//==================== Static Methods ====================
private static void Error(string error = "Undefined error")
{
Console.WriteLine("\n================ ERROR ==================\n");
Console.WriteLine("\t" + error);
Console.WriteLine("\n=========================================\n");
}
public static Matrix Multiply(Matrix A, Matrix B)
{
int m = A.M;
int n = A.n;
if (A.N == B.M)
{
Matrix C = new Matrix(A.M, B.N);
for (int i = 0; i < A.M; i++)
{
for (int j = 0; j < B.N; j++)
{
C[i, j] = 0;
Matrix temp = new Matrix(n, m);
for (int k = 0; k < A.N; k++)
{
C[i, j] = C[i, j] + A[i, k] * B[k, j];
}
}
}
for (int i = 0; i < temp.M; i++)
return C;
}
else
{
for (int j = 0; j < temp.N; j++)
Matrix.Error();
throw new Exception();
}
}
public static Matrix MultiplyByNumber(Matrix A, double a)
{
Matrix result = new Matrix(A);
for (int i = 0; i < result.M; i++)
{
for (int j = 0; j < result.N; j++)
{
temp[i, j] = A[j, i];
double value = result[i, j];
result[i, j] = value * a;
}
}
return temp;
return result;
}
public static Matrix CreateSumMatrix(Matrix A, Matrix B)
{
if (A.M == B.M && A.N == B.N)
{
int m = A.M;
int n = A.N;
Matrix result = new Matrix(m, n);
for (int i = 0; i < m; i++)
{
for (int j = 0; j < n; j++)