E126. Нахождение ранга матрицы

e-maxx algorithm original: C/C++ #algorithm #emaxx #linear-algebra #math #matrix

Der Aufgabentext wird für die gewählte Sprache aus dem Russischen übersetzt. Code bleibt unverändert.

Static

Источник: e-maxx.ru/algo, страница PDF 415.

Ранг матрицы - это наибольшее number линейно независимых строк/столбцов матрицы. Ранг определён не только для квадратных матриц; пусть матрица прямоугольна и имеет размер NxM. Также ранг матрицы можно определить как наибольший из порядков миноров матрицы, отличных от нуля. Заметим, что если матрица квадратная и её определитель отличен от нуля, то ранг равен N(=M), иначе он будет меньше. В общем случае, ранг матрицы не превосходит min(N,M).

Algorithmus

Искать ранг можно с помощью модифицированного метода Гаусса. Будем выполнять абсолютно те же самые операции, что и при решении системы или нахождении её определителя, но если на каком-либо шаге в i-ом столбце среди невыбранных до этого строк нет ненулевых, то мы этот шаг пропускаем, а ранг уменьшаем на единицу (изначально ранг полагаем равным max(N,M)). Иначе, если мы нашли на i-ом шаге строку с ненулевым elementом в i- ом столбце, то помечаем эту строку как выбранную, и выполняем обычные операции отнимания этой строки от остальных.

Implementierung

const double EPS = 1E-9;

int rank = max(n,m);

vector<char> line_used (n);

for (int i=0; i<m; ++i) {

int j;

for (j=0; j<n; ++j)

if (!line_used[j] && abs(a[j][i]) > EPS)

break;

if (j == n)

--rank;

else {

line_used[j] = true;

for (int p=i+1; p<m; ++p)

a[j][p] /= a[j][i];

for (int k=0; k<n; ++k)

if (k != j && abs (a[k][i]) > EPS)

for (int p=i+1; p<m; ++p)

a[k][p] -= a[j][p] * a[k][i];

} }

Вычисление определителя матрицы методом Гаусса

Пусть дана квадратная матрица A размером NxN. it is required вычислить её определитель.

Algorithmus

Воспользуемся идеями метода Гаусса решения систем линейных уравнений. Будем выполнять те же самые действия, что и при решении системы линейных уравнений, исключив только деление текущей строки на a[i][i] (точнее, само деление можно выполнять, но подразумевая, что number выносится за знак определителя). Тогда все операции, которые мы будем производить с матрицей, не будут изменять величину определителя матрицы, за исключением, быть может, знака (мы только обмениваем местами две строки, что меняет знак на противоположный, или прибавляем одну строку к другой, что не меняет величину определителя). Но матрица, к которой мы приходим после выполнения Algorithmusа Гаусса, является диагональной, и определитель её равен произведению elementов, стоящих на диагонали. Знак, как уже говорилось, будет определяться количеством обменов строк (если их нечётное, то знак определителя следует изменить на противоположный). Таким образом, мы можем с помощью Algorithmusа Гаусса вычислять определитель матрицы за O (N3). Осталось только заметить, что если в какой-то момент мы не найдём в текущем столбце ненулевого elementа, то Algorithmus следует остановить и вернуть 0.

Implementierung

const double EPS = 1E-9;

int n;

vector < vector<double> > a (n, vector<double> (n));

... чтение n и a ...

double det = 1;

for (int i=0; i<n; ++i) {

int k = i;

for (int j=i+1; j<n; ++j)

if (abs (a[j][i]) > abs (a[k][i]))

k = j;

if (abs (a[k][i]) < EPS) {

det = 0;

break;

}

swap (a[i], a[k]);

if (i != k)

det = -det;

det *= a[i][i];

for (int j=i+1; j<n; ++j)

a[i][j] /= a[i][i];

for (int j=0; j<n; ++j)

if (j != i && abs (a[j][i]) > EPS)

for (int k=i+1; k<n; ++k)

a[j][k] -= a[i][k] * a[j][i];

}

cout << det;

C# Lösung

Auto-Entwurf, vor dem Einreichen prüfen

using System;
using System.Collections.Generic;
using System.Linq;

public static class AlgorithmDraft
{
    // Auto-generated C# draft from the original e-maxx C/C++ listing. Review before production use.
    const double EPS = 1E-9;
    int rank = max(n,m);
    List<char> line_used (n);
    for (int i=0; i<m; ++i) {
            int j;
            for (j=0; j<n; ++j)
                    if (!line_used[j] && abs(a[j][i]) > EPS)
                            break;
            if (j == n)
                    --rank;
            else {
                    line_used[j] = true;
                    for (int p=i+1; p<m; ++p)
                            a[j][p] /= a[j][i];
                    for (int k=0; k<n; ++k)
                            if (k != j && abs (a[k][i]) > EPS)
                                    for (int p=i+1; p<m; ++p)
                                            a[k][p] -= a[j][p] * a[k][i];
            }
    }
    const double EPS = 1E-9;
    int n;
    vector < vector<double> > a (n, vector<double> (n));
    ... чтение n и a ...
    double det = 1;
    for (int i=0; i<n; ++i) {
            int k = i;
            for (int j=i+1; j<n; ++j)
                    if (abs (a[j][i]) > abs (a[k][i]))
                            k = j;
            if (abs (a[k][i]) < EPS) {
                    det = 0;
                    break;
            }
            swap (a[i], a[k]);
            if (i != k)
                    det = -det;
            det *= a[i][i];
            for (int j=i+1; j<n; ++j)
                    a[i][j] /= a[i][i];
            for (int j=0; j<n; ++j)
                    if (j != i && abs (a[j][i]) > EPS)
                            for (int k=i+1; k<n; ++k)
                                    a[j][k] -= a[i][k] * a[j][i];
    }
    Console.WriteLine( det;
}

C++ Lösung

zugeordnet/original

const double EPS = 1E-9;
int rank = max(n,m);
vector<char> line_used (n);
for (int i=0; i<m; ++i) {
        int j;
        for (j=0; j<n; ++j)
                if (!line_used[j] && abs(a[j][i]) > EPS)
                        break;
        if (j == n)
                --rank;
        else {
                line_used[j] = true;
                for (int p=i+1; p<m; ++p)
                        a[j][p] /= a[j][i];
                for (int k=0; k<n; ++k)
                        if (k != j && abs (a[k][i]) > EPS)
                                for (int p=i+1; p<m; ++p)
                                        a[k][p] -= a[j][p] * a[k][i];
        }
}
const double EPS = 1E-9;
int n;
vector < vector<double> > a (n, vector<double> (n));
... чтение n и a ...
double det = 1;
for (int i=0; i<n; ++i) {
        int k = i;
        for (int j=i+1; j<n; ++j)
                if (abs (a[j][i]) > abs (a[k][i]))
                        k = j;
        if (abs (a[k][i]) < EPS) {
                det = 0;
                break;
        }
        swap (a[i], a[k]);
        if (i != k)
                det = -det;
        det *= a[i][i];
        for (int j=i+1; j<n; ++j)
                a[i][j] /= a[i][i];
        for (int j=0; j<n; ++j)
                if (j != i && abs (a[j][i]) > EPS)
                        for (int k=i+1; k<n; ++k)
                                a[j][k] -= a[i][k] * a[j][i];
}
cout << det;

Java Lösung

Auto-Entwurf, vor dem Einreichen prüfen

import java.util.*;
import java.math.*;

public class AlgorithmDraft {
    // Auto-generated Java draft from the original e-maxx C/C++ listing. Review before production use.
    const double EPS = 1E-9;
    int rank = max(n,m);
    ArrayList<Character> line_used (n);
    for (int i=0; i<m; ++i) {
            int j;
            for (j=0; j<n; ++j)
                    if (!line_used[j] && abs(a[j][i]) > EPS)
                            break;
            if (j == n)
                    --rank;
            else {
                    line_used[j] = true;
                    for (int p=i+1; p<m; ++p)
                            a[j][p] /= a[j][i];
                    for (int k=0; k<n; ++k)
                            if (k != j && abs (a[k][i]) > EPS)
                                    for (int p=i+1; p<m; ++p)
                                            a[k][p] -= a[j][p] * a[k][i];
            }
    }
    const double EPS = 1E-9;
    int n;
    vector < vector<double> > a (n, vector<double> (n));
    ... чтение n и a ...
    double det = 1;
    for (int i=0; i<n; ++i) {
            int k = i;
            for (int j=i+1; j<n; ++j)
                    if (abs (a[j][i]) > abs (a[k][i]))
                            k = j;
            if (abs (a[k][i]) < EPS) {
                    det = 0;
                    break;
            }
            swap (a[i], a[k]);
            if (i != k)
                    det = -det;
            det *= a[i][i];
            for (int j=i+1; j<n; ++j)
                    a[i][j] /= a[i][i];
            for (int j=0; j<n; ++j)
                    if (j != i && abs (a[j][i]) > EPS)
                            for (int k=i+1; k<n; ++k)
                                    a[j][k] -= a[i][k] * a[j][i];
    }
    System.out.println( det;
}

Материал разбит как Algorithmusическая Aufgabe: изучить постановку, понять асимптотику и реализовать Algorithmus на выбранном языке.

Stellen zu dieser Aufgabe

aktive Stellen with overlapping task tags are angezeigt.

Alle Stellen

Es gibt noch keine aktiven Stellen.