Source code in the seed7 programming language

$ include "seed7_05.s7i";
  include "float.s7i";

const type: matrix is array array float;
 
const func string: str (in matrix: mat) is func
  result
    var string: stri is "";
  local
    var integer: row is 0;
    var integer: column is 0;
  begin
    for row range 1 to length(mat) do
      for column range 1 to length(mat[row]) do
        stri &:= str(mat[row][column]);
        if column < length(mat[row]) then
          stri &:= ", ";
        end if;
      end for;
      if row < length(mat) then
        stri &:= "\n";
      end if;
    end for;
  end func;

enable_output(matrix);

const func matrix: (in matrix: mat1) * (in matrix: mat2) is func
  result
    var matrix: product is matrix.value;
  local
    var integer: row is 0;
    var integer: column is 0;
    var integer: k is 0;
  begin
    product := length(mat1) times length(mat1) times 0.0;
    for row range 1 to length(mat1) do
      for column range 1 to length(mat1) do
        product[row][column] := 0.0;
        for k range 1 to length(mat1) do
          product[row][column] +:= mat1[row][k] * mat2[k][column];
        end for;
      end for;
    end for;
  end func;

const func matrix: (in var matrix: base) ** (in var integer: exponent) is func
  result
    var matrix: power is matrix.value;
  local
    var integer: row is 0;
    var integer: column is 0;
  begin
    if exponent < 0 then
      raise NUMERIC_ERROR;
    else
      if odd(exponent) then
        power := base;
      else
        # Create identity matrix
        power := length(base) times length(base) times 0.0;
        for row range 1 to length(base) do
          for column range 1 to length(base) do
            if row = column then
              power[row][column] := 1.0;
            end if;
          end for;
        end for;
      end if;
      exponent := exponent div 2;
      while exponent > 0 do
        base := base * base;
        if odd(exponent) then
          power := power * base;
        end if;
        exponent := exponent div 2;
      end while;
    end if;
  end func;

const proc: main is func
  local
    var matrix: m is [] (
      [] (4.0, 3.0),
      [] (2.0, 1.0));
    var integer: exponent is 0;
  begin
    for exponent range [] (0, 1, 2, 3, 5, 7, 11, 13, 17, 19, 23) do
      writeln("m ** " <& exponent <& " =");
      writeln(m ** exponent);
    end for;
  end func;