Метод Гаусса

[Elarium]

Необходимо преобразовать матрицу 4x4 к треугольному (ступенчатому) виду. Долго мозговал решение в общем виде и в голове появился следующий код:

Используется синтаксис Pascal

unit DetGauss;

interface

uses
  Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
  Dialogs, Grids, StdCtrls;

type
  TForm1 = class(TForm)
    StringGrid1: TStringGrid;
    Button1: TButton;
    Label1: TLabel;
    Button2: TButton;
    StringGrid2: TStringGrid;
    procedure Button2Click(Sender: TObject);
    procedure Button1Click(Sender: TObject);
  private
    { Private declarations }
  public
    { Public declarations }
  end;

var
  Form1: TForm1;

implementation

{$R *.dfm}

procedure TForm1.Button2Click(Sender: TObject);
begin
close;
end;

procedure TForm1.Button1Click(Sender: TObject);
var a: array [0..3,0..3] of real;
    i,j,k: integer;
    N:real;
begin

{Заполняем матрицу элементами}

 for i:=0 to 3 do
  for j:=0 to 3 do
   a[i,j]:=strtofloat(stringgrid1.cells[i,j]);

{Преобразование к треугольному виду}

k:=0;
repeat
for i:=k+1 to 3 do
begin
N:=a[i,k]/a[k,k];
for j:=k to 3 do
a[i,j]:=a[i,j]-a[k,j]*N;
end;
k:=k+1;
until (k=3);

 for i:=0 to 3 do
   for j:=0 to 3 do
   stringgrid2.Cells[i,j]:=floattostr(a[i,j]);
end;

end.
 

Расписывая ход работы на листке бумаги вроде как элементы матрицы должны сводится к треугольному виду, а получается следующая штука:



Если кто поможет разобраться, буду очень благодарен :)

[lim0n]

Всё у Вас нормально, только разберитесь со строками/столбцами в массивах и объектах типа TStringGrig. Обычно в массивах принято первым индексом обозначать строки, вторым - столбцы, но в свойстве Cells всё наоборот:

property Cells[ACol, ARow: Integer]: string;

И, естественно, не забывайте про ошибки округления при работе с числами с плавающей точкой (это к тому, что "треугольности" не видно).

[Elarium]

lim0n в сообщении #366311 писал(а):

Всё у Вас нормально, только разберитесь со строками/столбцами в массивах и объектах типа TStringGrig. Обычно в массивах принято первым индексом обозначать строки, вторым - столбцы, но в свойстве Cells всё наоборот:

property Cells[ACol, ARow: Integer]: string;

И, естественно, не забывайте про ошибки округления при работе с числами с плавающей точкой (это к тому, что "треугольности" не видно).

Слишком уж большие эти ошибки округляения (-1,88 6,15 и -5,05), всё таки что-то не так...

[lim0n]

Elarium в сообщении #366324 писал(а):

Слишком уж большие эти ошибки округляения

Эти числа в "научной" нотации, порядок $10^{-18}$, вы их целиком не видите.

[Circiter]

2Elarium
Хм, странно... Если перенести подвыражение /a[k,k], разместив его непосредственно после a[k,j], то результат становится более правдоподобным... :)

[ewert]

Circiter в сообщении #366351 писал(а):

то результат становится более правдоподобным... :)

Потому, что меняется порядок операций: раньше было сначала деление, потом умножение, а теперь -- наоборот.

Кстати (хоть и непринципиально): трудно понять, зачем понадобилось оформлять внешний цикл именно repeat'ом.

[Circiter]

А ещё результирующую матрицу можно тупо проквантовать по некоторому эпсилону, шоб красивше было...

Стандартный прием с выбором максимального элемента тоже конечно надо использовать, но в данном случае он, имхо, не поможет. :)

2ewert

Цитата:

раньше было сначала деление, потом умножение, а теперь -- наоборот

Да не... В справке к delphi написано, что все инфиксные операторы левоассоциативны. То есть, поменялось только делимое с a[i,k] на a[k,j]... А странно здесь то, что a[k,k] -- единица (то есть здесь какой-то более тонкий эффект).

[lim0n]

Если автор заглянет: проверьте работу программы, например, на нулевой матрице.

Цитата:

... левоассоциативны ...

т.е. выполняются слева направо

((a[k,j]*a[i,k])/a[k,k])

[Circiter]

2lim0n
Позвольте сам себя процитирую:

Circiter писал(а):

перенести подвыражение /a[k,k], разместив его непосредственно после a[k,j]

То есть: a[k,j]/a[k,k]*N. Этот вариант придает матрице треугольный вид с нулями под диагональю. Причем сначала выполняется деление (ассоциативность влево), да еще и на единицу. А раньше тоже сначала выполнялось деление, тоже на единицу, но делилось уже другое число. Похоже, что здесь не обошлось без оптимизатора. Правда это уже оффтопик попёр... :)

[lim0n]

Circiter в сообщении #366379 писал(а):

То есть: a[k,j]/a[k,k]*N.

В таком случае, т.к. само $N=a_{ik}/a_{kk}$, Вы вычисляете что-то не то.

[Circiter]

2lim0n
Так, ещё раз себя процитирую: :)

Circiter писал(а):

перенести подвыражение /a[k,k], разместив его непосредственно после a[k,j]

Под "переносом" я имел ввиду вырезание и вставку. :) То есть после переноса N:=a[i,k].

[lim0n]

Теперь ясно.