{$mode fpc}
{$modeswitch result}
unit dmatrix;

interface

type
generic TMatrix<TElement, TArray1, TArray2> = object
public type
  TMatrixPtr = ^TMatrix;
public var
  M: TArray2;
  constructor Init; overload;
  constructor Init(V: TElement); overload;
  function Width: Integer; inline;
  function Height: Integer; inline;
  function GetElement(I, J: Integer): TElement; inline;
  procedure SetElement(I, J: Integer; Element: TElement); inline;
  function Tr: TElement;
  function Det: TElement;
  function Mul(constref Other: TMatrix): TMatrix;
  procedure SwapColumns(I1, I2: Integer);
  procedure SwapRows(J1, J2: Integer);
  property Elements[I, J: Integer]: TElement read GetElement write SetElement; default;
end;

implementation

constructor TMatrix.Init;
begin
  Init(1);
end;

constructor TMatrix.Init(V: TElement);
var
  I, J: Integer;
begin
  for I := 0 to Width - 1 do
    for J := 0 to Height - 1 do
      if I = J then
        TArray1(M[I])[J] := V
      else
        TArray1(M[I])[J] := 0;
end;

function TMatrix.Width: Integer;
begin
  Result := SizeOf(TArray1) div SizeOf(TElement);
end;

function TMatrix.Height: Integer;
begin
  Result := SizeOf(TArray2) div (SizeOf(TArray1) * SizeOf(TElement));
end;

function TMatrix.GetElement(I, J: Integer): TElement;
begin
  Result := TArray1(M[I])[J];
end;

procedure TMatrix.SetElement(I, J: Integer; Element: TElement);
begin
  TArray1(M[I])[J] := Element;
end;

function TMatrix.Tr: TElement;
var
  I, J: Integer;
begin
  Result := 0;
  for I := 0 to Width - 1 do
    if I < Height then
      Result := Result + TArray1(M[I])[J];
end;

function TMatrix.Mul(constref Other: TMatrix): TMatrix;
var
  I, J, K: Integer;
begin
  // When FreePascal will support generic functions, parameter Other may
  // have other width and Height than Self
  Assert(Width = Other.Height);

  for I := 0 to Other.Width - 1 do
    for J := 0 to Height - 1 do begin
      Result[I, J] := 0;
      for K := 0 to Width - 1 do
        Result[I, J] := Result[I, J] + Self[K, J] * Other[I, K];
    end;
end;

procedure TMatrix.SwapColumns(I1, I2: Integer);
var
  J: Integer;
  E: TElement;
begin
  for J := 0 to Height - 1 do begin
    E := Elements[I1, J];
    Elements[I1, J] := Elements[I2, J];
    Elements[I2, J] := E;
  end;
end;

procedure TMatrix.SwapRows(J1, J2: Integer);
var
  I: Integer;
  E: TElement;
begin
  for I := 0 to Width - 1 do begin
    E := Elements[I, J1];
    Elements[I, J1] := Elements[I, J2];
    Elements[I, J2] := E;
  end;
end;

end.
