diff --git a/documentation/missing-features.text b/documentation/missing-features.text
index dc9c92be33e2b4aa59228f21bdcb0bac1266e502..2d7faf166b7d4e54262bf72e3d0b287e0edfc20d 100644
--- a/documentation/missing-features.text
+++ b/documentation/missing-features.text
@@ -11,8 +11,6 @@ to GSLL:
 
 GSL 1.9
 http://lists.gnu.org/archive/html/info-gsl/2007-02/msg00000.html
-** Added support for nonsymmetric eigensystems (Patrick Alken)
-** Added Mathieu functions (Lowell Johnson)
 ** Added a new BFGS2 minimisation method, requires substantially fewer
 function and gradient evaluations that the existing BFGS minimiser.
 ** Added updated Knuth generator, gsl_rng_knuthran2002, from 9th
@@ -25,7 +23,6 @@ gsl_matrix_isneg.
 
 GSL 1.10
 http://lists.gnu.org/archive/html/info-gsl/2007-09/msg00000.html
-** Added support for generalized eigensystems (Patrick Alken)
 ** Extended Cholesky routines to complex matrices (Patrick Alken)
 ** Added functions gsl_matrix_subrow and gsl_matrix_subcolumn
 ** Added function gsl_stats_correlation to compute Pearson correlation
diff --git a/documentation/status.text b/documentation/status.text
index 09d580a4cda9f4edf587efc87da0eda903ed411e..e582fd9b1011581e999e1d422998f8e6978acaea 100644
--- a/documentation/status.text
+++ b/documentation/status.text
@@ -52,3 +52,8 @@ Is this intentional?
 Wed Jan 28 2009:
 What are the _e functions in interpolation for, how are the results to
 be interpreted?
+
+=== SBCL problems? ===
+Mon Feb 16 2009
+(copy (mathieu-ce-array 0.0d0 (/ PI 2) 6) 'array)
+gives an error
diff --git a/gsll-tests.asd b/gsll-tests.asd
index c290b6d28ebc5cc71b6c93a962a2008cdd1dd53e..de55e2d3146f6ee759846846dc43b025c206adcb 100644
--- a/gsll-tests.asd
+++ b/gsll-tests.asd
@@ -1,6 +1,6 @@
 ;; Definition of GSLL system 
 ;; Liam Healy
-;; Time-stamp: <2009-02-11 22:15:11EST gsll-tests.asd>
+;; Time-stamp: <2009-02-16 17:57:09EST gsll-tests.asd>
 ;; $Id$
 
 (asdf:defsystem "gsll-tests"
@@ -90,6 +90,7 @@
 	     (:file "lognormal")
 	     (:file "lu")
 	     (:file "mathematical")
+	     (:file "mathieu")
 	     (:file "matrix-copy")
 	     (:file "matrix-copy-to-cl-and-back")
 	     (:file "matrix-div")
diff --git a/gsll.asd b/gsll.asd
index 22048c6cf2fff5fd1d457e08f340ae8b2154ce72..2b78987e739c14bcb3f399505fe850ad323c8c5a 100644
--- a/gsll.asd
+++ b/gsll.asd
@@ -1,6 +1,6 @@
 ;; Definition of GSLL system 
 ;; Liam Healy
-;; Time-stamp: <2009-02-16 14:36:41EST gsll.asd>
+;; Time-stamp: <2009-02-16 17:15:57EST gsll.asd>
 ;; $Id$
 
 (asdf:defsystem "gsll"
@@ -78,6 +78,7 @@
 	       (:file "lambert" :depends-on (return-structures))
 	       (:file "legendre" :depends-on (return-structures))
 	       (:file "logarithm" :depends-on (return-structures))
+	       (:file "mathieu" :depends-on (return-structures))
 	       (:file "power" :depends-on (return-structures))
 	       (:file "psi" :depends-on (return-structures))
 	       (:file "synchrotron" :depends-on (return-structures))
diff --git a/special-functions/mathieu.lisp b/special-functions/mathieu.lisp
new file mode 100644
index 0000000000000000000000000000000000000000..9ad7786019f52436c32ca088001596f1ac5860c0
--- /dev/null
+++ b/special-functions/mathieu.lisp
@@ -0,0 +1,280 @@
+;; Mathieu functions
+;; Liam Healy 2009-02-16 16:30:59EST mathieu.lisp
+;; Time-stamp: <2009-02-16 18:13:44EST mathieu.lisp>
+;; $Id: $
+
+(in-package :gsl)
+
+;;;;****************************************************************************
+;;;; Mathieu workspace
+;;;;****************************************************************************
+
+(defmobject mathieu "gsl_sf_mathieu"
+  ((n sizet) (qmax :double))
+  "workspace for Mathieu functions"
+  :gsl-version (1 9)
+  :documentation "Make a workspace needed for some Mathieu functions.")
+
+;;;;****************************************************************************
+;;;; Characteristic values
+;;;;****************************************************************************
+
+(defmfun mathieu-a (n q)
+  "gsl_sf_mathieu_a"
+  ((n :int) (q :double) (ret sf-result))
+  :gsl-version (1 9)
+  :documentation			; FDL
+  "Compute the characteristic value a_n(q) of the Mathieu function
+  ce_n(q,x) respectively.")
+
+(defmfun mathieu-b (n q)
+  "gsl_sf_mathieu_b"
+  ((n :int) (q :double) (ret sf-result))
+  :gsl-version (1 9)
+  :documentation			; FDL
+  "Compute the characteristic values b_n(q) of the Mathieu
+  function se_n(q,x), respectively.")
+
+(defmfun mathieu-a-array
+    (q &optional
+       size-or-array
+       (minimum-order 0)
+       (workspace
+	(make-mathieu (+ minimum-order (vdf-size size-or-array)) q))
+       &aux (size (vdf-size size-or-array)) (result-array (vdf size)))
+  "gsl_sf_mathieu_a_array"
+  ((minimum-order :int) ((+ minimum-order size) :int) (q :double)
+   ((mpointer workspace) :pointer) ((c-pointer result-array) :pointer))
+  :gsl-version (1 9)
+  :outputs (result-array)
+  :return (result-array)
+  :documentation			; FDL
+  "Compute a series of Mathieu characteristic values a_n(q) for n from
+    minimum-order to minimum-order + size inclusive, where size
+    is either the numerical value supplied in size-or-array, or the
+    the length of the vector supplied there.")
+
+(defmfun mathieu-b-array
+    (q &optional
+       size-or-array
+       (minimum-order 0)
+       (workspace
+	(make-mathieu (+ minimum-order (vdf-size size-or-array)) q))
+       &aux (size (vdf-size size-or-array)) (result-array (vdf size)))
+  "gsl_sf_mathieu_b_array"
+  ((minimum-order :int) ((+ minimum-order size) :int) (q :double)
+   ((mpointer workspace) :pointer) ((c-pointer result-array) :pointer))
+  :gsl-version (1 9)
+  :outputs (result-array)
+  :return (result-array)
+  :documentation			; FDL
+  "Compute a series of Mathieu characteristic values b_n(q) for n from
+    minimum-order to minimum-order + size inclusive, where size
+    is either the numerical value supplied in size-or-array, or the
+    the length of the vector supplied there.")
+
+;;;;****************************************************************************
+;;;; Angular Mathieu functions
+;;;;****************************************************************************
+
+(defmfun mathieu-ce (n q x)
+  "gsl_sf_mathieu_ce"
+  ((n :int) (q :double) (x :double) (ret sf-result))
+  :gsl-version (1 9)
+  :documentation			; FDL
+  "Compute the angular Mathieu functions ce_n(q,x).")
+
+(defmfun mathieu-se (n q x)
+  "gsl_sf_mathieu_se"
+  ((n :int) (q :double) (x :double) (ret sf-result))
+  :gsl-version (1 9)
+  :documentation			; FDL
+  "Compute the angular Mathieu functions se_n(q,x).")
+
+(defmfun mathieu-ce-array
+    (q x &optional
+       size-or-array
+       (minimum-order 0)
+       (workspace
+	(make-mathieu (+ minimum-order (vdf-size size-or-array)) q))
+       &aux (size (vdf-size size-or-array)) (result-array (vdf size)))
+  "gsl_sf_mathieu_ce_array"
+  ((minimum-order :int) ((+ minimum-order size) :int)
+   (q :double) (x :double)
+   ((mpointer workspace) :pointer) ((c-pointer result-array) :pointer))
+  :gsl-version (1 9)
+  :outputs (result-array)
+  :return (result-array)
+  :documentation			; FDL
+  "Compute a series of the angular Mathieu function ce_n(q) for n from
+    minimum-order to minimum-order + size inclusive, where size
+    is either the numerical value supplied in size-or-array, or the
+    the length of the vector supplied there.")
+
+(defmfun mathieu-se-array
+    (q x &optional
+       size-or-array
+       (minimum-order 0)
+       (workspace
+	(make-mathieu (+ minimum-order (vdf-size size-or-array)) q))
+       &aux (size (vdf-size size-or-array)) (result-array (vdf size)))
+  "gsl_sf_mathieu_se_array"
+  ((minimum-order :int) ((+ minimum-order size) :int)
+   (q :double) (x :double)
+   ((mpointer workspace) :pointer) ((c-pointer result-array) :pointer))
+  :gsl-version (1 9)
+  :outputs (result-array)
+  :return (result-array)
+  :documentation			; FDL
+  "Compute a series of the angular Mathieu function se_n(q) for n from
+    minimum-order to minimum-order + size inclusive, where size
+    is either the numerical value supplied in size-or-array, or the
+    the length of the vector supplied there.")
+
+;;;;****************************************************************************
+;;;; Radial Mathieu functions
+;;;;****************************************************************************
+
+(defmfun mathieu-Mc (j n q x)
+  "gsl_sf_mathieu_Mc"
+  ((j :int) (n :int) (q :double) (x :double) (ret sf-result))
+  :gsl-version (1 9)
+  :documentation			; FDL
+  "Compute the radial j-th kind Mathieu functions Mc_n^{(j)}(q,x) of
+   order n. The allowed values of j are 1 and 2. The functions for j =
+   3,4 can be computed as M_n^{(3)} = M_n^{(1)} + iM_n^{(2)} and
+   M_n^{(4)} = M_n^{(1)} - iM_n^{(2)}, where M_n^{(j)} = Mc_n^{(j)} or
+   Ms_n^{(j)}.")
+
+(defmfun mathieu-Ms (j n q x)
+  "gsl_sf_mathieu_Ms"
+  ((j :int) (n :int) (q :double) (x :double) (ret sf-result))
+  :gsl-version (1 9)
+  :documentation			; FDL
+  "Compute the radial j-th kind Mathieu functions Ms_n^{(j)}(q,x) of order n.
+   The allowed values of j are 1 and 2. The functions for j = 3,4 can
+   be computed as M_n^{(3)} = M_n^{(1)} + iM_n^{(2)} and M_n^{(4)} =
+   M_n^{(1)} - iM_n^{(2)}, where M_n^{(j)} = Mc_n^{(j)} or
+   Ms_n^{(j)}.")
+
+(defmfun mathieu-Mc-array
+    (j q x &optional
+       size-or-array
+       (minimum-order 0)
+       (workspace
+	(make-mathieu (+ minimum-order (vdf-size size-or-array)) q))
+       &aux (size (vdf-size size-or-array)) (result-array (vdf size)))
+  "gsl_sf_mathieu_Mc_array"
+  ((j :int) (minimum-order :int) ((+ minimum-order size) :int)
+   (q :double) (x :double)
+   ((mpointer workspace) :pointer) ((c-pointer result-array) :pointer))
+  :gsl-version (1 9)
+  :outputs (result-array)
+  :return (result-array)
+  :documentation			; FDL
+  "Compute a series of the radial Mathieu function of kind j for n from
+    minimum-order to minimum-order + size inclusive, where size
+    is either the numerical value supplied in size-or-array, or the
+    the length of the vector supplied there.")
+
+(defmfun mathieu-Ms-array
+    (j q x &optional
+       size-or-array
+       (minimum-order 0)
+       (workspace
+	(make-mathieu (+ minimum-order (vdf-size size-or-array)) q))
+       &aux (size (vdf-size size-or-array)) (result-array (vdf size)))
+  "gsl_sf_mathieu_Ms_array"
+  ((j :int) (minimum-order :int) ((+ minimum-order size) :int)
+   (q :double) (x :double)
+   ((mpointer workspace) :pointer) ((c-pointer result-array) :pointer))
+  :gsl-version (1 9)
+  :outputs (result-array)
+  :return (result-array)
+  :documentation			; FDL
+  "Compute a series of the radial Mathieu function of kind j for n from
+    minimum-order to minimum-order + size inclusive, where size
+    is either the numerical value supplied in size-or-array, or the
+    the length of the vector supplied there.")
+
+;;;;****************************************************************************
+;;;; Examples
+;;;;****************************************************************************
+
+;;; Tests from gsl-1.11/specfunc/test_mathieu.c
+(save-test mathieu
+	   (mathieu-ce 0 0.0d0 0.0d0)
+	   (mathieu-ce 0 0.0d0 (/ pi 2))
+	   (mathieu-ce 0 5.0d0 0.0d0)
+	   (mathieu-ce 0 5.0d0 (/ pi 2))
+	   (mathieu-ce 0 10.0d0 0.0d0)
+	   (mathieu-ce 0 10.0d0 (/ pi 2))
+	   (mathieu-ce 0 15.0d0 0.0d0)
+	   (mathieu-ce 0 15.0d0 (/ pi 2))
+	   (mathieu-ce 0 20.0d0 0.0d0)
+	   (mathieu-ce 0 20.0d0 (/ pi 2))
+	   (mathieu-ce 0 25.0d0 0.0d0)
+	   (mathieu-ce 0 25.0d0 (/ pi 2))
+	   (mathieu-ce 1 0.0d0 0.0d0)
+	   (mathieu-ce 1 5.0d0 0.0d0)
+	   (mathieu-ce 1 10.0d0 0.0d0)
+	   (mathieu-ce 1 15.0d0 0.0d0)
+	   (mathieu-ce 1 20.0d0 0.0d0)
+	   (mathieu-ce 1 25.0d0 0.0d0)
+	   (mathieu-se 1 0.0d0 (/ pi 2))
+	   (mathieu-se 1 5.0d0 (/ pi 2))
+	   (mathieu-se 1 10.0d0 (/ pi 2))
+	   (mathieu-se 1 15.0d0 (/ pi 2))
+	   (mathieu-se 1 20.0d0 (/ pi 2))
+	   (mathieu-se 1 25.0d0 (/ pi 2))
+	   (mathieu-ce 2 0.0d0 0.0d0)
+	   (mathieu-ce 2 0.0d0 (/ pi 2))
+	   (mathieu-ce 2 5.0d0 0.0d0)
+	   (mathieu-ce 2 5.0d0 (/ pi 2))
+	   (mathieu-ce 2 10.0d0 0.0d0)
+	   (mathieu-ce 2 10.0d0 (/ pi 2))
+	   (mathieu-ce 2 15.0d0 0.0d0)
+	   (mathieu-ce 2 15.0d0 (/ pi 2))
+	   (mathieu-ce 2 20.0d0 0.0d0)
+	   (mathieu-ce 2 20.0d0 (/ pi 2))
+	   (mathieu-ce 2 25.0d0 0.0d0)
+	   (mathieu-ce 2 25.0d0 (/ pi 2))
+	   (mathieu-ce 5 0.0d0 0.0d0)
+	   (mathieu-ce 5 5.0d0 0.0d0)
+	   (mathieu-ce 5 10.0d0 0.0d0)
+	   (mathieu-ce 5 15.0d0 0.0d0)
+	   (mathieu-ce 5 20.0d0 0.0d0)
+	   (mathieu-ce 5 25.0d0 0.0d0)
+	   (mathieu-se 5 0.0d0 (/ pi 2))
+	   (mathieu-se 5 5.0d0 (/ pi 2))
+	   (mathieu-se 5 10.0d0 (/ pi 2))
+	   (mathieu-se 5 15.0d0 (/ pi 2))
+	   (mathieu-se 5 20.0d0 (/ pi 2))
+	   (mathieu-se 5 25.0d0 (/ pi 2))
+	   (mathieu-ce 10 0.0d0 0.0d0)
+	   (mathieu-ce 10 0.0d0 (/ pi 2))
+	   (mathieu-ce 10 5.0d0 0.0d0)
+	   (mathieu-ce 10 5.0d0 (/ pi 2))
+	   (mathieu-ce 10 10.0d0 0.0d0)
+	   (mathieu-ce 10 10.0d0 (/ pi 2))
+	   (mathieu-ce 10 15.0d0 0.0d0)
+	   (mathieu-ce 10 15.0d0 (/ pi 2))
+	   (mathieu-ce 10 20.0d0 0.0d0)
+	   (mathieu-ce 10 20.0d0 (/ pi 2))
+	   (mathieu-ce 10 25.0d0 0.0d0)
+	   (mathieu-ce 10 25.0d0 (/ pi 2))
+	   (mathieu-ce 15 0.0d0 0.0d0)
+	   (mathieu-ce 15 5.0d0 0.0d0)
+	   (mathieu-ce 15 10.0d0 0.0d0)
+	   (mathieu-ce 15 15.0d0 0.0d0)
+	   (mathieu-ce 15 20.0d0 0.0d0)
+	   (mathieu-ce 15 25.0d0 0.0d0)
+	   (mathieu-se 15 0.0d0 (/ pi 2))
+	   (mathieu-se 15 5.0d0 (/ pi 2))
+	   (mathieu-se 15 10.0d0 (/ pi 2))
+	   (mathieu-se 15 15.0d0 (/ pi 2))
+	   (mathieu-se 15 20.0d0 (/ pi 2))
+	   (mathieu-se 15 25.0d0 (/ pi 2))
+	   (cl-array (mathieu-ce-array 0.0d0 (/ pi 2) 6))
+	   (cl-array (mathieu-ce-array 20.0d0 0.0d0 16))
+	   (cl-array (mathieu-se-array 20.0d0 (/ pi 2) 15 1)))
diff --git a/special-functions/return-structures.lisp b/special-functions/return-structures.lisp
index 61f2d35a619be8d5db796a6fe0ace810291589ce..d683e44a0a350dab762b68afd3884e682b2112cf 100644
--- a/special-functions/return-structures.lisp
+++ b/special-functions/return-structures.lisp
@@ -1,6 +1,6 @@
 ;; Structures returned by special functions.
 ;; Liam Healy, Mon Jan  1 2007 - 11:35
-;; Time-stamp: <2009-01-10 17:21:24EST return-structures.lisp>
+;; Time-stamp: <2009-02-16 16:49:45EST return-structures.lisp>
 ;; $Id$
 
 (in-package :gsl)
@@ -51,3 +51,9 @@
   (if (integerp size-or-array)
       (make-marray 'double-float :dimensions size-or-array)
       size-or-array))
+
+(defun vdf-size (size-or-array)
+  "Make or take a vector-double-float."
+  (if (integerp size-or-array)
+      size-or-array
+      (total-size size-or-array)))
diff --git a/tests/mathieu.lisp b/tests/mathieu.lisp
new file mode 100644
index 0000000000000000000000000000000000000000..fcbcc0560a1228bc63d41d981b3b21c085325fef
--- /dev/null
+++ b/tests/mathieu.lisp
@@ -0,0 +1,253 @@
+;; Regression test MATHIEU for GSLL, automatically generated
+
+(in-package :gsl)
+
+(LISP-UNIT:DEFINE-TEST MATHIEU
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.7071067811865475d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 0.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.7071067811865475d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 0.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.04480018165188903d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 5.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.334848674698019d0 5.927918932160465d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 5.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.007626517570935777d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 10.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.4686604707128563d0 6.522162679768934d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 10.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.0019325083152045943d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 15.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.5501081466866489d0 6.883863020442189d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 15.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 6.037438292241945d-4 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 20.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.609890857395926d0 7.149351588057966d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 20.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 2.1586301841465958d-4 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 25.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.657510298323475d0 7.360824387008136d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 0 25.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 1 0.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.25654287932236375d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 1 5.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.053598747747176455d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 1 10.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.015040066453826315d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 1 15.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.005051813764713131d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 1 20.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.0019110515066575852d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 1 25.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 1 0.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.3374338870223457d0 5.939399581144514d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 1 5.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.4687556641029376d0 6.522585423342775d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 1 10.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.5501150743575518d0 6.883893785481162d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 1 15.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.6098915926037722d0 7.149354853036681d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 1 20.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.6575103983745163d0 7.360824831324015d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 1 25.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 0.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 0.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.7352943084006845d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 5.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.7244881519676682d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 5.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.24588834929131909d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 10.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.9267592641263209d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 10.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.0787928278463933d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 15.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -1.019966226030262d0 4.529559953915294d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 15.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.028648943147074366d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 20.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -1.0752932287796875d0 4.77526120325894d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 20.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.011512866330887451d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 25.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -1.116278953295253d0 4.957274383411439d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 2 25.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 5 0.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.1248072506384799d0 4.995147631696639d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 5 5.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.2580199413082875d0 5.586730817112196d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 5 10.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.193432230413072d0 5.29990376213739d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 5 15.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.936575531422621d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 5 20.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.6106943100506989d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 5 25.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 5 0.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.906077930202355d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 5 5.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.8460384335355104d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 5 10.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.8379493400124841d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 5 15.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.8635431218533666d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 5 20.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 0.8992683245108409d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 5 25.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 0.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 0.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0259950270894385d0 4.556333208902423d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 5.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.975347487235964d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 5.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0538159921009345d0 4.679883112594638d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 10.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.9516453181789555d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 10.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0841063118392213d0 4.814399154181454d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 15.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.9285480638845388d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 15.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.1177886312593968d0 4.963978700353685d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 20.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.9057107845940974d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 20.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.1562399186322392d0 5.134736718624918d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 25.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.8826919105636902d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 10 25.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 15 0.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0112937325295657d0 4.491046346053754d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 15 5.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.022878282438181d0 4.54249208220761d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 15 10.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0347936522368735d0 4.595406953797336d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 15 15.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0470843441628868d0 4.649988590456674d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 15 20.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST 1.0598004418139373d0 4.706459408038987d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-CE 15 25.0d0 0.0d0)))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -1.0d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 15 0.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.9889607027406359d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 15 5.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.978142347183216d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 15 10.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.9675137031854539d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 15 15.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.9570452540612818d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 15 20.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST -0.9467086958780897d0 4.440892098500626d-16)
+                        (MULTIPLE-VALUE-LIST (MATHIEU-SE 15 25.0d0 (/ PI 2))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST
+                         #(0.7071067811865475d0 6.123233995736766d-17 -1.0d0
+                           -1.8369701987210297d-16 1.0d0
+                           3.061616997868383d-16))
+                        (MULTIPLE-VALUE-LIST
+                         (CL-ARRAY (MATHIEU-CE-ARRAY 0.0d0 (/ PI 2) 6))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST
+                         #(6.03743829223914d-4 0.005051813764712108d0
+                           0.02864894314707411d0 0.12611832163429212d0
+                           0.43023077677357696d0 0.9365755314225924d0
+                           1.2761355598483817d0 1.3076504092104844d0
+                           1.2172162968773408d0 1.154852006191031d0
+                           1.117788631259397d0 1.0936613176252825d0
+                           1.0767085541408374d0 1.0641859681403019d0
+                           1.054606963632144d0 1.0470843441628828d0))
+                        (MULTIPLE-VALUE-LIST
+                         (CL-ARRAY (MATHIEU-CE-ARRAY 20.0d0 0.0d0 16))))
+                       (LISP-UNIT::ASSERT-NUMERICAL-EQUAL
+                        (LIST
+                         #(1.6098915926037731d0 3.974439518389252d-16
+                           -1.0759417465195853d0 -4.3771386586551907d-16
+                           0.8635431218533666d0 4.514819072643059d-16
+                           -0.8254327649109297d0 -5.435816907711848d-16
+                           0.8843014096898285d0 6.616729635144726d-16
+                           -0.921662441765126d0 -7.791496473385796d-16
+                           0.9433025704899831d0 8.979725089944524d-16
+                           -0.9570452540612865d0))
+                        (MULTIPLE-VALUE-LIST
+                         (CL-ARRAY (MATHIEU-SE-ARRAY 20.0d0 (/ PI 2) 15 1)))))
+