 | |
 | |
 | |
 | ||
 |
 |
 |
 | ||
 | ||||
 |
|
Mathematics Packages > Sheafhom | |
|
| ||||
      | ||||
| ||||
|
|
most recent version : 2.1 |release date : March 17, 2005 | |
|
| ||||
|
|
Definition of Smith normal form. The Smith normal form (SNF) of an integer matrix is a way of breaking it down into basic pieces while respecting the integer structure. Let's begin with an example. International Allegro CL Professional Edition 7.0 [Windows] (Oct 20, 2004 22:05) Copyright (C) 1985-2004, Franz Inc., Oakland, CA, USA. All Rights Reserved. CL-USER(2): :cd \mark\math\shh2 c:\mark\math\shh2\ CL-USER(3): (load "make") [list of files being compiled and loaded] Enter a 3 by 3 matrix, listing the entries row by row. CL-USER(4): (defvar mm (input-csparse 3 3 '(5 -1 0 0 2 10 7 5 32))) MM CL-USER(5): mm 5 -1 0 0 2 10 7 5 32 This matrix is an instance of the class csparse. The name stands for "column sparse". A csparse stores a matrix by columns, and each column is a sparse vector. Now compute the SNF of the matrix and take a look at its internals.
CL-USER(6): (defvar snf-mm (make-snf mm))
SNF-MM
CL-USER(7): (describe snf-mm)
[SNF 3 by 3, [one unit] 2 [one 0]] is an instance of
#<STANDARD-CLASS SNF>:
The following slots have :INSTANCE allocation:
MAT
5 -1 0
0 2 10
7 5 32
D
-1 0 0
0 2 0
0 0 0
USE-P T
USE-Q T
P-LIST ((T 1 0 -2) (T 2 0 -5) (T 2 1 3) (SHEAFHOM2::P 1 2)
(T 2 1 5))
Q-LIST ((T 1 2 1) (T 0 1 -5) (SHEAFHOM2::P 0 1))
This say mm is a product of matrices
The P's and Q's are elementary matrices. For instance, P0, the 0th element of p-list, is a transvection matrix, built by taking the identity matrix and inserting -2 in the 1,0 position: 1 0 0 -2 1 0 0 0 1 For a given matrix A, the product P0 A alters the middle row of A by subtracting off two times the top row. P3 is the 2 by 2 permutation matrix 1 0 0 0 0 1 0 1 0 The product P3 A interchanges the middle and bottom rows of A. Multiplying by elementary matrices on the right performs the analogous column operations. The P's and Q's are isomorphisms representing changes of basis. The essential structure of mm is in the central matrix D, the diagonal matrix -1 0 0 0 2 0 0 0 0 Because there are only two non-zero diagonal entries, mm has rank 2. Viewing mm for the moment over the real numbers, rank 2 means the image of mm has only two dimensions instead of three. The kernel of mm is one-dimensional, generated by the following vector.
CL-USER(8): (shh2::ker snf-mm)
-1
-5
1
(We need the prefix "shh2::" because the current implementations of kernel, cokernel, and so on are provisional. We expect them to be made public in the next release.) Over the real numbers, the quotient of three-dimensional space by the image of mm is one-dimensional. Over the integers, there is more to the story, because the middle entry of D is 2 rather than +/- 1. The quotient is Z + (Z / 2 Z) including the torsion subgroup Z / 2 Z. Printing an SNF summarizes what's on the diagonal of D. CL-USER(9): snf-mm [SNF 3 by 3, [one unit] 2 [one 0]] A "unit" is an invertible element in a ring. Over the integers, the units are +/- 1. After the units, we list the torsion coefficients and the number of trailing zeroes. Operating on csparses. We'll need the next function later in the tutorial.
CL-USER(10): (defun csparse-incf (matrix i j)
"Adds 1 to the entry in the i,j position."
(let* ((entry (csparse-get matrix i j))
(val (if (null entry) 0 (sp-v entry)))
(new-val (1+ val))
(new-entry (make-sp i new-val)))
(csparse-set matrix new-entry j)))
CSPARSE-INCF
Let's walk through this function line by line, applying its steps to the csparse mm. CL-USER(11): mm 5 -1 0 0 2 10 7 5 32 We'll add a 1 to the 1,2 position, changing the 10 to 11. A csparse is stored by columns, and each column is a sparse vector. The "10" in column 2 must maintain its index in the vector (1) in addition to its value (10). Sheafhom wraps this data into an instance of the type sparse-elt. The function csparse-get extracts a sparse-elt. CL-USER(12): (csparse-get mm 1 2) (1 . 10) Over the integers, a sparse-elt is just a cons holding the index and value. For other rings, the representation could be more complicated. We extract the index of a sparse-elt with the macro sp-i, and the value with sp-v. CL-USER(13): (sp-v *) 10 To replace the 10 with 11, we need to create a new sparse-elt with index 1 and value 11. The macro make-sp creates a sparse-elt. CL-USER(14): (make-sp 1 11) (1 . 11) The function csparse-set puts the new sparse-elt into column 2 of the csparse. CL-USER(15): (csparse-set mm * 2) 5 -1 0 0 2 11 7 5 32 This change raises the rank of mm to 3. There's also a lot more torsion: the quotient by the image is cyclic of order 32. CL-USER(16): (make-snf mm) [SNF 3 by 3, [two units] 32 [zero 0s]] Sheafhom has many more operators. Check the documentation or the source files for more information. Operators on sparse-elts have "sp" in their names. Operators on sparse vectors (type sparse-v) have "v", csparses have "csparse", and general matrix operations have "m". An operator named "make-" is a basic constructor, and "input-" is a user-friendly constructor for use at the command prompt. "Alter" and "swap" provide elementary row/column operations. Modulo the dictionary. The funny article Homophonic quotients of free groups looks at the free group on the letters A through Z modulo equivalence of words that have the same pronunciation in French. For instance, soie = soi implies e = 1. The authors prove the quotient is trivial. Let's solve a related problem: find the quotient of the free abelian group on A through Z by the relation WWWW = 0, where WWWW runs through a complete English dictionary. If the quotient is too boring, we can also try subsets of the dictionary. We'll use this dictionary, which comes from a site that's of interest to Lispers for several reasons. :-) In the dictionary, no hyphens or other punctuation marks appear, just the letters a through z. Words of length one are excluded because they make combinatorial problems like this one too easy. We load the words into a list for flexibility.
CL-USER(17): (defconstant +words+
(let ((result '()))
(with-open-file (stream "WORD.LST" :direction :input)
(loop
(let ((symbol (read stream nil 'end-of-file)))
(if (eq symbol 'end-of-file)
(return (nreverse result))
(push (symbol-name symbol) result)))))))
+WORDS+
CL-USER(18): +words+
("AA" "AAH" "AAHED" "AAHING" "AAHS" "AAL" "AALII" "AALIIS" "AALS"
"AARDVARK" ...)
CL-USER(19): (length +words+)
173528
CL-USER(20): (nth 100000 +words+)
"NONJURY"
What's the range of word lengths? CL-USER(21): (remove-duplicates (sort (mapcar #'length +words+) #'<)) (2 3 4 5 6 7 8 9 10 11 ...) While we're at it, what is the longest word?
CL-USER(22): (last *)
(28)
CL-USER(23): (remove-if-not #'(lambda (w) (= (length w) 28)) +words+)
("ETHYLENEDIAMINETETRAACETATES")
ANTIDISESTABLISHMENTARIANISM is also 28 letters long, though it isn't in this dictionary for some reason. Back to the problem on free abelian groups. We'll create a 26 by 173528 csparse with one column for each word. The column for DAD is 1 0 0 2 0 ... 0 standing for one A and two D's. Here is a utility for converting A through Z to the row indices 0 through 25.
CL-USER(24): (defconstant +code-A+ (char-code #\A))
+CODE-A+
CL-USER(25): (defun az-to-int (ch)
"Converts ASCII A to 0, B to 1, ..., Z to 25."
(let ((result (- (char-code ch) +code-A+)))
(assert (<= 0 result 25) (ch)
"Not between A and Z: ~:C" ch)
result))
AZ-TO-INT
The next function takes a list of words and returns the csparse.
CL-USER(26): (defun az-list-to-csparse (words)
(let ((result (make-csparse-zero 26 (length words)))
(j 0))
(dolist (word words result)
(map nil #'(lambda (ch)
(csparse-incf result (az-to-int ch) j))
word)
(incf j))))
AZ-LIST-TO-CSPARSE
Let's apply this to the whole dictionary. CL-USER(27): (compile 'az-to-int) ;; to speed things up AZ-TO-INT NIL NIL CL-USER(28): (compile 'az-list-to-csparse) AZ-LIST-TO-CSPARSE NIL NIL CL-USER(29): (defvar words-all (az-list-to-csparse +words+)) WORDS-ALL ;; after a little while CL-USER(30): words-all [CSPARSE 26 by 173528] Finally, we find the SNF.
CL-USER(31): (let ((*show-stats* t))
(make-snf words-all nil nil t))
[SNF 26 by 173528, [twenty-six units] [zero 0s]]
Since the SNF has nothing but units down the diagonal, we've proved a theorem: the free abelian group on A through Z modulo the dictonary is trivial! We need to explain several items in the last form. When *show-stats* is true, make-snf pops up a line graph showing how the pattern of sparsity changed during the computation. The graph shows that the matrix started with 27.8% of the entries non-zero. The sparsity decreased slightly while the first 20 of the 26 rows were being reduced, then increased till the end. When *show-prog* is t, make-snf prints timing data to standard output. Another flag *show-csw* will be explained below. The last argument t to make-snf means the SNF was computed destructively. This is crucial for saving memory with large matrices. The entries of words-all were modified in place, with sparse-elts thrown to the garbage collector whenever possible. If the last argument had been nil (the default), words-all would have been preserved, and the snf would have been computed using a copy. As another way to save space, we gave the arguments "nil nil" to make-snf. These stop the program from storing the P's and Q's. Only the main D portion was recorded. As we've suggested, the trivial group is boring. What is the quotient of A through Z modulo the two-letter words? The quotient will be at least Z / 2 Z for parity reasons: no word with an odd number of letters can be equivalent to 0. But are there other relations?
CL-USER(32): (remove-if-not #'(lambda (w) (= (length w) 2)) +words+)
("AA" "AB" "AD" "AE" "AG" "AH" "AI" "AL" "AM" "AN" ...)
CL-USER(33): (length *)
96
CL-USER(34): (az-list-to-csparse **)
[CSPARSE 26 by 96]
With *show-csw* true, we see a movie of how the sparsity pattern changes as the SNF is computed. Csw stands for "csparse window".
CL-USER(35): (let ((*show-csw* t))
(make-snf * nil nil t))
[SNF 26 by 96, [twenty-one units] -2 [four 0s]]
The torsion coefficient 2 (up to sign) is the one we predicted. But the matrix has rank 22. Roughly speaking, there are four linearly independent "words" that are not combinations of two-letter words. The quotient is Z^4 + (Z / 2 Z). What about four-letter words?
CL-USER(36): (remove-if-not #'(lambda (w) (= (length w) 4)) +words+)
("AAHS" "AALS" "ABAS" "ABBA" "ABBE" "ABED" "ABET" "ABLE" "ABLY" "ABOS"
...)
CL-USER(37): (az-list-to-csparse **)
[CSPARSE 26 by 3919]
CL-USER(38): (make-snf * nil nil t)
[SNF 26 by 3919, [twenty-five units] 4 [zero 0s]]
We have torsion of order 4, again for parity reasons. Apart from that, the quotient by four-letter words is a small as possible. What about lengths greater than or equal to some cut-off value? The same methods show that the quotient by words of length 20 or more is trivial. But for length 21 or more, the quotient is Z.
CL-USER(39): (remove-if-not #'(lambda (w) (>= (length w) 21)) +words+)
("ACETYLCHOLINESTERASES" "ADRENOCORTICOSTEROIDS"
"ADRENOCORTICOTROPHINS" "ANTHROPOMORPHIZATIONS"
"ANTIAUTHORITARIANISMS" "ANTIFERROMAGNETICALLY"
"BUCKMINSTERFULLERENES" "CARBOXYMETHYLCELLULOSE"
"CARBOXYMETHYLCELLULOSES" "CLINICOPATHOLOGICALLY" ...)
CL-USER(40): (make-snf (az-list-to-csparse *) t nil t)
[SNF 26 by 118, [twenty-five units] [one 0]]
Note the arguments "t nil t", meaning we remember the P's but not the Q's, and work destructively. Because we remembered the P's, we can use shh::coker-section, which picks out a generator for the quotient. CL-USER(41): (shh2::coker-section *) The result is a column vector with all entries 0 except for a 1 standing for the letter Q. We've proved that CARBOXYMETHYLQCELLULOSE is not a word. Topology in dimension two. Sheafhom was written to solve problems in algebraic topology and derived categories. The file two-complex.lisp sets up a simplicial complex X. It's restricted to dimension two to make things easier. The method make-d turns the boundary maps of X into csparses. The method print-homology produces the homology of X with coefficients in Z, including the torsion. Here is a brief example. Consider a 3 by 3 grid of squares with vertices numbered as shown. 0 2 1 0 6 7 8 6 3 4 5 3 0 1 2 0 We list the nine squares. CL-USER(42): *klein-bottle-cells* ((6 7 2 0) (7 8 1 2) (8 6 0 1) (3 4 7 6) (4 5 8 7) (5 3 6 8) (0 1 4 3) (1 2 5 4) (2 0 3 5)) Now attach the left side of the grid to the right side by gluing 0 to 0, 3 to 3, and 6 to 6. The result is a cylindrical tube. Next, glue the bottom circle of the cylinder to the top circle, 0 to 0, 1 to 1, and 2 to 2. The trouble is that the order is reversed, 0-1-2-0 on the bottom and 0-2-1-0 on the top. After squinting, you'll see you can't do the gluing in three dimensions. You have to bend the tube into the fourth dimension and turn one end "upside down" to match up the circles in the proper order. The resulting space is a manifold of dimension two, the Klein bottle. The file two-complex.lisp lets you construct the Klein bottle this way. CL-USER(43): (defvar *klein-bottle* (make-two-complex *)) *KLEIN-BOTTLE* CL-USER(44): (print-homology *klein-bottle*) H_0 is free of rank 1 H_1 has free part of rank 1 and torsion coeffs (2) H_2 is free of rank 0 NIL | |
|
| ||||
 | ||||
 | ||||
|
|
Mark McConnell | |
|
| ||||
|
|
 |
|
|
| ||||
 | ||||