Commit c2a88130 authored by Liam Healy's avatar Liam Healy
Browse files

Householder: add tests from GSL

Added the standard set of linear algebra tests for the Householder
transformation; all pass.
parent df63f1b6
;; Definition of GSLL system
;; Liam Healy
;; Time-stamp: <2009-09-20 21:52:19EDT gsll-tests.asd>
;; Time-stamp: <2009-09-24 22:22:50EDT gsll-tests.asd>
(asdf:defsystem "gsll-tests"
:name "gsll-tests"
......@@ -70,6 +70,7 @@
(:file "hankel")
(:file "higher-moments")
(:file "histogram")
(:file "householder")
(:file "hypergeometric")
(:file "hypergeometric-randist")
(:file "index-max")
......
;; Householder Transformations
;; Liam Healy, Wed May 10 2006 - 10:03
;; Time-stamp: <2009-02-23 21:09:50EST householder.lisp>
;; Time-stamp: <2009-09-24 22:18:00EDT householder.lisp>
;; $Id$
;;; For householder-transform, it would be nice to be able to pick the
......@@ -91,3 +91,23 @@
stored in x and b is not modified. The matrix A is destroyed by
the Householder transformations. The solution is returned from the
function call.")
;;;;****************************************************************************
;;;; Examples and unit test, from linalg/test.c
;;;;****************************************************************************
(defun test-hh-solve-dim (matrix)
"Solve the linear equation using Householder with the supplied
matrix and a right-hand side vector which is the reciprocal of one
more than the index."
(householder-solve (copy matrix) (create-rhs-vector (dim0 matrix)) T))
(save-test householder
(test-hh-solve-dim *hilb2*)
(test-hh-solve-dim *hilb3*)
(test-hh-solve-dim *hilb4*)
(test-hh-solve-dim *hilb12*)
(test-hh-solve-dim *vander2*)
(test-hh-solve-dim *vander3*)
(test-hh-solve-dim *vander4*)
(test-hh-solve-dim *vander12*))
;; Regression test Householder for GSLL.
(in-package :gsl)
;;; Answers inserted from linalg/test.c
;;; GSL has #define GSL_DBL_EPSILON 2.2204460492503131e-16
;;; which is 2x what double-float-epsilon is.
(LISP-UNIT:DEFINE-TEST householder
(let ((lisp-unit:*epsilon* (* 2 16 double-float-epsilon)))
(LISP-UNIT:ASSERT-NUMERICAL-EQUAL
(LIST *hilb2-soln*)
(MULTIPLE-VALUE-LIST
(TEST-HH-SOLVE-DIM *HILB2*))))
(let ((lisp-unit:*epsilon* (* 2 128 double-float-epsilon)))
(LISP-UNIT:ASSERT-NUMERICAL-EQUAL
(LIST *hilb3-soln*)
(MULTIPLE-VALUE-LIST
(TEST-HH-SOLVE-DIM *HILB3*))))
(let ((lisp-unit:*epsilon* (* 2 2048 double-float-epsilon)))
(LISP-UNIT:ASSERT-NUMERICAL-EQUAL
(LIST *hilb4-soln*)
(MULTIPLE-VALUE-LIST (TEST-HH-SOLVE-DIM *hilb4*))))
(let ((lisp-unit:*epsilon* 0.5d0))
(LISP-UNIT:ASSERT-NUMERICAL-EQUAL
(LIST *hilb12-soln*)
(MULTIPLE-VALUE-LIST (TEST-HH-SOLVE-DIM *hilb12*))))
(let ((lisp-unit:*epsilon* (* 2 64 double-float-epsilon)))
(LISP-UNIT:ASSERT-NUMERICAL-EQUAL
(LIST *vander2-soln*)
(MULTIPLE-VALUE-LIST (TEST-HH-SOLVE-DIM *vander2*))))
(let ((lisp-unit:*epsilon* (* 2 64 double-float-epsilon)))
(LISP-UNIT:ASSERT-NUMERICAL-EQUAL
(LIST *vander3-soln*)
(MULTIPLE-VALUE-LIST (TEST-HH-SOLVE-DIM *vander3*))))
(let ((lisp-unit:*epsilon* (* 2 1024 double-float-epsilon)))
(LISP-UNIT:ASSERT-NUMERICAL-EQUAL
(LIST *vander4-soln*)
(MULTIPLE-VALUE-LIST (TEST-HH-SOLVE-DIM *vander4*))))
(let ((lisp-unit:*epsilon* 0.05d0))
(LISP-UNIT:ASSERT-NUMERICAL-EQUAL
(LIST *vander12-soln*)
(MULTIPLE-VALUE-LIST (TEST-HH-SOLVE-DIM *vander12*)))))
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment