Pascal code for matrix inversion

[finalist]

I have a symmetrical matrix A 8 x 8.
First I need to find matrix B = A[sup]T[/sup]*A;
Next I need to inverse B and find C = B[sup]-1[/sup];
Help please !

[FlamingClaw]

any code?

[finalist]

any code?

I am absolute beginner here, so I don't know have I some points to give them to you ....
If I have points, I will give you the points if you help me.
Thanks.

[finalist]

Please do not answer to this question, until I post an example here ....

[finalist]

I developed the multiplying of Matrices At*A, now I need to develope the Inverse of the result Matrice of multiplying ....

unit Unit1; 

{$mode objfpc}{$H+}

interface

uses
  Classes, SysUtils, FileUtil, LResources, Forms, Controls, Graphics, Dialogs,
  StdCtrls, Grids;

type
  T_Matrix_AT    = array[1..4, 1..8] of Double;
  T_Matrix_A     = array[1..8, 1..4] of Double;
  T_Matrix_MLTPL = array[1..4, 1..4] of Double;
  T_Matrix_INV   = array[1..4, 1..4] of Double;

  { TForm1 }

  TForm1 = class(TForm)
    Btn_Multiply_: TButton;
    Btn_Inverse_: TButton;
    StrGrid_A: TStringGrid;
    StrGrid_AT_A: TStringGrid;
    StrGrid_AT: TStringGrid;
    StrGrid_AT_A_INV: TStringGrid;
    procedure Btn_Multiply_Click(Sender: TObject);
    procedure Btn_Inverse_Click(Sender: TObject);
    procedure FormCreate(Sender: TObject);
  private
    _Matrix_AT:    T_Matrix_AT;
    _Matrix_A:     T_Matrix_A;
    _Matrix_MLTPL: T_Matrix_MLTPL;
    _Matrix_INV:   T_Matrix_INV;
    { private declarations }
  public
    { public declarations }
  end; 

var
  Form1: TForm1; 

implementation

{ TForm1 }

procedure TForm1.FormCreate(Sender: TObject);
var
  I:      Integer;
  J:      Integer;
begin
  _Matrix_A[1,1] := -0.993;
  _Matrix_AT[1,1] := -0.993;

  _Matrix_A[2,1] := -0.404;
  _Matrix_AT[1,2] := -0.404;

  _Matrix_A[3,1] := +0.877;
  _Matrix_AT[1,3] := +0.877;

  _Matrix_A[4,1] := +0.493;
  _Matrix_AT[1,4] := +0.493;

  _Matrix_A[5,1] := +0;
  _Matrix_AT[1,5] := +0;

  _Matrix_A[6,1] := +0;
  _Matrix_AT[1,6] := +0;

  _Matrix_A[7,1] := +0;
  _Matrix_AT[1,7] := +0;

  _Matrix_A[8,1] := +0;
  _Matrix_AT[1,8] := +0;

//................................
  _Matrix_A[1,2] := +0.119;
  _Matrix_AT[2,1] := +0.119;

  _Matrix_A[2,2] := -0.915;
  _Matrix_AT[2,2] := -0.915;

  _Matrix_A[3,2] := +0.480;
  _Matrix_AT[2,3] := +0.480;

  _Matrix_A[4,2] := -0.870;
  _Matrix_AT[2,4] := -0.870;

  _Matrix_A[5,2] := +0;
  _Matrix_AT[2,5] := +0;

  _Matrix_A[6,2] := +0;
  _Matrix_AT[2,6] := +0;

  _Matrix_A[7,2] := +0;
  _Matrix_AT[2,7] := +0;

  _Matrix_A[8,2] := +0;
  _Matrix_AT[2,8] := +0;
//................................
  _Matrix_A[1,3] := +0;
  _Matrix_AT[3,1] := +0;

  _Matrix_A[2,3] := -0;
  _Matrix_AT[3,2] := -0;

  _Matrix_A[3,3] := +0;
  _Matrix_AT[3,3] := +0;

  _Matrix_A[4,3] := -0.493;
  _Matrix_AT[3,4] := -0.493;

  _Matrix_A[5,3] := -0.790;
  _Matrix_AT[3,5] := -0.790;

  _Matrix_A[6,3] := -0.936;
  _Matrix_AT[3,6] := -0.936;

  _Matrix_A[7,3] := +0.163;
  _Matrix_AT[3,7] := +0.163;

  _Matrix_A[8,3] := +0.997;
  _Matrix_AT[3,8] := +0.997;

//................................

  _Matrix_A[1,4] := +0;
  _Matrix_AT[4,1] := +0;

  _Matrix_A[2,4] := -0;
  _Matrix_AT[4,2] := -0;

  _Matrix_A[3,4] := +0;
  _Matrix_AT[4,3] := +0;

  _Matrix_A[4,4] := +0.870;
  _Matrix_AT[4,4] := +0.870;

  _Matrix_A[5,4] := +0.613;
  _Matrix_AT[4,5] := +0.613;

  _Matrix_A[6,4] := -0.353;
  _Matrix_AT[4,6] := -0.353;

  _Matrix_A[7,4] := +0.987;
  _Matrix_AT[4,7] := +0.987;

  _Matrix_A[8,4] := -0.077;
  _Matrix_AT[4,8] := -0.077;
//................................

  for I := 1 to 8 do
  begin
    for J := 1 to 4 do
    begin
      StrGrid_AT.Cells[I, J] := FloatToSTr(_Matrix_AT[J, I]);
      StrGrid_A.Cells[J, I] := FloatToSTr(_Matrix_A[I, J]);
    end;
  end;
  StrGrid_AT.Cells[0,0] := 'At';
  StrGrid_A.Cells[0,0] := 'A';
  StrGrid_AT_A.Cells[0,0] := 'At * A';
  StrGrid_AT_A_INV.Cells[0,0] := 'INVERSED';
end;

procedure TForm1.Btn_Multiply_Click(Sender: TObject);
var
  I:      Integer;
  J:      Integer;
  K:      Integer;
  D:      Double;
begin
  for K := 1 to 4 do
  begin
    for J := 1 to 4 do
    begin
      D := 0;
      for I := 1 to 8 do
      begin
        D := D + _Matrix_AT[K,I]*_Matrix_A[I,J];
      end;
      StrGrid_AT_A.Cells[J,K] := FloatToStr(D);
      _Matrix_MLTPL[K,J] := D;
    end;
  end;
end;

procedure TForm1.Btn_Inverse_Click(Sender: TObject);
begin
(*
  // Inversing ....
  _Matrix_INV[..., ...] := ...;
*)
end;

initialization
  {$I Unit1.lrs}

end.

[FlamingClaw]

the inverse of a 2x2 matrix
http://www.maths.abdn.ac.uk/~igc/tch/eg1006/notes/node119.html
video about the inverse of matrix
http://www.5min.com/Video/Learn-about-Inverse-Matrix---part-1-99171123
program that shows the solution for inverse of 2x2,3x3,or more matrix
http://www.fileheap.com/authors/dekovsoft-com-18103.html

[finalist]

I found many links with samples, how to invert 2x2, 3x3 or #x# Matrix ....
These explanations are not solution for me, sorry !!!!
I need Mathematical formulas, that can be programmed on Pascal, or Pascal code ....

[FlamingClaw]

I read about this problem,some places says that the solution is : Gaus Jordan elimination
And when I search the web about it in pascal ,I found a right place,pascal sourcefiles....(sorry this problem is out of my range)...there is a man named Alan R. Miller who wrote a lot of pascal source about these problems http://infohost.nmt.edu/~armiller/pascal.htm
Hope this helps....

[pogson]

In the old days SWAG (Software Archive Group) was active and their collection is still out there.

see http://www.bsdg.org/swag/MATH/0080.PAS.html

That unit has a routine for inversion but I note that it does not check for division by zero or do pivotting to optimize the result but it may be useful as a start.

The mathematics of inversion is that reducing to the unit matrix of the subject matrix accomplishes the reverse process, turning the unit matrix into the inverse of the subject matrix, so you can just augment the subject matrix with the unit matrix of the same size and crunch away with the Gaussian elimination method.

[finalist]

pogson,
Sorry but:
procedure matinv(a:matrix;var ainv:matrix;n:integer);
does not invert correctly, the result matrix ainv is equal to a matrix;

[pogson]

pogson,
Sorry but:
procedure matinv(a:matrix;var ainv:matrix;n:integer);
does not invert correctly, the result matrix ainv is equal to a matrix;

Works for me. Maybe you inverted the identity matrix. ;-)

[B][I]sh-3.2$ ./test
1.000  2.000  4.000  
2.000  3.000  5.000  
3.000  4.000  4.000  

-4.000  4.000  -1.000  
3.500  -4.000  1.500  
-0.500  1.000  -0.500  

1.000  0.000  0.000  
0.000  1.000  0.000  
-0.000  0.000  1.000  
sh-3.2$ cat test.pas
program test;
uses algebra;
var a,a1,c:matrix;i,j:integer;
begin
 a[1,1]:=1.0;
a[2,1]:=2.0;
a[3,1]:=3.0;
a[1,2]:=2.0;
a[2,2]:=3.0;
a[3,2]:=4.0;
a[1,3]:=4.0;
a[2,3]:=5.0;
a[3,3]:=4.0;

matinv(a,a1,3);
for i:= 1 to 3 do begin for j:= 1 to 3 do begin write(a[i,j]:5:3,'  ') end; writeln end;
writeln;
for i:= 1 to 3 do begin for j:= 1 to 3 do begin write(a1[i,j]:5:3,'  ') end; writeln end;
writeln;
matmult(a,a1,c,3);

for i:= 1 to 3 do begin for j:= 1 to 3 do begin write(c[i,j]:5:3,'  ') end; writeln end;
end.
[/I][/B]

[finalist]

pogson,
Sorry it was a hidden problem of another procedure before conversion ....
The procedure before conversion did requare the matrix as var parameter and did change it inside itself, so I did send to conversion changed matrix and the result was confused ....
I did correct the problem, and send the matrix as normal parameter, it does not change inside previous procedure and after that goes to conversion in right content ....
Thanks !