[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Axiom-developer] 20070927.01.tpd.patch
From: |
daly |
Subject: |
[Axiom-developer] 20070927.01.tpd.patch |
Date: |
Thu, 27 Sep 2007 22:16:22 -0500 |
This patch adds Martin Rubey's implementation of Gunter Rote's Pfaffian
algorithm as well as some documentation and regression test cases.
========================================================================
diff --git a/changelog b/changelog
index d085c48..c254ce6 100644
--- a/changelog
+++ b/changelog
@@ -1,3 +1,5 @@
+20070927 tpd src/input/Makefile pfaffian regression test added
+20070927 mxr src/input/pfaffian.input regression test added
20070920 tpd src/input/Makefile add bug101.input regression test
20070920 tpd src/input/bug101.input test laplace(log(z),z,w) bug 101
20070920 wxh src/algebra/laplace.spad fix laplace(log(z),z,w) bug 101
diff --git a/src/input/Makefile.pamphlet b/src/input/Makefile.pamphlet
index b2d5fe9..b9e5374 100644
--- a/src/input/Makefile.pamphlet
+++ b/src/input/Makefile.pamphlet
@@ -338,9 +338,9 @@ REGRES= algaggr.regress algbrbf.regress algfacob.regress
alist.regress \
oct.regress ode.regress odpol.regress op1.regress \
opalg.regress operator.regress op.regress ovar.regress \
padic.regress parabola.regress pascal1.regress pascal.regress \
- patch51.regress \
+ patch51.regress page.regress \
patmatch.regress pat.regress perman.regress perm.regress \
- pfr1.regress pfr.regress page.regress pmint.regress \
+ pfaffian.regress pfr1.regress pfr.regress pmint.regress \
poly1.regress polycoer.regress poly.regress psgenfcn.regress \
quat1.regress quat.regress r20abugs.regress r20bugs.regress \
r21bugsbig.regress r21bugs.regress radff.regress radix.regress \
@@ -592,8 +592,9 @@ FILES= ${OUT}/algaggr.input ${OUT}/algbrbf.input
${OUT}/algfacob.input \
${OUT}/pascal1.input \
${OUT}/pascal.input \
${OUT}/patch51.input \
- ${OUT}/patmatch.input ${OUT}/perman.input \
- ${OUT}/perm.input ${OUT}/pfr.input ${OUT}/pfr1.input \
+ ${OUT}/patmatch.input ${OUT}/perman.input \
+ ${OUT}/perm.input ${OUT}/pfaffian.input \
+ ${OUT}/pfr.input ${OUT}/pfr1.input \
${OUT}/pinch.input ${OUT}/plotfile.input ${OUT}/pollevel.input \
${OUT}/pmint.input ${OUT}/polycoer.input \
${OUT}/poly1.input ${OUT}/psgenfcn.input \
@@ -869,7 +870,8 @@ DOCFILES= \
${DOC}/pat.input.dvi ${DOC}/patch51.input.dvi \
${DOC}/patmatch.input.dvi \
${DOC}/pdecomp0.as.dvi ${DOC}/perman.input.dvi \
- ${DOC}/perm.input.dvi ${DOC}/pfr1.input.dvi \
+ ${DOC}/perm.input.dvi ${DOC}/pfaffian.input.dvi \
+ ${DOC}/pfr1.input.dvi \
${DOC}/pfr.input.dvi ${DOC}/pinch.input.dvi \
${DOC}/plotfile.input.dvi ${DOC}/plotlist.input.dvi \
${DOC}/pollevel.input.dvi ${DOC}/poly1.input.dvi \
diff --git a/src/input/pfaffian.input.pamphlet
b/src/input/pfaffian.input.pamphlet
new file mode 100755
index 0000000..ad4d7f0
--- /dev/null
+++ b/src/input/pfaffian.input.pamphlet
@@ -0,0 +1,483 @@
+\documentclass{article}
+\usepackage{amsfonts}
+\usepackage{axiom}
+\begin{document}
+\title{\$SPAD/src/input pfaffian}
+\author{Timothy Daly, Gunter Rote, Martin Rubey}
+\maketitle
+\begin{abstract}
+In mathematics, the determinant of a skew-symmetric matrix can always
+be written as the square of a polynomial in the matrix entries. This
+polynomial is called the Pfaffian of the matrix. The Pfaffian is
+nonvanishing only for $2n \times 2n$ skew-symmetric matrices, in which
+case it is a polynomial of degree $n$.
+\end{abstract}
+\eject
+\tableofcontents
+\eject
+\section{Examples}
+$$
+{\rm Pfaffian}\left[
+\begin{array}{cc}
+0 & a\\
+-a & 0
+\end{array}
+\right] = a
+$$
+
+$$
+{\rm Pfaffian}\left[
+\begin{array}{cccc}
+0 & a & b & c\\
+-a & 0 & d & e\\
+-b & -d & 0 & f\\
+-c & -e & -f & 0
+\end{array}
+\right] = af -be + dc
+$$
+
+$$
+{\rm Pfaffian}\left[
+\begin{array}{cccc}
+\begin{array}{cc}
+0 & \lambda_1\\
+-\lambda_1 & 0
+\end{array} & 0 & \cdots & 0\\
+0 &
+\begin{array}{cc}
+0 & \lambda_2\\
+-\lambda_2 & 0
+\end{array} & & 0\\
+\vdots & & \ddots & \vdots\\
+0 & 0 & \cdots &
+\begin{array}{cc}
+0 & \lambda_n\\
+-\lambda_n & 0
+\end{array}
+\end{array}
+\right] = \lambda_1\lambda_2\cdots\lambda_n
+$$
+\section{Formal definition}
+
+Let $A = \{a_{i,j}\}$ be a $2n \times 2n$ skew-symmetric matrix.
+The Pfaffian of A is defined by the equation
+
+$$Pf(A) = \frac{1}{s^n n!}\sum_{\sigma\in{}S_{2n}}sgn(\sigma)
+\prod_{i=1}^n{a_{\sigma(2i-1),\sigma(2i)}}$$
+
+where $S_{2n}$ is the symmetric group and $sgn(\sigma)$
+is the signature of $\sigma$.
+
+One can make use of the skew-symmetry of $A$ to avoid summing over all
+possible permutations. Let $\Pi$ be the set of all partitions of
+$\{1, 2, \ldots, 2n\}$ into pairs without regard to order.
+There are $(2n - 1)!!$ such partitions. An element
+$\alpha \in \Pi$, can be written as
+
+$$\alpha=\{(i_1,j_1),(i_2,j_2),\cdots,(i_n,j_n)\}$$
+
+with $i_k < j_k$ and $i_1 < i_2 < \cdots < i_n$. Let
+
+$$
+\pi=
+\left[
+\begin{array}{cccccc}
+1 & 2 & 3 & 4 & \cdots & 2n \\
+i_1 & j_1 & i_2 & j_2 & \cdots & j_{n}
+\end{array}
+\right]
+$$
+
+be a corresponding permutation. This depends only on the partition $\alpha$
+and not on the particular choice of $\Pi$. Given a partition $\alpha$ as above
+define
+
+$$A_\alpha =sgn(\pi)a_{i_1,j_1}a_{i_2,j_2}\cdots a_{i_n,j_n}$$
+
+The Pfaffian of $A$ is then given by
+
+$$Pf(A)=\sum_{\alpha\in\Pi} A_\alpha$$
+
+The Pfaffian of a $n\times n$ skew-symmetric matrix for n odd is
+defined to be zero.
+
+\subsection{Alternative definition}
+
+One can associate to any skew-symmetric $2n \times 2n$
+matrix $A=\{a_{ij}\}$ a bivector
+
+$$\omega=\sum_{i<j} a_{ij}~e^i\wedge e^j$$
+
+where $\{e^1, e^2, \ldots, e^{2n}\}$ is the standard basis of
+$\mathbb{R}^{2n}$. The Pfaffian is then defined by the equation
+
+$$\frac{1}{n!}\omega^n = \mbox{Pf}(A)~e^1\wedge e^2\wedge\cdots\wedge e^{2n}$$
+
+here $\omega^n$ denotes the wedge product of $n$ copies of
+$\omega^n$ with itself.
+
+\subsection{Derivation from Determinant}
+
+The Pfaffian can be derived from the determinant for a skew-symmetric
+matrix $A$ as follows. Using Laplace's formula we can write the
+determinant as
+
+$$\det(A)=(-1)^{p+1}a_{p1}A_{p1} +
+(-1)^{p+2}a_{p2}A_{p2}+\cdots+(-1)^{n+p}a_{pn}A_{pn}$$
+
+where $A_{pi}$ is the $p,i$ minor matrix of the matrix $A$. We further
+use Laplace's formula to note that
+
+$$\det(A[A_{ij}]) = \vert A \vert^n$$
+
+since this determinant is that of an $n \times n$ matrix whose only
+non-zero elements are the diagonals (each with value det(A)) and
+$[A_{ij}]$ is a matrix whose $i,j$th component is the corresponding $i,j$
+minor matrix. In this way, following a proof by Parameswaran, we can
+write the determinant as,
+
+$$\det(A)\equiv\Delta_n=
+\left|
+\begin{array}{cccc}
+a_{11}&a_{12}&\cdots&a_{1n}\\
+a_{21}&a_{22}&\cdots&a_{2n}\\
+\vdots&&&\vdots\\
+a_{n1}&a_{n2}&\cdots&a_{nn}
+\end{array}
+\right|
+$$
+
+The minor of
+$$
+\left|
+\begin{array}{cc}
+a_{11}&a_{12}\\
+a_{21}&a_{22}
+\end{array}
+\right|
+$$
+would be $\Delta_{n - 2}$. With this notation
+
+$$\left|
+\begin{array}{cccc}
+1&0&\cdots&0\\
+0&1&\cdots&0\\
+a_{31}&a_{32}&\cdots&a_{3n}\\
+\cdots&\cdots&\cdots&\cdots\\
+a_{n1}&a_{n2}&\cdots&a_{nn}
+\end{array}
+\right|
+\times
+\left|
+\begin{array}{cccc}
+A_{11}&A_{21}&\cdots&A_{n1}\\
+A_{12}&A_{22}&\cdots&A_{n2}\\
+A_{13}&A_{23}&\cdots&A_{n3}\\
+\cdots&\cdots&\cdots&\cdots\\
+A_{1n}&A_{2n}&\cdots&A_{nn}
+\end{array}
+\right| =
+\left|
+\begin{array}{cc}
+A_{11}&A_{21}\\
+A_{12}&a_{22}
+\end{array}
+\right|
+\Delta_n^{n-2}
+$$
+
+Thus
+
+$$\Delta_{n-2}\Delta_n^{n-1}=
+\left|
+\begin{array}{cc}
+A_{11}&A_{21}\\
+A_{12}&a_{22}
+\end{array}
+\right|
+\Delta_n^{n-2}
+$$
+
+Of course, it was only arbitrarily that we chose to remove the first
+two rows, and more generically we can write
+
+$$A_{rr}A_{ss} - A_{rs}A_{sr} = \Delta_{rs,rs}\Delta_n$$
+
+where $\Delta_{rs,rs}$ is the determinant of the original matrix with
+the rows $r$ and $s$, as well as the columns $r$ and $s$ removed. The
+equation above simplifies in the skew-symmetric case to
+
+$$A_{rs} = \sqrt{\Delta_{rs,rs}\Delta_n}$$
+
+We now plug this back into the original formula for the determinant,
+
+$$\Delta_n=
+(-1)^{p+1}a_{p1}\sqrt{\Delta_{p1,p1}\Delta_n} +
+(-1)^{p+2}a_{p2}\sqrt{\Delta_{p2,p2}\Delta_n} +
+\cdots +
+(-1)^{n+p}a_{pn}\sqrt{\Delta_{pn,pn}\Delta_n}
+$$
+
+or with slight manipulation,
+
+$$\sqrt{\Delta_n}=(-1)^{p+1}
+\left( \ a_{p1}\sqrt{\Delta_{p1,p1}} -
+a_{p2}\sqrt{\Delta_{p2,p2}} +
+\cdots +
+(-1)^{n-1}a_{pn}\sqrt{\Delta_{pn,pn}} \
+\right)
+$$
+
+The determinant is thus the square of the right hand side, and so we
+identify the right hand side as the Pfaffian.
+
+\section{Identities}
+
+For a $2n \times 2n$ skew-symmetric matrix $A$ and an arbitrary
+$2n \times 2n$ matrix B,
+\begin{itemize}
+\item $Pf(A)^2 = \det(A)$
+\item $Pf(BAB^T) = \det(B)Pf(A)$
+\item $Pf(\lambda{}A) = \lambda^nPf(A)$
+\item $Pf(A^T) = ( - 1)^nPf(A)$
+\item For a block-diagonal matrix
+$$A_1\oplus A_2=
+\left[
+\begin{array}{cc}
+A_1 & 0 \\
+0 & A_2
+\end{array}
+\right]
+$$
+$$Pf(A1\oplus A2) = Pf(A_11)Pf(A_2)$$
+\item For an arbitrary $n \times n$ matrix $M$:
+$$\mbox{Pf}
+\left[
+\begin{array}{cc}
+0 & M \\
+-M^T & 0
+\end{array}
+\right] =
+(-1)^{n(n-1)/2}\det M
+$$
+\end{itemize}
+\section{Applications}
+
+The Pfaffian is an invariant polynomial of a skew-symmetric matrix
+(note that it is not invariant under a general change of basis but
+rather under a proper orthogonal transformation). As such, it is
+important in the theory of characteristic classes. In particular, it
+can be used to define the Euler class of a Riemannian manifold which
+is used in the generalized Gauss-Bonnet theorem.
+
+The number of perfect matchings in a planar graph turns out to be the
+absolute value of a Pfaffian, hence is polynomial time
+computable. This is surprising given that for a general graph, the
+problem is very difficult (so called \#P-complete). This result is used
+to calculate the partition function of Ising models of spin glasses in
+physics, respectively of Markov random fields in machine learning
+(Globerson and Jaakkola, 2007), where the underlying graph is
+planar. Recently it is also used to derive efficient algorithms for
+some otherwise seemingly intractable problems, including the efficient
+simulation of certain types of restricted quantum computation.
+
+The calculation of the number of possible ways to tile a standard
+chessboard or 8-by-8 checkerboard with 32 dominoes is a simple example
+of a problem which may be solved through the use of the Pfaffian
+technique. There are 12,988,816 possible ways to tile a chessboard in
+this manner. Specifically, 12988816 is the number of possible ways to
+cover an 8-by-8 square with 32 1-by-2 rectangles. 12988816 is a square
+number: $12988816 = 3604^2$). Note that 12988816 can be written in the
+form: $2\times 1802^2 + 2\times 1802^2$,
+where all the numbers have a digital root of 2.
+
+More generally, the number of ways to cover a $2n \times 2n$ square with
+$2n^2$ dominoes (as calculated independently by Temperley and M.E. Fisher and
+Kasteleyn in 1961) is given by
+
+$$
+\prod_{j=1}^N
+\prod_{k=1}^N
+\left ( 4\cos^2 \frac{\pi j}{2n + 1} +
+4\cos^2 \frac{\pi k}{2n + 1} \right )
+$$
+
+This technique can be applied in many mathematics-related subjects,
+for example, in the classical, 2-dimensional computation of the
+dimer-dimer correlator function in quantum mechanics.
+
+\section{History}
+
+The term Pfaffian was introduced by Arthur Cayley, who used the term
+in 1852: ``The permutants of this class (from their connection with the
+researches of Pfaff on differential equations) I shall term
+Pfaffians.'' The term honors German mathematician Johann Friedrich
+Pfaff.
+
+\section{Axiom code}
+I have hacked together an algorithm to compute a Pfaffian, using an algorithm
+of Gunter Rote. Currently it's only an .input script, but if it's useful for
+somebody else than myself, we could make it a little more professional.
+
+Martin
+
+<<*>>=
+)spool pfaffian.output
+)set message test on
+)set message auto off
+)clear all
+
+--S 1 of 12
+B0 n == matrix [[(if i=j+1 and odd? j then -1 else _
+ if i=j-1 and odd? i then 1 else 0) _
+ for j in 1..n] for i in 1..n]
+--R
+--R Type:
Void
+--E 1
+
+--S 2 of 12
+PfChar(lambda, A) ==
+ n := nrows A
+ (n = 2) => lambda^2 + A.(1,2)
+ M := subMatrix(A, 3, n, 3, n)
+ r := subMatrix(A, 1, 1, 3, n)
+ s := subMatrix(A, 3, n, 2, 2)
+
+ p := PfChar(lambda, M)
+ d := degree(p, lambda)
+
+ B := B0(n-2)
+ C := r*B
+ g := [(C*s).(1,1), A.(1,2), 1]
+ if d >= 4 then
+ B := M*B
+ for i in 4..d by 2 repeat
+ C := C*B
+ g := cons((C*s).(1,1), g)
+ g := reverse! g
+
+ res := 0
+ for i in 0..d by 2 for j in 2..d+2 repeat
+ c := coefficient(p, lambda, i)
+ for e in first(g, j) for k in 2..-d by -2 repeat
+ res := res + c * e * lambda^(k+i)
+
+ res
+--R
+--R Type:
Void
+--E 2
+
+--S 3 of 12
+pfaffian A == eval(PfChar(l, A), l=0)
+--R
+--R Type:
Void
+--E 3
+
+--S 4 of 12
+m:Matrix(Integer):=[[0,15],[-15,0]]
+--R
+--R + 0 15+
+--R (4) | |
+--R +- 15 0 +
+--R Type: Matrix
Integer
+--E 4
+
+--S 5 of 12
+pfaffian m
+--R
+--R (5) 15
+--R Type: Polynomial
Integer
+--E 5
+
+--S 6 of 12
+m1:Matrix(Polynomial(Integer)):=[[0,a,b,c],[-a,0,d,e],[-b,-d,0,f],[-c,-e,-f,0]]
+--R
+--R
+--R + 0 a b c+
+--R | |
+--R |- a 0 d e|
+--R (6) | |
+--R |- b - d 0 f|
+--R | |
+--R +- c - e - f 0+
+--R Type: Matrix Polynomial
Integer
+--E 6
+
+--S 7 of 12
+pfaffian m1
+--R
+--R Compiling function B0 with type PositiveInteger -> Matrix Integer
+--R
+--R (7) a f - b e + c d
+--R Type: Polynomial
Integer
+--E 7
+
+--S 8 of 12
+(a,b,c,d,e,f):=(3,5,7,11,13,17)
+--R
+--R (8) 17
+--R Type:
PositiveInteger
+--E 8
+
+--S 9 of 12
+m1
+--R
+--R
+--R + 0 a b c+
+--R | |
+--R |- a 0 d e|
+--R (9) | |
+--R |- b - d 0 f|
+--R | |
+--R +- c - e - f 0+
+--R Type: Matrix Polynomial
Integer
+--E 9
+
+--S 10 of 12
+a*f-b*e+d*c
+--R
+--R
+--R (10) 63
+--R Type:
PositiveInteger
+--E 10
+
+--S 11 of 12
+n:=pfaffian m1
+--R
+--R (11) a f - b e + c d
+--R Type: Polynomial
Integer
+--E 9
+
+--S 12 of 12
+eval(n,['a,'b,'c,'d,'e,'f]::List(Symbol),[a,b,c,d,e,f])
+--R
+--R (12) 63
+--R Type: Polynomial
Integer
+--E 12
+)spool
+)lisp (bye)
+@
+\eject
+\begin{thebibliography}{99}
+\bibitem{1} {\bf ``http://en.wikipedia.org/wiki/Pfaffian''}
+\bibitem{2} ``The statistics of dimers on a lattice, Part I'', Physica,
+vol.27, 1961, pp.1209-25, P.W. Kasteleyn.
+\bibitem{3} Propp, James (2004),
+``Lambda-determinants and domino-tilings'', arXiv:math.CO/0406301.
+\bibitem{4} Globerson, Amir and Tommi Jaakkola (2007),
+``Approximate inference using planar graph decomposition'',
+Advances in Neural Information Processing Systems 19, MIT Press.
+\bibitem{5} ``The Games and Puzzles Journal'',
+No.14, 1996, pp.204-5, Robin J. Chapman, University of Exeter
+\bibitem{6} ``Domino Tilings and Products of Fibonacci and Pell numbers'',
+Journal of Integer Sequences, Vol.5, 2002, Article 02.1.2,
+James A. Sellers, The Pennsylvania State University
+\bibitem{7} ``The Penguin Dictionary of Curious and Interesting Numbers'',
+revised ed., 1997, ISBN 0-14-026149-4, David Wells, p.182.
+\bibitem{8} ``A Treatise on the Theory of Determinants'',
+1882, Macmillan and Co., Thomas Muir, Online
+\bibitem{9} ``Skew-Symmetric Determinants'',
+The American Mathematical Monthly, vol. 61, no.2., 1954, p.116,
+S. Parameswaran Online-Subscription
+\end{thebibliography}
+\end{document}
[Prev in Thread] |
Current Thread |
[Next in Thread] |
- [Axiom-developer] 20070927.01.tpd.patch,
daly <=