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

Initial commit

parents
.vs
obj
bin
using System;
using System.Collections.Generic;
using System.ComponentModel;
using System.Text;
namespace SVD
{
static class HouseHolder
{
static bool output1 = true; //пошаговый вывод для метода House
static bool output2 = false; //пошаговый вывод для метода HouseRow
static bool output3 = false; //пошаговый вывод для метода HouseCol
static bool output4 = true; //пошаговый вывод для метода Diag
//Вычисление вектора Хаусхолдера (работает корректно, проверено на бумаге)
public static Matrix House(Matrix vector1) //стр. 183
{
Matrix vector = new Matrix(vector1);
if (vector.N != 1)
{
throw new Exception("Not a vector!");
}
int n = vector.M;
float q = Matrix.VectorNorm(vector);
//>>>
if (output1)
{
Console.WriteLine("========= Вычисление вектора Хаусхолдера =========\n");
Console.WriteLine(" Переданный вектор:");
vector.Output();
}
//<<<
if (q != 0 )
{
float el = vector[0, 0];
float beta = el + Math.Sign(el) * q;
//>>
if (output1)
{
Console.WriteLine(" 1. Берем нулевой элемент:");
Console.WriteLine(" vector[0] == " + el);
Console.WriteLine();
Console.WriteLine(" 2. Вычисляем beta:");
Console.WriteLine(" beta = vector[0] + sign(vector[0]) * норма_вектора = {0} + {1} * {2} = {3}", el, Math.Sign(el), q, beta);
Console.WriteLine();
}
//<<
for (int i = 1; i < n; i++)
{
vector[i, 0] = vector[i, 0] / beta;
}
//>>
if (output1)
{
Console.WriteLine(" 3. Делим каждый элемент кроме нулевого на beta и получаем следующий вектор:");
vector.Output();
}
//<<
}
vector[0, 0] = 1;
//>>
if (output1)
{
Console.WriteLine(" 4. В конце устанавливаем нулевой элемент равным единице:");
vector.Output();
Console.WriteLine(" Метод Хаусхолдера завершен. Возвращение вектора...");
Console.WriteLine("====================================================\n");
}
//<<
return vector;
}
//Умножение на матрицу Хаусхолдера слева
public static Matrix HouseRow(Matrix A1, Matrix v) //стр. 184
{
Matrix A = new Matrix(A1);
Matrix AT = Matrix.Transp(A);
Matrix vT = Matrix.Transp(v);
//>>
if (output2)
{
Console.WriteLine("========= Умножение на матрицу Хаусхолдера слева =========\n");
Console.WriteLine(" Переданная матрица A:");
A.Output();
Console.WriteLine(" Переданный вектор v:");
v.Output();
Console.WriteLine(" 1. Транспонируем v:");
vT.Output();
}
//<<
float temp = Matrix.Multiply(vT, v)[0, 0]; //берем значение из полученной матрицы 1 на 1
float beta = -2 / temp;
//>>
if (output2)
{
Console.WriteLine(" 2. Вычисляем beta:\n");
Console.WriteLine(" beta = -2 / (vT * v) = -2 / {0} = {1}", temp, beta);
Console.WriteLine();
}
//<<
Matrix w = Matrix.Multiply(AT, v);
w = Matrix.MultiplyByNumber(w, beta);
//>>
if (output2)
{
Console.WriteLine(" 3. Транспонируем A:");
AT.Output();
Console.WriteLine(" 4. Вычисляем AT * v:");
Matrix.Multiply(AT, v).Output();
Console.WriteLine(" 5. Вычисляем w:\n");
Console.WriteLine(" w = beta * (AT * v) =");
w.Output();
}
//<<
Matrix tempA = Matrix.Multiply(v, Matrix.Transp(w));
A = Matrix.Addition(A, tempA);
//>>
if (output2)
{
Console.WriteLine(" 6. Транспонируем w:");
Matrix.Transp(w).Output();
Console.WriteLine(" 7. Вычисляем v * wT:");
tempA.Output();
Console.WriteLine(" 8. Вычисляем A:\n");
Console.WriteLine(" A = A + v * wT =");
A.Output();
Console.WriteLine(" Метод завершен. Возвращение матрицы...");
Console.WriteLine("====================================================\n");
}
//<<
return A;
}
//Умножение на матрицу Хаусхолдера справа
public static Matrix HouseCol(Matrix A1, Matrix v) //стр. 184
{
Matrix A = new Matrix(A1);
Matrix vT = Matrix.Transp(v);
//>>
if (output3)
{
Console.WriteLine("========= Умножение на матрицу Хаусхолдера справа =========\n");
Console.WriteLine(" Переданная матрица A:");
A.Output();
Console.WriteLine(" Переданный вектор v:");
v.Output();
Console.WriteLine(" 1. Транспонируем v:");
vT.Output();
}
//<<
float temp = Matrix.Multiply(vT, v)[0, 0];
float beta = -2 / temp;
//>>
if (output3)
{
Console.WriteLine(" 2. Вычисляем beta:\n");
Console.WriteLine(" beta = -2 / (vT * v) = -2 / {0} = {1}", temp, beta);
Console.WriteLine();
}
//<<
Matrix w = Matrix.Multiply(A, v);
w = Matrix.MultiplyByNumber(w, beta);
//>>
if (output3)
{
Console.WriteLine(" 3. Вычисляем A * v:");
Matrix.Multiply(A, v).Output();
Console.WriteLine(" 4. Вычисляем w:\n");
Console.WriteLine(" w = beta * (A * v) =");
w.Output();
}
//<<
Matrix tempA = Matrix.Multiply(w, vT);
A = Matrix.Addition(A, tempA);
//>>
if (output3)
{
Console.WriteLine(" 5. Вычисляем w * vT:");
Matrix.Multiply(w, vT).Output();
Console.WriteLine(" 6. Вычисляем A:\n");
Console.WriteLine(" A = A + w * vT =");
A.Output();
Console.WriteLine(" Метод завершен. Возвращение матрицы...");
Console.WriteLine("====================================================\n");
}
//<<
return A;
}
//Двухдиагонализация Хаусхолдера
public static Matrix Diag(Matrix A1)
{
Matrix A = new Matrix(A1);
int m = A.M;
int n = A.N;
Matrix v = new Matrix(m, 1);
//>>
if (output4)
{
Console.WriteLine(" Исходная матрица:");
A.Output();
}
//<<
for (int j = 0; j < n; j++)
{
for (int i = 0; i < v.M; i++)
{
v[i, 0] = 0;
}
//>>
Console.WriteLine(" Вектор v:");
v.Output();
//<<
Matrix subMatrix1 = Matrix.SubMatrix(A, j, m - 1, j, j);
Matrix house1 = HouseHolder.House(subMatrix1);
int t = j;
for (int i = j; i < m; i++)
{
v[i, 0] = house1[i - t, 0];
}
//>>
if (output4)
{
Console.WriteLine("================= Начало цикла при j == {0} ==================\n", j);
Console.WriteLine(" 1. Выделяем из матрицы A вектор, состоящий из элементов столбца {0} в диапозоне индексов {1}-{2}.", j, j, m - 1);
subMatrix1.Output();
Console.WriteLine(" 2. Вычисляем вектор Хаусхолдера:");
house1.Output();
Console.WriteLine(" 3. Присваиваем элементам вектора v с индексами из диапозона {0}-{1} соответствующие значения элементов вектора Хаусхолдера. Все остальные элементы полагаем равными нулю.", j, m - 1);
Console.WriteLine("\n Вектор v:");
v.Output();
}
//<<
Matrix subMatrix2 = Matrix.SubMatrix(A, j, m - 1, j, n - 1);
Matrix subMatrix3 = Matrix.SubMatrix(v, j, m - 1, 0, 0);
Matrix rowHouse = HouseHolder.HouseRow(subMatrix2, subMatrix3);
for (int i = j; i < m; i++)
{
for (int k = j; k < n; k++)
{
A[i, k] = rowHouse[i - t, k - t];
}
}
//>>
if (output4)
{
Console.WriteLine(" 4. Выделяем из матрицы A подматрицу в диапозоне {0}-{1} по строкам и {2}-{3} по столбцам:", j, m - 1, j, n - 1);
subMatrix2.Output();
Console.WriteLine(" 5. Выделяем из вектора v вектор в диапозоне {0}-{1}:", j, m - 1);
subMatrix3.Output();
Console.WriteLine(" 6. Умножим на матрицу Хаусхолдера слева и получим новую матрицу:");
rowHouse.Output();
Console.WriteLine(" 7. Присвоим значения элементов новой матрицы соответствующим значениям элементов матрицы А, лежащих в диапозоне {0}-{1} строк и {2}-{3} столбцов. Тогда матрица А будет иметь следующий вид:", j, m - 1, j, n - 1);
A.Output();
}
//<<
for (int i = j + 1; i < m; i++)
{
A[i, j] = v[i, 0];
}
//>>
if (output4)
{
Console.WriteLine(" 8. Присвоим значения элементов вектора v в диапозоне {0}-{1} значениям элементов матрицы А в диапозоне строк {2}-{3} столбца {4}", j + 1, m - 1, j + 1, m - 1, j);
Console.WriteLine("\n Вектор v:");
v.Output();
Console.WriteLine("\n Матрица A:");
A.Output();
}
//<<
if (j <= n - 3)
{
Matrix subMatrix4 = Matrix.SubMatrix(A, j, j, j + 1, n - 1);
Matrix house2 = HouseHolder.House(Matrix.Transp(subMatrix4));
for (int i = j + 1; i < n; i++)
{
v[i, 0] = house2[i - t - 1, 0];
}
//>>
if (output4)
{
Console.WriteLine(" 9. Условие j меньше или равно n - 3 выполняется, т.к. {0} меньше или равно {1}. Поэтому выполняем следующие действия.\n", j, n - 2);
Console.WriteLine(" 10. Выделяем из матрицы А вектор в строке {0} в диапозоне столбцов {1}-{2}:", j, j + 1, n - 1);
subMatrix4.Output();
Console.WriteLine(" 11. Транспонируем:");
Matrix.Transp(subMatrix4).Output();
Console.WriteLine(" 12. Вычисляем от этого вектора вектор Хаусхолдера:");
house2.Output();
Console.WriteLine(" 13. Присваиваем элементы вектора Хаусахолдера элементам вектора v в диапозоне {0}-{1}:", j + 1, n - 1);
Console.WriteLine("\n Вектор v:");
v.Output();
}
//<<
Matrix subMatrix5 = Matrix.SubMatrix(A, j, m - 1, j + 1, n - 1);
Matrix subMatrix6 = Matrix.SubMatrix(v, j + 1, n - 1, 0, 0);
Matrix colHouse = HouseHolder.HouseCol(subMatrix5, subMatrix6);
for (int i = j; i < m; i++)
{
for (int k = j + 1; k < n; k++)
{
A[i, k] = colHouse[i - t, k - t - 1];
}
}
//>>
if (output4)
{
Console.WriteLine(" 14. Выделяем из матрицы А подматрицу в диапозоне {0}-{1} по строкам и {2}-{3} по столбцам:", j, m - 1, j + 1, n - 1);
subMatrix5.Output();
Console.WriteLine(" 15. Выделяем из вектора v вектор в диапозоне {0}-{1}:", j + 1, n - 1);
subMatrix6.Output();
Console.WriteLine(" 16. Умножим на матрицу Хаусхолдера справа и получим новую матрицу:");
colHouse.Output();
Console.WriteLine(" 17. Присвоим значения элементов новой матрицы соответствующим значениям элементов матрицы А в диапозоне {0}-{1} по строкам и {2}-{3} по столбцам:", j, m - 1, j + 1, n - 1);
Console.WriteLine("\n Матрица A:");
A.Output();
}
//<<
Matrix subMatrix7 = Matrix.SubMatrix(v, j + 2, n - 1, 0, 0);
Matrix subMatrix8 = Matrix.Transp(subMatrix7);
for (int i = j + 2; i < n; i++)
{
A[j, i] = subMatrix8[0, i - t - 2];
}
//>>
if (output4)
{
Console.WriteLine(" 18. Выделим из вектора v вектор в диапозоне {0}-{1}:", j + 2, n - 1);
Console.WriteLine("\n Вектор v:");
v.Output();
Console.WriteLine("\n Выделенный вектор:");
subMatrix7.Output();
Console.WriteLine(" 19. Транспонируем его:");
subMatrix8.Output();
Console.WriteLine(" 20. Присваиваем значения его элементов элементам матрицы А в строке {0} с диапозоном {1}-{2} по столбцам:", j, j + 2, n - 1);
Console.WriteLine(" 21. Цикл при j == {0} завершен. Полученная матрица А:", j);
A.Output();
Console.WriteLine("============================================================");
}
//<<
}
else
{
Console.WriteLine(" 9. Условие j меньше или равно n - 2 НЕ выполняется, т.к. {0} больше {1}. Поэтому заканчиваем цикл при j == {2}.\n", j, n - 2, j);
}
}
//>>
Console.WriteLine("------------------------------------------------------");
Console.WriteLine(" Метод завершен ");
Console.WriteLine(" Полученная матрица ");
Console.WriteLine("------------------------------------------------------");
A.Output();
//<<
return A;
}
}
}
using System;
using System.Collections.Generic;
using System.Numerics;
using System.Runtime.Intrinsics;
namespace SVD
{
class Matrix
{
//Fields
protected float[,] data;
protected int m;
protected int n;
//Properties
public int M { get => this.m; }
public int N { get => this.n; }
//Constructor
public Matrix(int m, int n)
{
this.m = m;
this.n = n;
this.Initialize();
}
public Matrix(int m, int n, float value)
{
this.m = m;
this.n = n;
this.FillMatrix(value);
}
public Matrix(Matrix matrix)
{
m = matrix.m;
n = matrix.n;
data = new float[m, n];
for (int i = 0; i < m; i++)
{
for (int j = 0; j < n; j++)
{
data[i, j] = matrix.data[i, j];
}
}
}
//Overloading
public float this[int index1, int index2]
{
get
{
return data[index1, index2];
}
set
{
data[index1, index2] = value;
}
}
//Methods
protected void Initialize()
{
var rand = new Random();
data = new float[m, n];
for (int i = 0; i < m; i++)
{
for (int j = 0; j < n; j++)
{
data[i, j] = rand.Next() % 10 + 1;
}
}
}
protected void FillMatrix(float value)
{
data = new float[m, n];
for (int i = 0; i < m; i++)
{
for (int j = 0; j < n; j++)
{
data[i, j] = value;
}
}
}
public void Output()
{
Console.WriteLine();
for (int i = 0; i < m; i++)
{
Console.Write("{0, 15}", "║");
for (int j = 0; j < n; j++)
{
Console.Write("{0, 15}", data[i, j]);
}
Console.WriteLine("{0, 15}", "║");
}
Console.WriteLine();
}
public static Matrix Multiply(Matrix A, Matrix B)
{
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;
for (int k = 0; k < A.N; k++)
{
C[i, j] = C[i, j] + A[i, k] * B[k, j];
}
}
}
return C;
}
else
{
Console.WriteLine("================ ERROR ==================");
Console.WriteLine("=========== Incorrect input =============");
Console.WriteLine("=========================================");
Console.WriteLine();
throw new Exception();
}
}
public static Matrix MultiplyByNumber(Matrix A1, float a)
{
Matrix A = new Matrix(A1);
for (int i = 0; i < A.M; i++)
{
for (int j = 0; j < A.N; j++)
{
float value = A[i, j];
A[i, j] = value * a;
}
}
return A;
}
public static Matrix Addition(Matrix A, Matrix B)
{
if (A.M == B.M && A.N == B.N)
{
int m = A.M;
int n = A.N;
Matrix C = new Matrix(m, n);
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;
}
}
return C;
}
else
{
throw new Exception();
}
}
public static Matrix SubMatrix(Matrix A, int a1, int a2, int b1, int b2)
{
Matrix B = new Matrix(a2 - a1 + 1, b2 - b1 + 1);
int a = a1;
int b = b1;
for (int i = a1; i <= a2; i++)
{
for (int j = b1; j <= b2; j++)
{
B[i - a, j - b] = A[i, j];