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}