Obliczanie wyznacznika macierzy

Deti

No.. w końcu skończyłem pisać tą manianę :) Funkcja oblicza wyznacznik przez metodę klasyczną, czyli

(-1)^(wiersz*kolumna)*wyznacznik naprzeciwległy + ... itd każda kolumna

.. a wszystko przez rekurencję.

Najpierw utwórzmy sobie nowy typ

type
  TDynamicIntegerArray = array of array of Integer;

.. aby móc odpalić funkcję tworzymy sobie jakąś tablice (macierz):

var
  tab: TDynamicIntegerArray;

Teraz można przypisać sobie do niej jakieś wartości:

SetLength(tab, 3, 3);
  tab[0][0] := 1;
  tab[0][1] := 2;
  tab[0][2] := 0;

  tab[1][0] := 2;
  tab[1][1] := 3;
  tab[1][2] := 0;

  tab[2][0] := 1;
  tab[2][1] := -1;
  tab[2][2] := 1;

.. i jak już mamy to gotowe - można odpalić funkcję:

ShowMessage('Wyznacznik wynosi: ' + IntToStr(MatrixDet(tab)));

A teraz sama funkcja :)

function MatrixDet(Matrix : TDynamicIntegerArray): Integer;

(*
============================================================================
Matrix determinant function [Copyright (c) HAKGERSoft 2000 - 2004]

This code is owned by HAKGERSoft, any modifications without HAKGERSoft
permision is prohibited!

Author:
  DetoX [detox@satanbsd.org]

============================================================================
*)

var
  MXArr: TDynamicIntegerArray;       // Recursion array buffer
  MXInpArrX, MXInpArrY, MXArrX, MXArrY: Integer;

begin
  Result := 0;         // Zero determinant start !
  if High(Matrix) = 0 then           // Simple [1][1] matrix
    Result := Matrix[0][0]
  else if High(Matrix) = 1 then      // Simple [2][2] matrix
    Result := (Matrix[0][0] * Matrix[1][1]) - (Matrix[1][0] * Matrix[0][1])
  else
  begin                              // [>2][>2] array
    for MXInpArrX := 0 to High(Matrix) do
    begin
      SetLength(MXArr, High(Matrix), High(Matrix));   // Set inner array
      for MXInpArrY := 0 to High(Matrix) do
      begin
        MXArrX := 0;  // Zero X[D] array counter - for loop
        MXArrY := 0;  // Zero Y[D] array counter - for loop
        repeat
          if (MXArrX <> MXInpArrX) then   // Must be different X dimension
          begin
            if MXInpArrY + 1 > High(Matrix) then   // Check array size
              Break;
            MXArr[MXInpArrY][MXArrY] := Matrix[MXInpArrY + 1][MXArrX];  // Set MXArr Values
            Inc(MXArrY);  // Next Y
          end;
          Inc(MXArrX);  // Next X
        until MXArrX = High(Matrix) + 1;
      end;

      // Recursion determinant calculate
      if not Odd(MXInpArrX + 1) then
        Result := Result + (Matrix[0][MXInpArrX] * MatrixDet(MXArr))
      else
        Result := Result - (Matrix[0][MXInpArrX] * MatrixDet(MXArr));
    end;
  end;
end;

No i wszystko. Na pewno są lepsze sposoby na to obliczanie, ale chciałem sam takie coś napisać.. no i się udało :) - może komuś się przyda...

PS: sorry za komentarze po angielsku - ale jeszcze gdzie indziej to wysyłałem

5 komentarzy

wszystko fajnie ale czy ktos moze mi powiedziec w jakiej kolejności są reprezentowane kolumny i wiersze??
albo czy może mi ktoś powiedzieć jak przerobić żeby było tak: tab[kolumna][wiersz]:=1;

Jeżeli elementy macierzy są rzeczywiste to najlepiej skorzystać z eliminacji Gaussa lub rozłożyć macierz na iloczyn LU
W eliminacji Gaussa trzeba użyć częściowego wyboru elementu głównego (zamiana wierszy) a w rozkładzie LU
trzeba ponadto obliczyć macierz permutacji dla ustalenia znaku wyznacznika
Jeżeli elementy macierzy są liczbami całkowitymi to trzeba zastosować rekurencyjne rozwinięcie Laplace lub
ewentualnie wygenerować permutacje W tym przypadku metoda eliminacji Gaussa nie sprawdza się ponieważ
dzielenie nie jest wykonalne w całej dziedzinie (można wejść w zbiór liczb wymiernych)

lepiej jest chyba napisać eliminację gaussa a potem ponożyć elementy na przekątnej ;-)

No nie wiem - u mnie wszystko śmiga...

dziwny error chyba w lini: "InnerArray[y][k] := Matrix[y+1][x];" - "Access violation at address 00401B59 in module "Project1.exe". Read of address 0000001E."