Integer Smith Normal Form Package

Victor Alvarez Solano, Jose Andres Armario Sampalo, Maria Dolores Frau Garcia and Pedro Real Jurado
Departamento de Matematica Aplicada I, University of Seville (Spain)
6 November 2006

valvarez@us.es

In spite of the fact that some packages (even in Mathematica) for computing the integer Smith Normal Form of a matrix over the integers are already available, some of them eventually fail to give the diagonal matrix properly sorted. This notebook provides an alternate method for computing the Smith Normal Form for matrix over the integers. Notice that the point of this notebook is correctness better than efficiency, so that previous algorithms will eventually run faster than ours.

Reference

Title

Integer Smith Normal Form

Author

Victor Alvarez Solano, Jose Andres Armario Sampalo, Maria Dolores Frau Garcia and Pedro Real Jurado

Summary

This package is designed to find out the Integer Smith Normal Form (the diagonal matrix, and the change basis matrices) of a matrix over the integers.

Copyright

Rights reserved by Victor Alvarez Solano et al, Departamento de Matematica Aplicada I, University of Seville (Spain).

Notebook Version

2.0

Mathematica Version

4.0

History

A Mathematica notebook for computing the Smith Normal Form is already known (see Jabon, 1994, for instance). Unfortunately, this package fails when working with some matrices (see example 1 below). This source of difficulties has provided inspiration for this work. The present package is used in other notebooks provided by the authors, http://library.wolfram.com/infocenter/MathSource/6384/ and http://library.wolfram.com/infocenter/MathSource/6516/

Keywords

Smith Normal Form.

Limitation

This package works only on matrices over the integers (it is not designed to work on polynomials).

Discussion

References: 
D. Jabon and A. Pascoletti. Smith Normal Forms. http://library.wolfram.com/infocenter/MathSource/682/

Requirements

Context`IntegerSmithNormalForm`

Interface

This part declares the publicly visible functions, options, and values.

Set up the package context, including public imports

Off[General::spell];
Off[General::spell1];
BeginPackage["`IntegerSmithNormalForm`"]

$ContextPath

Usage messages for the exported functions and the context itself

IntegerSmithNormalForm::usage="IntegerSmithNormalForm[A] gives, for a matrix over the integers A, a list {D,{P,Q}}, D being the integer Smith Normal Form related to A, and P,Q being the change basis matrices, so that D=P·A·Q. Notice that the matrices P and Q are not uniquely determined."

Implementation

This part contains the actual definitions and any auxiliary functions that should not be visible outside.

Begin the private context (implementation part)

Begin["`Private`"]
{$Context,$ContextPath}

Definition of auxiliary functions and local (static) variables

Most of auxiliary functions uses the maps `chi` to be defined by the user, so we have decided to include these maps `chi` as an argument of the main command. This way, the auxiliary functions are incorporated to the body of the principal function, so that there is no need for adding the maps `chi` as arguments anymore.

Main program

M[i_, j_, t_, k_] := 
        ReplacePart[IdentityMatrix[k], 
          t, {i, j}];
T[i_, j_, k_] := 
        Module[{m, n}, m = IdentityMatrix[k]; n = m[[i]]; 
          m[[i]] = m[[j]]; m[[j]] = n; 
          m];
Completar[C_, m_, n_] := 
        Module[{k, i, j, ii}, i = m - Length[C]; k = {}; 
          Do[k = Append[k, Table[0, {ii, 1, n}]], {j, 1, i}]; 
          Do[k = Append[k, 
                Join[Table[0, {ii, 1, n - Length[C[[1]]]}], 
                  C[[j]]]], {j, 1, Length[C]}]; 
          k];
IntegerSmithNormalForm[A_] := 
        Module[{k, C, coli, colj, P, Q, m, n, H, i, it, control}, 
          m = Length[A]; n = Length[A[[1]]]; P = {IdentityMatrix[m]}; 
          Q = {IdentityMatrix[n]}; C = A; i = 1; 
          Fin = Max[Abs[Take[C, {i, m}, {i, n}]]] > 
              0; 
While[
            Fin, {coli, 
                colj} = Position[
                  Completar[Take[Abs[C], {i, m}, {i, n}], m, n], 
                  Min[Select[
                      Flatten[
                        Take[Abs[C], {i, m}, {i, n}]], #1 > 
                          0 &]]][[1]]; P = Prepend[P, T[i, coli, m]]; 
            Q = Append[Q, T[i, colj, n]]; it = i + 1; 
            C = First[P] . C . Last[Q]; 
            control = True; 
While[it ≤ n, 
              If[IntegerQ[C[[i, it]]/C[[i, i]]], it = it + 1, 
                Q = Append[Q, 
                    M[i, it, -Quotient[C[[i, it]], C[[i, i]]], 
                        n] . T[it, i, n]]; C = C . Last[Q]; control = False; 
                it = n + 1]
]; 
            If[control, 
              Do[Q = Append[Q, 
                    M[i, it, -C[[i, it]]/C[[i, i]], n]]; 
                C = C . Last[Q], {it, i + 1, n}]; 
              it = i + 1; 
While[it ≤ m, 
                If[IntegerQ[C[[it, i]]/C[[i, i]]], it = it + 1, 
                  P = Prepend[P, 
                      T[it, i, m] . 
                        M[it, i, -Quotient[C[[it, i]], 
                              C[[i, i]]], m]]; C = First[P] . C; 
                  control = False; it = m + 1]
]; 
              If[control, 
                Do[P = Prepend[P, 
                      M[it, i, -C[[it, i]]/C[[i, i]], m]]; 
                  C = First[P] . C, {it, i + 1, m}]; 
Catch[
                  Do[
Do[
                      If[IntegerQ[C[[coli, colj]]/C[[i, i]]], , 
                        P = Prepend[P, M[i, coli, 1, m]]; C = First[P] . C; 
                        control = False; Throw[False]], {colj, i + 1, 
                        n}], {coli, i + 1, m}]; Throw[True]]; 
                If[ control, 
                  P = Prepend[P, M[i, i, Sign[C[[i, i]]], m]]; 
                  C = First[P] . C; i = i + 1; 
                  Fin = Max[Abs[Take[C, {i, m}, {i, n}]]] > 
                      0]
];
];
]; 
P = Fold[Dot, IdentityMatrix[m], P]; 
          Q = Fold[Dot, IdentityMatrix[n], Q]; {C, {P, 
              Q}}
];

End the private context

End[ ]

Epilog

This section ends the package.

End the package context

EndPackage[ ]

Examples, Tests

Examples, tests for the use of the package can go here.

Example 1

Let consider the matrix (({{11, 0, 0}, {0, 0, 0}, {-9, 0, 0}, {4, 2, 0}, {-4, 0, 0}, {0, 0, 2}})).

A=(({{11, 0, 0}, {0, 0, 0}, {-9, 0, 0}, {4, 2, 0}, {-4, 0, 0}, {0, 0, 2}}));{F, {P, Q}} = IntegerSmithNormalForm[A];Print[F // MatrixForm];Print[P // MatrixForm];Print[Q // MatrixForm];

The diagonal matrix F is obtained as the product P·A·Q

P . A . Q == F