Посчитать матрицы по заданой формуле - Free Pascal

{$mode objfpc}
 
const
  size = 4;
 
type
  TMatrix = array[1 .. size, 1 .. size] of double;
 
operator * (const a, b: TMatrix) m: TMatrix;
var i, j, k: integer;
begin
  for i := 1 to size do
    for j := 1 to size do begin
      m[i, j] := 0;
      for k := 1 to size do
        m[i, j] := m[i, j] + a[i, k] * b[k, j]
    end;
end;
 
operator + (const a, b: TMatrix) m: TMatrix;
var i, j: integer;
begin
  for i := 1 to size do
    for j := 1 to size do m[i, j] := a[i, j] * b[i, j];
end;
 
operator ** (const a: TMatrix; const pow: integer) m: TMatrix;
var i, j: Integer;
begin
  if pow = 0 then begin
    for i := 1 to size do
      for j := 1 to size do
        m[i, j] := Byte(i = j);
    exit
  end;
 
  m := a;
  for i := 1 to pred(pow) do
    m := m * a;
end;
 
procedure matrixInput(const s : string; var m : TMatrix);
var i, j : Integer;
begin
  writeln('Input matrix ' + s);
  for i := 1 to size do
    for j := 1 to size do
    begin
      write(s+'[', i, ',', j, '] = ');
      readln(m[i, j]);
    end;
end;
 
procedure matrixPrint(a: TMatrix);
var i, j: integer;
begin
  for i := 1 to size do begin
    for j := 1 to size do
      write(a[i, j]:9:2);
    writeln
  end
end;
 
var
  A, B, C : TMatrix;
  Res: TMatrix;
 
begin
  matrixInput('A', A);
  matrixInput('B', B);
  matrixInput('C', C);
  writeln;
  Res := A**2 + B**2 + A*B*C;
  matrixPrint(Res)
end.