NAME=20
DESCRIPTION=20
SYNOPSIS=20
- Constructor=20 Methods And Such=20
- Matrix=20 Row, Column and Element operations=20
- Matrix=20 Operations=20
- Arithmetic=20 Operations=20
- Boolean=20 Matrix Operations=20
- Eigensystems=20
- Miscellaneous=20
OVERLOADED=20 OPERATORS=20
- SYNOPSIS=20
- DESCRIPTION=20
SEE=20 ALSO=20
VERSION=20
AUTHORS=20
CREDITS=20
COPYRIGHT=20
LICENSE=20 AGREEMENT

        for(my $i =3D 0;$i =
< $rows; $i++ ){
                $big =3D abs($matrix->[0][$i][0]);
                for(my $j=3D0;$j < $cols; $j++ ){
                        $big =3D $big < abs($matrix->[0][$i][$j])=20
                                ? abs($matrix->[0][$i][$j]) : $big;
                }
                next unless $big;
                # now $big is biggest element in row
                for(my $j =3D 0;$j < $cols; $j++ ){
                        $matrix->[0][$i][$j]  /=3D $big;
                }
                # now all elements are between [-1,1] and=20
                # dependence can be easily checked
        }

        for(my $i =3D 0;$i < $rows; $i++ ){
                for(my $k =3D 0; $k < $rows; $k++ ){
                        next if ($i =3D=3D $k);
                        print "i,k =3D $i,$k\n";
                        $tmprow =3D $matrix->row($i+1) + =
$matrix->row($k+1);
                        print "tmprow=3D $tmprow\n";
                        print "tmprow has " . =
$tmprow->_count_zeroes_row . " zeroes.\n";
                        print "i+1 has " . =
$matrix->row($i+1)->_count_zeroes_row . " zeroes\n";
                        print "k+1 has " . =
$matrix->row($k+1)->_count_zeroes_row . " zeroes\n";

                        $tmprow =3D $tmprow->_normalize_row;     =20
                        if( $tmprow->norm_sum =3D=3D 0 ){
                                for(my $j=3D0;$j<$cols;$j++){
                                        $matrix->[0][$k][$j] =3D 0;
                                }
                        } elsif ( $tmprow->_count_zeroes_row > =
$matrix->row($i+1)->_count_zeroes_row ){
                                print "tmprow has more zeroes, replacing =
row " . ($i+1) . "\n";
                                ## assign row  =20
                                $matrix->assign_row($i+1,$tmprow);
                        } elsif ( $tmprow->_count_zeroes_row > =
$matrix->row($k+1)->_count_zeroes_row ){
                                print "tmprow has more zeroes, replacing =
row " . ($k+1) . "\n";
                                $matrix->assign_row($k+1,$tmprow);
                        }      =20
                                       =20
                }
        }
        return $matrix;
}

NAME

Math::MatrixReal - Matrix of Reals

Implements the data type "matrix of reals" (and consequently also = "vector of=20 reals").

DESCRIPTION

Implements the data type "matrix of reals", which can be used almost = like any=20 other basic Perl type thanks to OPERATOR OVERLOADING, = i.e.,

  $product =3D $matrix1 * $matrix2;

does what you would like it to do (a matrix multiplication).

Also features many important operations and methods: matrix norm, = matrix=20 transposition, matrix inverse, determinant of a matrix, order and = numerical=20 condition of a matrix, scalar product of vectors, vector product of = vectors,=20 vector length, projection of row and column vectors, a comfortable way = for=20 reading in a matrix from a file, the keyboard or your code, and many = more.

Allows to solve linear equation systems using an efficient algorithm = known as=20 "L-R-decomposition" and several approximative (iterative) methods.

Features an implementation of Kleene's algorithm to compute the = minimal costs=20 for all paths in a graph with weighted edges (the "weights" being the = costs=20 associated with each edge).

SYNOPSIS

Constructor Methods And = Such

use Math::MatrixReal;=20
Makes the methods and overloaded operators of this module available = to your=20 program.
$new_matrix =3D new Math::MatrixReal($rows,$columns); =
The matrix object constructor method. A new matrix of size $rows by = $columns will be created, with the value 0.0 for all=20 elements.

Note that this method is implicitly called by many of the other = methods in=20 this module.
$new_matrix =3D=20 $some_matrix->new($rows,$columns);=20
Another way of calling the matrix object constructor method.

Matrix "$some_matrix" is not changed by this in any = way.
$new_matrix =3D $matrix->new_from_cols( [=20 $column_vector|$array_ref|$string, ... ] )=20
Creates a new matrix given a reference to an array of any of the=20 following:
- column vectors ( n by 1 Math::MatrixReal matrices )=20
- references to arrays=20
- strings properly formatted to create a column with = Math::MatrixReal's=20 new_from_string command
You may mix and match these as you wish. However, all must be of = the same=20 dimension--no padding happens automatically. Example:
```
        =
my $matrix =3D Math::MatrixReal->new_from_cols( [ [1,2], [3,4] ] );
        print $matrix;
```
will print
```
        [  1.000000000000E+00  =
3.000000000000E+00 ]
        [  2.000000000000E+00  4.000000000000E+00 ]
```
new_from_rows( [ $row_vector|$array_ref|$string, ... ] = )=20
Creates a new matrix given a reference to an array of any of the=20 following:
- row vectors ( 1 by n Math::MatrixReal matrices )=20
- references to arrays=20
- strings properly formatted to create a row with = Math::MatrixReal's=20 new_from_string command
You may mix and match these as you wish. However, all must be of = the same=20 dimension--no padding happens automatically. Example:
```
        =
my $matrix =3D Math::MatrixReal->new_from_rows( [ [1,2], [3,4] ] );
        print $matrix;
```
will print
```
        [  1.000000000000E+00  =
2.000000000000E+00 ]
        [  3.000000000000E+00  4.000000000000E+00 ]
```

$new_matrix =3D Math::MatrixReal->new_diag( $array_ref =
);

=20

This method allows you to create a diagonal matrix by only = specifying the=20 diagonal elements. Example:

        $matrix =3D =
Math::MatrixReal->new_diag( [ 1,2,3,4 ] );
        print $matrix;

will print

        [  1.000000000000E+00  =
0.000000000000E+00  0.000000000000E+00  0.000000000000E+00 ]
        [  0.000000000000E+00  2.000000000000E+00  0.000000000000E+00  =
0.000000000000E+00 ]
        [  0.000000000000E+00  0.000000000000E+00  3.000000000000E+00  =
0.000000000000E+00 ]
        [  0.000000000000E+00  0.000000000000E+00  0.000000000000E+00  =
4.000000000000E+00 ]

$new_matrix =3D=20 Math::MatrixReal->new_from_string($string);=20
This method allows you to read in a matrix from a string (for = instance,=20 from the keyboard, from a file or from your code).

The syntax is simple: each row must start with "[ " = and end=20 with " ]\n" ("\n" being the newline = character and=20 " " a space or tab) and contain one or more numbers, all=20 separated from each other by spaces or tabs.

Additional spaces or tabs can be added at will, but no = comments.

Examples:
```
  $string =3D "[ 1 2 3 ]\n[ 2 2 -1 ]\n[ 1 1 1 =
]\n";
  $matrix =3D Math::MatrixReal->new_from_string($string);
  print "$matrix";
```
By the way, this prints
```
  [  1.000000000000E+00  =
2.000000000000E+00  3.000000000000E+00 ]
  [  2.000000000000E+00  2.000000000000E+00 -1.000000000000E+00 ]
  [  1.000000000000E+00  1.000000000000E+00  1.000000000000E+00 ]
```
But you can also do this in a much more comfortable way using the=20 shell-like "here-document" syntax:
```
  $matrix =3D =
Math::MatrixReal->new_from_string(<<'MATRIX');
  [  1  0  0  0  0  0  1  ]
  [  0  1  0  0  0  0  0  ]
  [  0  0  1  0  0  0  0  ]
  [  0  0  0  1  0  0  0  ]
  [  0  0  0  0  1  0  0  ]
  [  0  0  0  0  0  1  0  ]
  [  1  0  0  0  0  0 -1  ]
  MATRIX
```
You can even use variables in the matrix:
```
  $c1 =3D   2  /  =
3;
  $c2 =3D  -2  /  5;
  $c3 =3D  26  /  9;

  $matrix =3D Math::MatrixReal->new_from_string(<<"MATRIX");

      [   3    2    0   ]
      [   0    3    2   ]
      [  $c1  $c2  $c3  ]

  MATRIX
```
(Remember that you may use spaces and tabs to format the matrix to = your=20 taste)

Note that this method uses exactly the same representation for a = matrix as=20 the "stringify" operator "": this means that you can convert any = matrix into a=20 string with $string =3D "$matrix"; and read it back in = later (for=20 instance from a file!).

Note however that you may suffer a precision loss in this process = because=20 only 13 digits are supported in the mantissa when printed!!

If the string you supply (or someone else supplies) does not obey = the=20 syntax mentioned above, an exception is raised, which can be caught by = "eval"=20 as follows:
```
  print "Please enter your matrix (in one line): =
";
  $string =3D <STDIN>;
  $string =3D~ s/\\n/\n/g;
  eval { $matrix =3D Math::MatrixReal->new_from_string($string); };
  if ($@)
  {
      print "$@";
      # ...
      # (error handling)
  }
  else
  {
      # continue...
  }
```
or as follows:
```
  eval { $matrix =3D =
Math::MatrixReal->new_from_string(<<"MATRIX"); };
  [   3    2    0   ]
  [   0    3    2   ]
  [  $c1  $c2  $c3  ]
  MATRIX
  if ($@)
  # ...
```
Actually, the method shown above for reading a matrix from the = keyboard is=20 a little awkward, since you have to enter a lot of "\n"'s for the=20 newlines.

A better way is shown in this piece of code:
```
  while (1)
  {
      print "\nPlease enter your matrix ";
      print "(multiple lines, <ctrl-D> =3D done):\n";
      eval { $new_matrix =3D
          Math::MatrixReal->new_from_string(join('',<STDIN>)); =
};
      if ($@)
      {
          $@ =3D~ s/\s+at\b.*?$//;
          print "${@}Please try again.\n";
      }
      else { last; }
  }
```
Possible error messages of the "new_from_string()" method = are:
```
  Math::MatrixReal::new_from_string(): syntax error in =
input string
  Math::MatrixReal::new_from_string(): empty input string
```
If the input string has rows with varying numbers of columns, the = following=20 warning will be printed to STDERR:
```
  =
Math::MatrixReal::new_from_string(): missing elements will be set to =
zero!
```
If everything is okay, the method returns an object reference to = the (newly=20 allocated) matrix containing the elements you specified.
$new_matrix =3D $some_matrix->shadow();=20
Returns an object reference to a NEW but EMPTY matrix = (filled=20 with zero's) of the SAME SIZE as matrix=20 "$some_matrix".

Matrix "$some_matrix" is not changed by this in any = way.
$matrix1->copy($matrix2);=20
Copies the contents of matrix "$matrix2" to an = ALREADY=20 EXISTING matrix "$matrix1" (which must have the same = size as=20 matrix "$matrix2"!).

Matrix "$matrix2" is not changed by this in any = way.
$twin_matrix =3D $some_matrix->clone();=20
Returns an object reference to a NEW matrix of the SAME = SIZE=20 as matrix "$some_matrix". The contents of matrix=20 "$some_matrix" have ALREADY BEEN COPIED to the new = matrix=20 "$twin_matrix". This is the method that the operator = "=3D" is=20 overloaded to when you type $a =3D $b, when = $a and=20 $b are matrices.

Matrix "$some_matrix" is not changed by this in any=20 way.

Matrix Row, = Column and=20 Element operations

$row_vector =3D $matrix->row($row);=20
This is a projection method which returns an object reference to a=20 NEW matrix (which in fact is a (row) vector since it has only = one row)=20 to which row number "$row" of matrix = "$matrix" has=20 already been copied.

Matrix "$matrix" is not changed by this in any = way.
$column_vector =3D $matrix->column($column);=20
This is a projection method which returns an object reference to a=20 NEW matrix (which in fact is a (column) vector since it has = only one=20 column) to which column number "$column" of matrix=20 "$matrix" has already been copied.

Matrix "$matrix" is not changed by this in any = way.
$matrix->assign($row,$column,$value);=20
Explicitly assigns a value "$value" to a single = element of the=20 matrix "$matrix", located in row "$row" and = column=20 "$column", thereby replacing the value previously stored=20 there.
$value =3D = $matrix->element($row,$column);=20
Returns the value of a specific element of the matrix=20 "$matrix", located in row "$row" and column=20 "$column".
$new_matrix =3D $matrix->each( \&function );=20
Creates a new matrix by evaluating a code reference on each element = of the=20 given matrix. The function is passed the element, the row index and = the column=20 index, in that order. The value the function returns ( or the value of = the=20 last executed statement ) is the value given to the corresponding = element in=20 $new_matrix.

Example:
```
        # add 1 to every element in the matrix
        $matrix =3D $matrix->each ( sub { (shift) + 1 } );
```
Example:
```
        my $cofactor =3D $matrix->each( sub { =
my(undef,$i,$j) =3D @_;
                ($i+$j) % 2 =3D=3D 0 ? =
$matrix->minor($i,$j)->det()
                : -1*$matrix->minor($i,$j)->det();
                } );
```
This code needs some explanation. For each element of $matrix, it = throws=20 away the actual value and stores the row and column indexes in $i and = $j. Then=20 it sets element [$i,$j] in $cofactor to the determinant of=20 $matrix->minor($i,$j) if it is an "even" element, or=20 -1*$matrix->minor($i,$j) if it is an "odd" = element.
$new_matrix =3D $matrix->each_diag( \&function = );=20
Creates a new matrix by evaluating a code reference on each = diagonal=20 element of the given matrix. The function is passed the element, the = row index=20 and the column index, in that order. The value the function returns ( = or the=20 value of the last executed statement ) is the value given to the = corresponding=20 element in $new_matrix.
$matrix->swap_col( $col1, $col2 );=20
This method takes two one-based column numbers and swaps the values = of each=20 element in each column. $matrix->swap_col(2,3) would = replace=20 column 2 in $matrix with column 3, and replace column 3 with column = 2.
$matrix->swap_row( $row1, $row2 );=20
This method takes two one-based row numbers and swaps the values of = each=20 element in each row. $matrix->swap_row(2,3) would = replace row=20 2 in $matrix with row 3, and replace row 3 with row 2.

Matrix Operations

$det =3D $matrix->det();=20
Returns the determinant of the matrix, without going through the = rigamarole=20 of computing a LR decomposition. This method should be much faster = than LR=20 decomposition if the matrix is diagonal or triangular. Otherwise, it = is just a=20 wrapper for $matrix->decompose_LR->det_LR. If the=20 determinant is zero, there is no inverse and vice-versa. Only = quadratic=20 matrices have determinants.
$inverse =3D $matrix->inverse();=20
Returns the inverse of a matrix, without going through the = rigamarole of=20 computing a LR decomposition. If no inverse exists, undef is returned = and an=20 error is printed via carp(). This is nothing but a = wrapper for=20 $matrix->decompose_LR->invert_LR.
($rows,$columns) =3D $matrix->dim();=20
Returns a list of two items, representing the number of rows and = columns=20 the given matrix "$matrix" contains.
$norm_one =3D $matrix->norm_one();=20
Returns the "one"-norm of the given matrix = "$matrix".

The "one"-norm is defined as follows:

For each column, the sum of the absolute values of the elements in = the=20 different rows of that column is calculated. Finally, the maximum of = these=20 sums is returned.

Note that the "one"-norm and the "maximum"-norm are mathematically=20 equivalent, although for the same matrix they usually yield a = different=20 value.

Therefore, you should only compare values that have been calculated = using=20 the same norm!

Throughout this package, the "one"-norm is (arbitrarily) used for = all=20 comparisons, for the sake of uniformity and comparability, except for = the=20 iterative methods "solve_GSM()", "solve_SSM()" and "solve_RM()" which = use=20 either norm depending on the matrix itself.
$norm_max =3D $matrix->norm_max();=20
Returns the "maximum"-norm of the given matrix $matrix.

The "maximum"-norm is defined as follows:

For each row, the sum of the absolute values of the elements in the = different columns of that row is calculated. Finally, the maximum of = these=20 sums is returned.

Note that the "maximum"-norm and the "one"-norm are mathematically=20 equivalent, although for the same matrix they usually yield a = different=20 value.

Therefore, you should only compare values that have been calculated = using=20 the same norm!

Throughout this package, the "one"-norm is (arbitrarily) used for = all=20 comparisons, for the sake of uniformity and comparability, except for = the=20 iterative methods "solve_GSM()", "solve_SSM()" and "solve_RM()" which = use=20 either norm depending on the matrix itself.
$norm_sum =3D $matrix->norm_sum();=20
This is a very simple norm which is defined as the sum of the = absolute=20 values of every element.
$p_norm =3D $matrix->norm_p($n);>=20
This function returns the "p-norm" of a vector. The argument $n = must be a=20 number greater than or equal to 1 or the string "Inf". The p-norm is = defined=20 as (sum(x_i^p))^(1/p). In words, it raised each element to the p-th = power,=20 adds them up, and then takes the p-th root of that number. If the = string "Inf"=20 is passed, the "infinity-norm" is computed, which is really the limit = of the=20 p-norm as p goes to infinity. It is defined as the maximum element of = the=20 vector. Also, note that the familiar Euclidean distance between two = vectors is=20 just a special case of a p-norm, when p is equal to 2.

Example: $a =3D Math::MatrixReal->new_from_cols([[1,2,3]]); $p1 = =3D=20 $a->norm_p(1); $p2 =3D $a->norm_p(2); $p3 =3D $a->norm_p(3); = $pinf =3D=20 $a->norm_p("Inf");
```
        print "(1,2,3,Inf) =
norm:\n$p1\n$p2\n$p3\n$pinf\n";

        $i1 =3D $a->new_from_rows([[1,0]]);
        $i2 =3D $a->new_from_rows([[0,1]]);

        # this should be sqrt(2) since it is the same as the=20
        # hypotenuse of a 1 by 1 right triangle

        $dist  =3D ($i1-$i2)->norm_p(2);
        print "Distance is $dist, which should be " . sqrt(2) . =
"\n";
```
Output:
```
        (1,2,3,Inf) norm:
        6
        3.74165738677394139
        3.30192724889462668
        3

        Distance is 1.41421356237309505, which should be =
1.41421356237309505
```
$frob_norm =3D = $matrix->norm_frobenius();=20
This norm is similar to that of a p-norm where p is 2, except it = acts on a=20 matrix, not a vector. Each element of the matrix is squared, = this is=20 added up, and then a square root is taken.
$matrix->spectral_radius();=20
Returns the maximum value of the absolute value of all eigenvalues. = Currently this computes all eigenvalues, then sifts through = them to=20 find the largest in absolute value. Needless to say, this is very = inefficient,=20 and in the future an algorithm that computes only the largest = eigenvalue may=20 be implemented.
$matrix1->transpose($matrix2);=20
Calculates the transposed matrix of matrix $matrix2 and stores the = result=20 in matrix "$matrix1" (which must already exist and have = the same=20 size as matrix "$matrix2"!).

This operation can also be carried out "in-place", i.e., input and = output=20 matrix may be identical.

Transposition is a symmetry operation: imagine you rotate the = matrix along=20 the axis of its main diagonal (going through elements (1,1), (2,2), = (3,3) and=20 so on) by 180 degrees.

Another way of looking at it is to say that rows and columns are = swapped.=20 In fact the contents of element (i,j) are swapped with = those of=20 element (j,i).

Note that (especially for vectors) it makes a big difference if you = have a=20 row vector, like this:
```
  [ -1 0 1 ]
```
or a column vector, like this:
```
  [ -1 ]
  [  0 ]
  [  1 ]
```
the one vector being the transposed of the other!

This is especially true for the matrix product of two = vectors:
```
               [ -1 ]
  [ -1 0 1 ] * [  0 ]  =3D  [ 2 ] ,  whereas
               [  1 ]

                             *     [ -1  0  1 ]
  [ -1 ]                                            [  1  0 -1 ]
  [  0 ] * [ -1 0 1 ]  =3D  [ -1 ]   [  1  0 -1 ]  =3D  [  0  0  0 ]
  [  1 ]                  [  0 ]   [  0  0  0 ]     [ -1  0  1 ]
                          [  1 ]   [ -1  0  1 ]
```
So be careful about what you really mean!

Hint: throughout this module, whenever a vector is explicitly = required for=20 input, a COLUMN vector is expected!
$trace =3D $matrix->trace();=20
This returns the trace of the matrix, which is defined as the sum = of the=20 diagonal elements. The matrix must be quadratic.
$minor =3D $matrix->minor($row,$col);=20
Returns the minor matrix corresponding to $row and $col. $matrix = must be=20 quadratic. If $matrix is n rows by n cols, the minor of $row and $col = will be=20 an (n-1) by (n-1) matrix. The minor is defined as crossing out the row = and the=20 col specified and returning the remaining rows and columns as a = matrix. This=20 method is used by cofactor().
$cofactor =3D $matrix->cofactor();=20
The cofactor matrix is constructed as follows:

For each element, cross out the row and column that it sits in. = Now, take=20 the determinant of the matrix that is left in the other rows and = columns.=20 Multiply the determinant by (-1)^(i+j), where i is the row index, and = j is the=20 column index. Replace the given element with this value.

The cofactor matrix can be used to find the inverse of the matrix. = One=20 formula for the inverse of a matrix is the cofactor matrix transposed = divided=20 by the original determinant of the matrix.

The following two inverses should be exactly the same:
```
     =
   my $inverse1 =3D $matrix->inverse;
        my $inverse2 =3D ~($matrix->cofactor)->each( sub { =
(shift)/$matrix->det() } );
```
Caveat: Although the cofactor matrix is simple algorithm to compute = the=20 inverse of a matrix, and can be used with pencil and paper for small = matrices,=20 it is comically slower than the native inverse() = function. Here=20 is a small benchmark:
```
        # $matrix1 is 15x15
        $det =3D $matrix1->det;
        timethese( 10,
        {'inverse' =3D> sub { $matrix1->inverse(); },
          'cofactor' =3D> sub { (~$matrix1->cofactor)->each ( =
sub { (shift)/$det; } ) }
        } );


        Benchmark: timing 10 iterations of LR, cofactor, inverse...
        inverse:  1 wallclock secs ( 0.56 usr +  0.00 sys =3D  0.56 CPU) =
@ 17.86/s (n=3D10)
        cofactor: 36 wallclock secs (36.62 usr +  0.01 sys =3D 36.63 =
CPU) @  0.27/s (n=3D10)
```
$adjoint =3D $matrix->adjoint();=20
The adjoint is just the transpose of the cofactor matrix. This = method is=20 just an alias for ~($matrix->cofactor).

Arithmetic Operations

$matrix1->add($matrix2,$matrix3);=20
Calculates the sum of matrix "$matrix2" and matrix=20 "$matrix3" and stores the result in matrix=20 "$matrix1" (which must already exist and have the same = size as=20 matrix "$matrix2" and matrix = "$matrix3"!).

This operation can also be carried out "in-place", i.e., the output = and one=20 (or both) of the input matrices may be identical.
$matrix1->subtract($matrix2,$matrix3);=20
Calculates the difference of matrix "$matrix2" minus = matrix=20 "$matrix3" and stores the result in matrix=20 "$matrix1" (which must already exist and have the same = size as=20 matrix "$matrix2" and matrix = "$matrix3"!).

This operation can also be carried out "in-place", i.e., the output = and one=20 (or both) of the input matrices may be identical.

Note that this operation is the same as=20 $matrix1->add($matrix2,-$matrix3);, although the = latter is a=20 little less efficient.
$matrix1->multiply_scalar($matrix2,$scalar);=20
Calculates the product of matrix "$matrix2" and the = number=20 "$scalar" (i.e., multiplies each element of matrix=20 "$matrix2" with the factor "$scalar") and = stores the=20 result in matrix "$matrix1" (which must already exist and = have=20 the same size as matrix "$matrix2"!).

This operation can also be carried out "in-place", i.e., input and = output=20 matrix may be identical.
$product_matrix =3D $matrix1->multiply($matrix2);=20
Calculates the product of matrix "$matrix1" and matrix = "$matrix2" and returns an object reference to a new = matrix=20 "$product_matrix" in which the result of this operation = has been=20 stored.

Note that the dimensions of the two matrices = "$matrix1" and=20 "$matrix2" (i.e., their numbers of rows and columns) must = harmonize in the following way (example):
```
                     =
     [ 2 2 ]
                          [ 2 2 ]
                          [ 2 2 ]

              [ 1 1 1 ]   [ * * ]
              [ 1 1 1 ]   [ * * ]
              [ 1 1 1 ]   [ * * ]
              [ 1 1 1 ]   [ * * ]
```
I.e., the number of columns of matrix "$matrix1" has = to be the=20 same as the number of rows of matrix "$matrix2".

The number of rows and columns of the resulting matrix=20 "$product_matrix" is determined by the number of rows of = matrix=20 "$matrix1" and the number of columns of matrix=20 "$matrix2", respectively.
$matrix1->negate($matrix2);=20
Calculates the negative of matrix "$matrix2" (i.e., = multiplies=20 all elements with "-1") and stores the result in matrix=20 "$matrix1" (which must already exist and have the same = size as=20 matrix "$matrix2"!).

This operation can also be carried out "in-place", i.e., input and = output=20 matrix may be identical.
$matrix_to_power =3D $matrix1->exponent($integer); =
Raises the matrix to the $integer power. Obviously,=20 $integer must be an integer. If it is zero, the identity = matrix=20 is returned. If a negative integer is given, the inverse will be = computed (if=20 it exists) and then raised the the absolute value of = $integer.=20 The matrix must be quadratic.

Boolean Matrix = Operations

$matrix->is_quadratic();=20
Returns a boolean value indicating if the given matrix is quadratic = (also=20 know as "square" or "n by n"). A matrix is quadratic if it has the = same number=20 of rows as it does columns.
$matrix->is_square();=20
This is an alias for is_quadratic().
$matrix->is_symmetric();=20
Returns a boolean value indicating if the given matrix is = symmetric. By=20 definition, a matrix is symmetric if and only if=20 (M[i,j]=3DM[j,i]). This is = equivalent=20 to ($matrix =3D=3D ~$matrix) but without memory = allocation. Only=20 quadratic matrices can be symmetric.

Notes: A symmetric matrix always has real eigenvalues/eigenvectors. = A=20 matrix plus its transpose is always symmetric.
$matrix->is_skew_symmetric();=20
Returns a boolean value indicating if the given matrix is skew = symmetric.=20 By definition, a matrix is symmetric if and only if=20 (M[i,j]=3D-M[j,i]). This is = equivalent=20 to ($matrix =3D=3D -(~$matrix)) but without memory = allocation. Only=20 quadratic matrices can be skew symmetric.
$matrix->is_diagonal();=20
Returns a boolean value indicating if the given matrix is diagonal, = i.e.=20 all of the nonzero elements are on the main diagonal. Only quadratic = matrices=20 can be diagonal.
$matrix->is_tridiagonal();=20
Returns a boolean value indicating if the given matrix is = tridiagonal, i.e.=20 all of the nonzero elements are on the main diagonal or the diagonals = above=20 and below the main diagonal. Only quadratic matrices can be = tridiagonal.
$matrix->is_upper_triangular();=20
Returns a boolean value indicating if the given matrix is upper = triangular,=20 i.e. all of the nonzero elements not on the main diagonal are above = it. Only=20 quadratic matrices can be upper triangular. Note: diagonal matrices = are both=20 upper and lower triangular.
$matrix->is_lower_triangular();=20
Returns a boolean value indicating if the given matrix is lower = triangular,=20 i.e. all of the nonzero elements not on the main diagonal are below = it. Only=20 quadratic matrices can be lower triangular. Note: diagonal matrices = are both=20 upper and lower triangular.
$matrix->is_orthogonal();=20
Returns a boolean value indicating if the given matrix is = orthogonal. An=20 orthogonal matrix is has the property that the transpose equals the = inverse of=20 the matrix. Instead of computing each and comparing them, this method=20 multiplies the matrix by it's transpose, and returns true if this = turns out to=20 be the identity matrix, false otherwise. Only quadratic matrices can=20 orthogonal.
$matrix->is_binary();=20
Returns a boolean value indicating if the given matrix is binary. A = matrix=20 is binary if it contains only zeroes or ones.
$matrix->is_gramian();=20
Returns a boolean value indicating if the give matrix is Gramian. A = matrix=20 $A is Gramian if and only if there exists a square matrix = $B such that $A =3D ~$B*$B. This is = equivalent to=20 checking if $A is symmetric and has all nonnegative = eigenvalues,=20 which is what Math::MatrixReal uses to check for this property.
$matrix->is_LR();=20
Returns a boolean value indicating if the matrix is an LR = decomposition=20 matrix.
$matrix->is_positive();=20
Returns a boolean value indicating if the matrix contains only = positive=20 entries. Note that a zero entry is not positive and will cause=20 is_positive() to return false.
$matrix->is_negative();=20
Returns a boolean value indicating if the matrix contains only = negative=20 entries. Note that a zero entry is not negative and will cause=20 is_negative() to return false.
$matrix->is_periodic($k);=20
Returns a boolean value indicating if the matrix is periodic with = period=20 $k. This is true if $matrix ** ($k+1) =3D=3D $matrix. = When $k=20 =3D=3D 1, this reduces down to the is_idempotent() = function.
$matrix->is_idempotent();=20
Returns a boolean value indicating if the matrix is idempotent, = which is=20 defined as the square of the matrix being equal to the original = matrix, i.e=20 $matrix ** 2 =3D=3D $matrix.
$matrix->is_row_vector();=20
Returns a boolean value indicating if the matrix is a row vector. A = row=20 vector is a matrix which is 1xn. Note that the 1x1 matrix is both a = row and=20 column vector.
$matrix->is_col_vector();=20
Returns a boolean value indicating if the matrix is a col vector. A = col=20 vector is a matrix which is nx1. Note that the 1x1 matrix is both a = row and=20 column vector.

Eigensystems

($l, $V) =3D $matrix->sym_diagonalize();=20
This method performs the diagonalization of the quadratic = symmetric=20 matrix M stored in $matrix. On output, l is a column = vector=20 containing all the eigenvalues of M and V is an = orthogonal=20 matrix which columns are the corresponding normalized eigenvectors. = The=20 primary property of an eigenvalue l and an eigenvector x = is of=20 course that: M * x =3D l * x.

The method uses a Householder reduction to tridiagonal form = followed by a=20 QL algoritm with implicit shifts on this tridiagonal. (The tridiagonal = matrix=20 is kept internally in a compact form in this routine to save memory.) = In fact,=20 this routine wraps the householder() and tri_diagonalize() methods = described=20 below when their intermediate results are not desired. The overall = algorithmic=20 complexity of this technique is O(N^3). According to several books, = the=20 coefficient hidden by the 'O' is one of the best possible for general=20 (symmetric) matrixes.
($T, $Q) =3D $matrix->householder();=20
This method performs the Householder algorithm which reduces the = n=20 by n real symmetric matrix M contained in $matrix = to=20 tridiagonal form. On output, T is a symmetric tridiagonal = matrix (only=20 diagonal and off-diagonal elements are non-zero) and Q is an=20 orthogonal matrix performing the tranformation between M = and=20 T ($M =3D=3D $Q * $T * ~$Q).
($l, $V) =3D $T->tri_diagonalize([$Q]);=20
This method diagonalizes the symmetric tridiagonal matrix T. = On=20 output, $l and $V are similar to the output values described for=20 sym_diagonalize().

The optional argument $Q corresponds to an orthogonal = transformation matrix=20 Q that should be used additionally during V = (eigenvectors)=20 computation. It should be supplied if the desired eigenvectors = correspond to a=20 more general symmetric matrix M previously reduced by the = householder()=20 method, not a mere tridiagonal. If T is really a tridiagonal = matrix,=20 Q can be omitted (it will be internally created in fact as an = identity=20 matrix). The method uses a QL algorithm (with implicit shifts).
$l =3D $matrix->sym_eigenvalues();=20
This method computes the eigenvalues of the quadratic = symmetric=20 matrix M stored in $matrix. On output, l is a column = vector=20 containing all the eigenvalues of M. Eigenvectors are not = computed (on=20 the contrary of sym_diagonalize()) and this method is = more=20 efficient (even though it uses a similar algorithm with two phases). = However,=20 understand that the algorithmic complexity of this technique is still = also=20 O(N^3). But the coefficient hidden by the 'O' is better by a factor = of...,=20 well, see your benchmark, it's wiser.

This routine wraps the householder_tridiagonal() and = tri_eigenvalues()=20 methods described below when the intermediate tridiagonal matrix is = not=20 needed.
$T =3D $matrix->householder_tridiagonal();=20
This method performs the Householder algorithm which reduces the = n=20 by n real symmetric matrix M contained in $matrix = to=20 tridiagonal form. On output, T is the obtained symmetric = tridiagonal=20 matrix (only diagonal and off-diagonal elements are non-zero). The = operation=20 is similar to the householder() method, but potentially a little more=20 efficient as the transformation matrix is not computed.
$l =3D $T->tri_eigenvalues();=20
This method computesthe eigenvalues of the symmetric tridiagonal = matrix=20 T. On output, $l is a vector containing the eigenvalues = (similar to=20 sym_eigenvalues()). This method is much more efficient = than=20 tri_diagonalize() when eigenvectors are not needed.

Miscellaneous

$matrix->zero();=20
Assigns a zero to every element of the matrix = "$matrix", i.e.,=20 erases all values previously stored there, thereby effectively = transforming=20 the matrix into a "zero"-matrix or "null"-matrix, the neutral element = of the=20 addition operation in a Ring.

(For instance the (quadratic) matrices with "n" rows and columns = and matrix=20 addition and multiplication form a Ring. Most prominent characteristic = of a=20 Ring is that multiplication is not commutative, i.e., in general,=20 "matrix1 * matrix2" is not the same as "matrix2 *=20 matrix1"!)
$matrix->one();=20
Assigns one's to the elements on the main diagonal (elements (1,1), = (2,2),=20 (3,3) and so on) of matrix "$matrix" and zero's to all = others,=20 thereby erasing all values previously stored there and transforming = the matrix=20 into a "one"-matrix, the neutral element of the multiplication = operation in a=20 Ring.

(If the matrix is quadratic (which this method doesn't require, = though),=20 then multiplying this matrix with itself yields this same matrix = again, and=20 multiplying it with some other matrix leaves that other matrix = unchanged!)
$latex_string =3D $matrix->as_latex( align=3D> "c", = format =3D>=20 "%s", name =3D> "" );=20
This function returns the matrix as a LaTeX string. It takes a hash = as an=20 argument which is used to control the style of the output. The hash = element=20 align may be "c","l" or "r", corresponding to center, = left and=20 right, respectively. The format element is a format = string that=20 is given to sprintf to control the style of number = format, such a=20 floating point or scientific notation. The name element = can be=20 used so that a LaTeX string of "$name =3D " is prepended to the = string.

Example:
```
        my $a =3D =
Math::MatrixReal->new_from_cols([[ 1.234, 5.678, 9.1011],[1,2,3]] );
        print $a->as_latex( ( format =3D> "%.2f", align =3D> =
"l",name =3D> "A" ) );
```
Output: $A =3D $ $ \left( \begin{array}{ll} 1.23&1.00 \\ = 5.68&2.00 \\=20 9.10&3.00 \end{array} \right) $
$yacas_string =3D $matrix->as_yacas( format =3D> "%s", = name =3D>=20 "", semi =3D> 0 );=20
This function returns the matrix as a string that can be read by = Yacas. It=20 takes a hash as an an argument which controls the style of the output. = The=20 format element is a format string that is given to=20 sprintf to control the style of number format, such a = floating=20 point or scientific notation. The name element can be = used so=20 that "$name =3D " is prepended to the string. The <semi> element = can be=20 set to 1 to that a semicolon is appended (so Matlab does not print out = the=20 matrix.)

Example:
```
        $a =3D =
Math::MatrixReal->new_from_cols([[ 1.234, 5.678, 9.1011],[1,2,3]] );
        print $a->as_yacas( ( format =3D> "%.2f", align =3D> =
"l",name =3D> "A" ) );
```
Output:
```
        A :=3D =
{{1.23,1.00},{5.68,2.00},{9.10,3.00}}
```
$matlab_string =3D $matrix->as_matlab( format =3D> = "%s", name=20 =3D> "", semi =3D> 0 );=20
This function returns the matrix as a string that can be read by = Matlab. It=20 takes a hash as an an argument which controls the style of the output. = The=20 format element is a format string that is given to=20 sprintf to control the style of number format, such a = floating=20 point or scientific notation. The name element can be = used so=20 that "$name =3D " is prepended to the string. The <semi> element = can be=20 set to 1 to that a semicolon is appended (so Matlab does not print out = the=20 matrix.)

Example:
```
        my $a =3D =
Math::MatrixReal->new_from_rows([[ 1.234, 5.678, 9.1011],[1,2,3]] );
        print $a->as_matlab( ( format =3D> "%.3f", name =3D> =
"A",semi =3D> 1 ) );
```
Output: A =3D [ 1.234 5.678 9.101; 1.000 2.000 3.000];
$scilab_string =3D $matrix->as_scilab( format =3D> = "%s", name=20 =3D> "", semi =3D> 0 );=20
This function is just an alias for as_matlab(), since = both=20 Scilab and Matlab have the same matrix format.
$minimum =3D = Math::MatrixReal::min($number1,$number2);=20
Returns the minimum of the two numbers "number1" and=20 "number2".
$maximum =3D = Math::MatrixReal::max($number1,$number2);=20
Returns the maximum of the two numbers "number1" and=20 "number2".
$minimal_cost_matrix =3D $cost_matrix->kleene();=20
Copies the matrix "$cost_matrix" (which has to be = quadratic!)=20 to a new matrix of the same size (i.e., "clones" the input matrix) and = applies=20 Kleene's algorithm to it.

See Math::Kleene(3)<= /A> for=20 more details about this algorithm!

The method returns an object reference to the new matrix.

Matrix "$cost_matrix" is not changed by this method in = any=20 way.
($norm_matrix,$norm_vector) =3D = $matrix->normalize($vector);=20
This method is used to improve the numerical stability when solving = linear=20 equation systems.

Suppose you have a matrix "A" and a vector "b" and you want to find = out a=20 vector "x" so that A * x =3D b, i.e., the vector "x" = which solves=20 the equation system represented by the matrix "A" and the vector = "b".

Applying this method to the pair (A,b) yields a pair (A',b') where = each row=20 has been divided by (the absolute value of) the greatest coefficient = appearing=20 in that row. So this coefficient becomes equal to "1" (or "-1") in the = new=20 pair (A',b') (all others become smaller than one and greater than = minus=20 one).

Note that this operation does not change the equation system itself = because=20 the same division is carried out on either side of the equation = sign!

The method requires a quadratic (!) matrix "$matrix" = and a=20 vector "$vector" for input (the vector must be a column = vector=20 with the same number of rows as the input matrix) and returns a list = of two=20 items which are object references to a new matrix and a new vector, in = this=20 order.

The output matrix and vector are clones of the input matrix and = vector to=20 which the operation explained above has been applied.

The input matrix and vector are not changed by this in any way.

Example of how this method can affect the result of the methods to = solve=20 equation systems (explained immediately below following this = method):

Consider the following little program:
```
  #!perl -w

  use Math::MatrixReal qw(new_from_string);

  $A =3D Math::MatrixReal->new_from_string(<<"MATRIX");
  [  1   2   3  ]
  [  5   7  11  ]
  [ 23  19  13  ]
  MATRIX

  $b =3D Math::MatrixReal->new_from_string(<<"MATRIX");
  [   0   ]
  [   1   ]
  [  29   ]
  MATRIX

  $LR =3D $A->decompose_LR();
  if (($dim,$x,$B) =3D $LR->solve_LR($b))
  {
      $test =3D $A * $x;
      print "x =3D \n$x";
      print "A * x =3D \n$test";
  }

  ($A_,$b_) =3D $A->normalize($b);

  $LR =3D $A_->decompose_LR();
  if (($dim,$x,$B) =3D $LR->solve_LR($b_))
  {
      $test =3D $A * $x;
      print "x =3D \n$x";
      print "A * x =3D \n$test";
  }
```
This will print:
```
  x =3D
  [  1.000000000000E+00 ]
  [  1.000000000000E+00 ]
  [ -1.000000000000E+00 ]
  A * x =3D
  [  4.440892098501E-16 ]
  [  1.000000000000E+00 ]
  [  2.900000000000E+01 ]
  x =3D
  [  1.000000000000E+00 ]
  [  1.000000000000E+00 ]
  [ -1.000000000000E+00 ]
  A * x =3D
  [  0.000000000000E+00 ]
  [  1.000000000000E+00 ]
  [  2.900000000000E+01 ]
```
You can see that in the second example (where "normalize()" has = been used),=20 the result is "better", i.e., more accurate!
$LR_matrix =3D $matrix->decompose_LR();=20
This method is needed to solve linear equation systems.

Suppose you have a matrix "A" and a vector "b" and you want to find = out a=20 vector "x" so that A * x =3D b, i.e., the vector "x" = which solves=20 the equation system represented by the matrix "A" and the vector = "b".

You might also have a matrix "A" and a whole bunch of different = vectors=20 "b1".."bk" for which you need to find vectors "x1".."xk" so that = A * xi=20 =3D bi, for i=3D1..k.

Using Gaussian transformations (multiplying a row or column with a = factor,=20 swapping two rows or two columns and adding a multiple of one row or = column to=20 another), it is possible to decompose any matrix "A" into two = triangular=20 matrices, called "L" and "R" (for "Left" and "Right").

"L" has one's on the main diagonal (the elements (1,1), (2,2), = (3,3) and so=20 so), non-zero values to the left and below of the main diagonal and = all zero's=20 in the upper right half of the matrix.

"R" has non-zero values on the main diagonal as well as to the = right and=20 above of the main diagonal and all zero's in the lower left half of = the=20 matrix, as follows:
```
          [ 1 0 0 0 0 ]      [ x x x x x ]
          [ x 1 0 0 0 ]      [ 0 x x x x ]
      L =3D [ x x 1 0 0 ]  R =3D [ 0 0 x x x ]
          [ x x x 1 0 ]      [ 0 0 0 x x ]
          [ x x x x 1 ]      [ 0 0 0 0 x ]
```
Note that "L * R" is equivalent to matrix "A" in the = sense=20 that L * R * x =3D b <=3D=3D> A * x =3D b for all = vectors "x",=20 leaving out of account permutations of the rows and columns (these are = taken=20 care of "magically" by this module!) and numerical errors.

Trick:

Because we know that "L" has one's on its main diagonal, we can = store both=20 matrices together in the same array without information loss! = I.e.,
```
                 [ R R R R R ]
                 [ L R R R R ]
            LR =3D [ L L R R R ]
                 [ L L L R R ]
                 [ L L L L R ]
```
Beware, though, that "LR" and "L * R" are not the = same!!!

Note also that for the same reason, you cannot apply the method=20 "normalize()" to an "LR" decomposition matrix. Trying to do so will = yield=20 meaningless rubbish!

(You need to apply "normalize()" to each pair (Ai,bi) BEFORE = decomposing the matrix "Ai'"!)

Now what does all this help us in solving linear equation = systems?

It helps us because a triangular matrix is the next best thing that = can=20 happen to us besides a diagonal matrix (a matrix that has non-zero = values only=20 on its main diagonal - in which case the solution is trivial, simply = divide=20 "b[i]" by "A[i,i]" to get = "x[i]"!).

To find the solution to our problem "A * x =3D b", we = divide=20 this problem in parts: instead of solving A * x =3D b = directly, we=20 first decompose "A" into "L" and "R" and then solve "L * y =3D = b"=20 and finally "R * x =3D y" (motto: divide and rule!).

From the illustration above it is clear that solving "L * y = =3D=20 b" and "R * x =3D y" is straightforward: we = immediately know=20 that y[1] =3D b[1]. We then deduce swiftly that
```
  =
y[2] =3D b[2] - L[2,1] * y[1]
```
(and we know "y[1]" by now!), that
```
  y[3] =3D =
b[3] - L[3,1] * y[1] - L[3,2] * y[2]
```
and so on.

Having effortlessly calculated the vector "y", we now proceed to = calculate=20 the vector "x" in a similar fashion: we see immediately that = x[n] =3D y[n]=20 / R[n,n]. It follows that
```
  x[n-1] =3D ( y[n-1] - =
R[n-1,n] * x[n] ) / R[n-1,n-1]
```
and
```
  x[n-2] =3D ( y[n-2] - R[n-2,n-1] * x[n-1] - R[n-2,n] =
* x[n] )
           / R[n-2,n-2]
```
and so on.

You can see that - especially when you have many vectors "b1".."bk" = for=20 which you are searching solutions to A * xi =3D bi - this = scheme is=20 much more efficient than a straightforward, "brute force" = approach.

This method requires a quadratic matrix as its input matrix.

If you don't have that many equations, fill up with zero's (i.e., = do=20 nothing to fill the superfluous rows if it's a "fresh" matrix, i.e., a = matrix=20 that has been created with "new()" or "shadow()").

The method returns an object reference to a new matrix containing = the=20 matrices "L" and "R".

The input matrix is not changed by this method in any way.

Note that you can "copy()" or "clone()" the result of this method = without=20 losing its "magical" properties (for instance concerning the hidden=20 permutations of its rows and columns).

However, as soon as you are applying any method that alters the = contents of=20 the matrix, its "magical" properties are stripped off, and the matrix=20 immediately reverts to an "ordinary" matrix (with the values it just = happens=20 to contain at that moment, be they meaningful as an ordinary matrix or = not!).
($dimension,$x_vector,$base_matrix) =3D=20 $LR_matrix->solve_LR($b_vector);=20
Use this method to actually solve an equation system.

Matrix "$LR_matrix" must be a (quadratic) matrix = returned by=20 the method "decompose_LR()", the LR decomposition matrix of the matrix = "A" of=20 your equation system A * x =3D b.

The input vector "$b_vector" is the vector "b" in your = equation system A * x =3D b, which must be a column = vector and have=20 the same number of rows as the input matrix = "$LR_matrix".

The method returns a list of three items if a solution exists or an = empty=20 list otherwise (!).

Therefore, you should always use this method like this:
```
  =
if ( ($dim,$x_vec,$base) =3D $LR->solve_LR($b_vec) )
  {
      # do something with the solution...
  }
  else
  {
      # do something with the fact that there is no solution...
  }
```
The three items returned are: the dimension = "$dimension" of=20 the solution space (which is zero if only one solution exists, one if = the=20 solution is a straight line, two if the solution is a plane, and so = on), the=20 solution vector "$x_vector" (which is the vector "x" of = your=20 equation system A * x =3D b) and a matrix=20 "$base_matrix" representing a base of the solution space = (a set=20 of vectors which put up the solution space like the spokes of an=20 umbrella).

Only the first "$dimension" columns of this base = matrix=20 actually contain entries, the remaining columns are all zero.

Now what is all this stuff with that "base" good for?

The output vector "x" is ALWAYS a solution of your equation = system=20 A * x =3D b.

But also any vector "$vector"
```
  $vector =3D =
$x_vector->clone();

  $machine_infinity =3D 1E+99; # or something like that

  for ( $i =3D 1; $i <=3D $dimension; $i++ )
  {
      $vector +=3D rand($machine_infinity) * =
$base_matrix->column($i);
  }
```
is a solution to your problem A * x =3D b, i.e., if=20 "$A_matrix" contains your matrix "A", then
```
  =
print abs( $A_matrix * $vector - $b_vector ), "\n";
```
should print a number around 1E-16 or so!

By the way, note that you can actually calculate those vectors=20 "$vector" a little more efficient as follows:
```
  =
$rand_vector =3D $x_vector->shadow();

  $machine_infinity =3D 1E+99; # or something like that

  for ( $i =3D 1; $i <=3D $dimension; $i++ )
  {
      $rand_vector->assign($i,1, rand($machine_infinity) );
  }

  $vector =3D $x_vector + ( $base_matrix * $rand_vector );
```
Note that the input matrix and vector are not changed by this = method in any=20 way.
$inverse_matrix =3D $LR_matrix->invert_LR();=20
Use this method to calculate the inverse of a given matrix=20 "$LR_matrix", which must be a (quadratic) matrix returned = by the=20 method "decompose_LR()".

The method returns an object reference to a new matrix of the same = size as=20 the input matrix containing the inverse of the matrix that you = initially fed=20 into "decompose_LR()" IF THE INVERSE EXISTS, or an empty list=20 otherwise.

Therefore, you should always use this method in the following = way:
```
  if ( $inverse_matrix =3D $LR->invert_LR() )
  {
      # do something with the inverse matrix...
  }
  else
  {
      # do something with the fact that there is no inverse matrix...
  }
```
Note that by definition (disregarding numerical errors), the = product of the=20 initial matrix and its inverse (or vice-versa) is always a matrix = containing=20 one's on the main diagonal (elements (1,1), (2,2), (3,3) and so on) = and zero's=20 elsewhere.

The input matrix is not changed by this method in any way.
$condition =3D = $matrix->condition($inverse_matrix);=20
In fact this method is just a shortcut for
```
  abs($matrix) * =
abs($inverse_matrix)
```
Both input matrices must be quadratic and have the same size, and = the=20 result is meaningful only if one of them is the inverse of the other = (for=20 instance, as returned by the method "invert_LR()").

The number returned is a measure of the "condition" of the given = matrix=20 "$matrix", i.e., a measure of the numerical stability of = the=20 matrix.

This number is always positive, and the smaller its value, the = better the=20 condition of the matrix (the better the stability of all subsequent=20 computations carried out using this matrix).

Numerical stability means for example that if
```
  abs( =
$vec_correct - $vec_with_error ) < $epsilon
```
holds, there must be a "$delta" which doesn't depend = on the=20 vector "$vec_correct" (nor = "$vec_with_error", by the=20 way) so that
```
  abs( $matrix * $vec_correct - $matrix * =
$vec_with_error ) < $delta
```
also holds.
$determinant =3D $LR_matrix->det_LR();=20
Calculates the determinant of a matrix, whose LR decomposition = matrix=20 "$LR_matrix" must be given (which must be a (quadratic) = matrix=20 returned by the method "decompose_LR()").

In fact the determinant is a by-product of the LR decomposition: It = is (in=20 principle, that is, except for the sign) simply the product of the = elements on=20 the main diagonal (elements (1,1), (2,2), (3,3) and so on) of the LR=20 decomposition matrix.

(The sign is taken care of "magically" by this module)
$order =3D $LR_matrix->order_LR();=20
Calculates the order (called "Rang" in German) of a matrix, whose = LR=20 decomposition matrix "$LR_matrix" must be given (which = must be a=20 (quadratic) matrix returned by the method "decompose_LR()").

This number is a measure of the number of linear independent row = and column=20 vectors (=3D number of linear independent equations in the case of a = matrix=20 representing an equation system) of the matrix that was initially fed = into=20 "decompose_LR()".

If "n" is the number of rows and columns of the (quadratic!) = matrix, then=20 "n - order" is the dimension of the solution space of the associated = equation=20 system.
$rank =3D $LR_matrix->rank_LR();=20
This is an alias for the order_LR() function. The = "order" is=20 usually called the "rank" in the United States.
$scalar_product =3D = $vector1->scalar_product($vector2);=20
Returns the scalar product of vector "$vector1" and = vector=20 "$vector2".

Both vectors must be column vectors (i.e., a matrix having several = rows but=20 only one column).

This is a (more efficient!) shortcut for
```
  $temp           =
=3D ~$vector1 * $vector2;
  $scalar_product =3D  $temp->element(1,1);
```
or the sum i=3D1..n of the products vector1[i] * = vector2[i].

Provided none of the two input vectors is the null vector, then the = two=20 vectors are orthogonal, i.e., have an angle of 90 degrees between = them,=20 exactly when their scalar product is zero, and vice-versa.
$vector_product =3D = $vector1->vector_product($vector2);=20
Returns the vector product of vector "$vector1" and = vector=20 "$vector2".

Both vectors must be column vectors (i.e., a matrix having several = rows but=20 only one column).

Currently, the vector product is only defined for 3 dimensions = (i.e.,=20 vectors with 3 rows); all other vectors trigger an error message.

In 3 dimensions, the vector product of two vectors "x" and "y" is = defined=20 as
```
              |  x[1]  y[1]  e[1]  |
  determinant |  x[2]  y[2]  e[2]  |
              |  x[3]  y[3]  e[3]  |
```
where the "x[i]" and "y[i]" are the = components of=20 the two vectors "x" and "y", respectively, and the "e[i]" = are=20 unity vectors (i.e., vectors with a length equal to one) with a one in = row "i"=20 and zero's elsewhere (this means that you have numbers and vectors as = elements=20 in this matrix!).

This determinant evaluates to the rather simple formula
```
  =
z[1] =3D x[2] * y[3] - x[3] * y[2]
  z[2] =3D x[3] * y[1] - x[1] * y[3]
  z[3] =3D x[1] * y[2] - x[2] * y[1]
```
A characteristic property of the vector product is that the = resulting=20 vector is orthogonal to both of the input vectors (if neither of both = is the=20 null vector, otherwise this is trivial), i.e., the scalar product of = each of=20 the input vectors with the resulting vector is always zero.

$length =3D $vector->length();=20

This is actually a shortcut for

  $length =3D sqrt( =
$vector->scalar_product($vector) );

and returns the length of a given (column!) vector=20 "$vector".

Note that the "length" calculated by this method is in fact the = "two"-norm=20 of a vector "$vector"!

The general definition for norms of vectors is the = following:

  sub vector_norm
  {
      croak "Usage: \$norm =3D \$vector->vector_norm(\$n);"
        if (@_ !=3D 2);

      my($vector,$n) =3D @_;
      my($rows,$cols) =3D ($vector->[1],$vector->[2]);
      my($k,$comp,$sum);

      croak "Math::MatrixReal::vector_norm(): vector is not a column =
vector"
        unless ($cols =3D=3D 1);

      croak "Math::MatrixReal::vector_norm(): norm index must be > 0"
        unless ($n > 0);

      croak "Math::MatrixReal::vector_norm(): norm index must be =
integer"
        unless ($n =3D=3D int($n));

      $sum =3D 0;
      for ( $k =3D 0; $k < $rows; $k++ )
      {
          $comp =3D abs( $vector->[0][$k][0] );
          $sum +=3D $comp ** $n;
      }
      return( $sum ** (1 / $n) );
  }

Note that the case "n =3D 1" is the "one"-norm for matrices applied = to a=20 vector, the case "n =3D 2" is the euclidian norm or length of a = vector, and if=20 "n" goes to infinity, you have the "infinity"- or "maximum"-norm for = matrices=20 applied to a vector!

$xn_vector =3D=20 = $matrix->solve_GSM($x0_vector,$b_vector,$epsilon);= =20
$xn_vector =3D=20 = $matrix->solve_SSM($x0_vector,$b_vector,$epsilon);= =20
$xn_vector =3D=20 = $matrix->solve_RM($x0_vector,$b_vector,$weight,$epsilon);==20
In some cases it might not be practical or desirable to solve an = equation=20 system "A * x =3D b" using an analytical algorithm like = the=20 "decompose_LR()" and "solve_LR()" method pair.

In fact in some cases, due to the numerical properties (the = "condition") of=20 the matrix "A", the numerical error of the obtained result can be = greater than=20 by using an approximative (iterative) algorithm like one of the three=20 implemented here.

All three methods, GSM ("Global Step Method" or = "Gesamtschrittverfahren"),=20 SSM ("Single Step Method" or "Einzelschrittverfahren") and RM = ("Relaxation=20 Method" or "Relaxationsverfahren"), are fix-point iterations, that is, = can be=20 described by an iteration function "x(t+1) =3D Phi( x(t) = )" which=20 has the property:
```
  Phi(x)  =3D  x    <=3D=3D>    A * x  =
=3D  b
```
We can define "Phi(x)" as follows:
```
  Phi(x)  =
:=3D  ( En - A ) * x  +  b
```
where "En" is a matrix of the same size as "A" ("n" rows and = columns) with=20 one's on its main diagonal and zero's elsewhere.

This function has the required property.

Proof:
```
           A * x        =3D   b

  <=3D=3D>  -( A * x )      =3D  -b

  <=3D=3D>  -( A * x ) + x  =3D  -b + x

  <=3D=3D>  -( A * x ) + x + b  =3D  x

  <=3D=3D>  x - ( A * x ) + b  =3D  x

  <=3D=3D>  ( En - A ) * x + b  =3D  x
```
This last step is true because
```
  x[i] - ( a[i,1] x[1] + ... =
+ a[i,i] x[i] + ... + a[i,n] x[n] ) + b[i]
```
is the same as
```
  ( -a[i,1] x[1] + ... + (1 - a[i,i]) x[i] + =
... + -a[i,n] x[n] ) + b[i]
```
qed

Note that actually solving the equation system "A * x =3D = b"=20 means to calculate
```
        a[i,1] x[1] + ... + a[i,i] x[i] + =
... + a[i,n] x[n]  =3D  b[i]

  <=3D=3D>  a[i,i] x[i]  =3D
        b[i]
        - ( a[i,1] x[1] + ... + a[i,i] x[i] + ... + a[i,n] x[n] )
        + a[i,i] x[i]

  <=3D=3D>  x[i]  =3D
        ( b[i]
            - ( a[i,1] x[1] + ... + a[i,i] x[i] + ... + a[i,n] x[n] )
            + a[i,i] x[i]
        ) / a[i,i]

  <=3D=3D>  x[i]  =3D
        ( b[i] -
            ( a[i,1] x[1] + ... + a[i,i-1] x[i-1] +
              a[i,i+1] x[i+1] + ... + a[i,n] x[n] )
        ) / a[i,i]
```
There is one major restriction, though: a fix-point iteration is = guaranteed=20 to converge only if the first derivative of the iteration function has = an=20 absolute value less than one in an area around the point = "x(*)"=20 for which "Phi( x(*) ) =3D x(*)" is to be true, and if = the start=20 vector "x(0)" lies within that area!

This is best verified grafically, which unfortunately is impossible = to do=20 in this textual documentation!

See literature on Numerical Analysis for details!

In our case, this restriction translates to the following three=20 conditions:

There must exist a norm so that the norm of the matrix of the = iteration=20 function, ( En - A ), has a value less than one, the = matrix "A"=20 may not have any zero value on its main diagonal and the initial = vector=20 "x(0)" must be "good enough", i.e., "close enough" to the = solution "x(*)".

(Remember school math: the first derivative of a straight line = given by=20 "y =3D a * x + b" is "a"!)

The three methods expect a (quadratic!) matrix = "$matrix" as=20 their first argument, a start vector "$x0_vector", a = vector=20 "$b_vector" (which is the vector "b" in your equation = system=20 "A * x =3D b"), in the case of the "Relaxation Method" = ("RM"), a=20 real number "$weight" best between zero and two, and = finally an=20 error limit (real number) "$epsilon".

(Note that the weight "$weight" used by the = "Relaxation=20 Method" ("RM") is NOT checked to lie within any reasonable = range!)

The three methods first test the first two conditions of the three=20 conditions listed above and return an empty list if these conditions = are not=20 fulfilled.

Therefore, you should always test their return value using some = code=20 like:
```
  if ( $xn_vector =3D =
$A_matrix->solve_GSM($x0_vector,$b_vector,1E-12) )
  {
      # do something with the solution...
  }
  else
  {
      # do something with the fact that there is no solution...
  }
```
Otherwise, they iterate until abs( Phi(x) - x ) <=20 epsilon.

(Beware that theoretically, infinite loops might result if the = starting=20 vector is too far "off" the solution! In practice, this shouldn't be a = problem. Anyway, you can always press <ctrl-C> if you think that = the=20 iteration takes too long!)

The difference between the three methods is the following:

In the "Global Step Method" ("GSM"), the new vector = "x(t+1)"=20 (called "y" here) is calculated from the vector "x(t)" = (called=20 "x" here) according to the formula:
```
  y[i] =3D
  ( b[i]
      - ( a[i,1] x[1] + ... + a[i,i-1] x[i-1] +
          a[i,i+1] x[i+1] + ... + a[i,n] x[n] )
  ) / a[i,i]
```
In the "Single Step Method" ("SSM"), the components of the vector=20 "x(t+1)" which have already been calculated are used to = calculate=20 the remaining components, i.e.
```
  y[i] =3D
  ( b[i]
      - ( a[i,1] y[1] + ... + a[i,i-1] y[i-1] +  # note the "y[]"!
          a[i,i+1] x[i+1] + ... + a[i,n] x[n] )  # note the "x[]"!
  ) / a[i,i]
```
In the "Relaxation method" ("RM"), the components of the vector=20 "x(t+1)" are calculated by "mixing" old and new value = (like cold=20 and hot water), and the weight "$weight" determines the=20 "aperture" of both the "hot water tap" as well as of the "cold water = tap",=20 according to the formula:
```
  y[i] =3D
  ( b[i]
      - ( a[i,1] y[1] + ... + a[i,i-1] y[i-1] +  # note the "y[]"!
          a[i,i+1] x[i+1] + ... + a[i,n] x[n] )  # note the "x[]"!
  ) / a[i,i]
  y[i] =3D weight * y[i] + (1 - weight) * x[i]
```
Note that the weight "$weight" should be greater than = zero and=20 less than two (!).

The three methods are supposed to be of different efficiency.=20 Experiment!

Remember that in most cases, it is probably advantageous to first=20 "normalize()" your equation system prior to solving it!

OVERLOADED OPERATORS

SYNOPSIS

Unary operators:=20
"-", "~", "abs", = test,=20 "!", '""'
Binary (arithmetic) operators:=20
"+", "-", "*", = "**",=20 "+=3D", "-=3D", "*=3D",=20 "**=3D"
Binary (relational) operators:=20
"=3D=3D", "!=3D", "<",=20 "<=3D", ">", = ">=3D"

"eq", "ne", "lt", = "le",=20 "gt", "ge"

Note that the latter ("eq", "ne", ... ) = are just=20 synonyms of the former ("=3D=3D", "!=3D", = ... ), defined=20 for convenience only.

DESCRIPTION

'-'=20

Unary minus

Returns the negative of the given matrix, i.e., the matrix with all = elements multiplied with the factor "-1".

Example:

    $matrix =3D -$matrix;

'~'=20

Transposition

Returns the transposed of the given matrix.

Examples:

    $temp =3D ~$vector * $vector;
    $length =3D sqrt( $temp->element(1,1) );

    if (~$matrix =3D=3D $matrix) { # matrix is symmetric ... }

abs=20

Norm

Returns the "one"-Norm of the given matrix.

Example:

    $error =3D abs( $A * $x - $b );

test=20

Boolean test

Tests wether there is at least one non-zero element in the = matrix.

Example:

    if ($xn_vector) { # result of iteration is not =
zero ... }

'!'=20

Negated boolean test

Tests wether the matrix contains only zero's.

Examples:

    if (! $b_vector) { # heterogenous equation =
system ... }
    else             { # homogenous equation system ... }

    unless ($x_vector) { # not the null-vector! }

'""""'=20

"Stringify" operator

Converts the given matrix into a string.

Uses scientific representation to keep precision loss to a minimum = in case=20 you want to read this string back in again later with = "new_from_string()".

Uses a 13-digit mantissa and a 20-character field for each element = so that=20 lines will wrap nicely on an 80-column screen.

Examples:

    $matrix =3D =
Math::MatrixReal->new_from_string(<<"MATRIX");
    [ 1  0 ]
    [ 0 -1 ]
    MATRIX
    print "$matrix";

    [  1.000000000000E+00  0.000000000000E+00 ]
    [  0.000000000000E+00 -1.000000000000E+00 ]

    $string =3D "$matrix";
    $test =3D Math::MatrixReal->new_from_string($string);
    if ($test =3D=3D $matrix) { print ":-)\n"; } else { print ":-(\n"; =
}

'+'=20

Addition

Returns the sum of the two given matrices.

Examples:

    $matrix_S =3D $matrix_A + $matrix_B;

    $matrix_A +=3D $matrix_B;

'-'=20

Subtraction

Returns the difference of the two given matrices.

Examples:

    $matrix_D =3D $matrix_A - $matrix_B;

    $matrix_A -=3D $matrix_B;

Note that this is the same as:

    $matrix_S =3D $matrix_A =
+ -$matrix_B;

    $matrix_A +=3D -$matrix_B;

(The latter are less efficient, though)

'*'=20

Multiplication

Returns the matrix product of the two given matrices or the product = of the=20 given matrix and scalar factor.

Examples:

    $matrix_P =3D $matrix_A * $matrix_B;

    $matrix_A *=3D $matrix_B;

    $vector_b =3D $matrix_A * $vector_x;

    $matrix_B =3D -1 * $matrix_A;

    $matrix_B =3D $matrix_A * -1;

    $matrix_A *=3D -1;

'**'=20

Exponentiation

Returns the matrix raised to an integer power. If 0 is passed, the = identity=20 matrix is returned. If a negative integer is passed, it computes the = inverse=20 (if it exists) and then raised the inverse to the absolute value of = the=20 integer. The matrix must be quadratic.

Examples:

    $matrix2 =3D $matrix ** 2;

    $matrix **=3D 2;

    $inv2 =3D $matrix ** -2;

    $ident =3D $matrix ** 0;

'=3D=3D'=20

Equality

Tests two matrices for equality.

Example:

    if ( $A * $x =3D=3D $b ) { print "EUREKA!\n"; =
}

Note that in most cases, due to numerical errors (due to the finite = precision of computer arithmetics), it is a bad idea to compare two = matrices=20 or vectors this way.

Better use the norm of the difference of the two matrices you want = to=20 compare and compare that norm with a small number, like this:

 =
   if ( abs( $A * $x - $b ) < 1E-12 ) { print "BINGO!\n"; }

'!=3D'=20

Inequality

Tests two matrices for inequality.

Example:

    while ($x0_vector !=3D $xn_vector) { # proceed =
with iteration ... }

(Stops when the iteration becomes stationary)

Note that (just like with the '=3D=3D' operator), it is usually a = bad idea to=20 compare matrices or vectors this way. Compare the norm of the = difference of=20 the two matrices with a small number instead.

'<'=20

Less than

Examples:

    if ( $matrix1 < $matrix2 ) { # ... }

    if ( $vector < $epsilon ) { # ... }

    if ( 1E-12 < $vector ) { # ... }

    if ( $A * $x - $b < 1E-12 ) { # ... }

These are just shortcuts for saying:

    if ( abs($matrix1) =
< abs($matrix2) ) { # ... }

    if ( abs($vector) < abs($epsilon) ) { # ... }

    if ( abs(1E-12) < abs($vector) ) { # ... }

    if ( abs( $A * $x - $b ) < abs(1E-12) ) { # ... }

Uses the "one"-norm for matrices and Perl's built-in "abs()" for=20 scalars.

'<=3D'=20

Less than or equal

As with the '<' operator, this is just a shortcut for the same=20 expression with "abs()" around all arguments.

Example:

    if ( $A * $x - $b <=3D 1E-12 ) { # ... =
}

which in fact is the same as:

    if ( abs( $A * $x - $b ) =
<=3D abs(1E-12) ) { # ... }

Uses the "one"-norm for matrices and Perl's built-in "abs()" for=20 scalars.

'">'>'=20

Greater than

As with the '<' and '<=3D' operator, this

    if ( =
$xn - $x0 > 1E-12 ) { # ... }

is just a shortcut for:

    if ( abs( $xn - $x0 ) > =
abs(1E-12) ) { # ... }

Uses the "one"-norm for matrices and Perl's built-in "abs()" for=20 scalars.

=3D'">'>=3D'=20

Greater than or equal

As with the '<', '<=3D' and '>' operator, the = following

    if ( $LR >=3D $A ) { # ... }

is simply a shortcut for:

    if ( abs($LR) >=3D abs($A) =
) { # ... }

Uses the "one"-norm for matrices and Perl's built-in "abs()" for=20 scalars.

VERSION

This man page documents Math::MatrixReal version 1.9.

The latest version can always be found at http://leto.net/code/Math-= MatrixReal/

AUTHORS

Steffen Beyer <sb@engelschall.com>, Rodolphe Ortalo=20 <ortalo@laas.fr>, Jonathan Leto <jonathan@leto.net>.

Currently maintained by Jonathan Leto, send all bugs/patches to = me.

CREDITS

Many thanks to Prof. Pahlings for stoking the fire of my enthusiasm = for=20 Algebra and Linear Algebra at the university (RWTH Aachen, Germany), and = to=20 Prof. Esser and his assistant, Mr. Jarausch, for their fascinating = lectures in=20 Numerical Analysis!

COPYRIGHT

LICENSE AGREEMENT

This package is free software; you can redistribute it and/or modify = it under=20 the same terms as Perl itself.

NAME

DESCRIPTION

SYNOPSIS

Constructor Methods And = Such

Matrix Row, = Column and=20 Element operations

Matrix Operations

Arithmetic Operations

Boolean Matrix = Operations

Eigensystems

Miscellaneous

OVERLOADED OPERATORS

SYNOPSIS

DESCRIPTION

SEE ALSO

VERSION

AUTHORS

CREDITS

COPYRIGHT

LICENSE AGREEMENT