Skip to content
Find file
Fetching contributors…
Cannot retrieve contributors at this time
400 lines (350 sloc) 11.2 KB
;;; -*- Mode: LISP; Package: (SCREAMS :USE CL :COLON-MODE :EXTERNAL); Base: 10; Syntax: Ansi-common-lisp -*-
;;; LaHaShem HaAretz U'Mloah
(in-package :screamer-user)
(eval-when (:compile-toplevel :load-toplevel)
(require :iterate))
(screamer:define-screamer-package :screams (:use :iterate))
(in-package :screams)
(defvar *infinity* 1d38)
(defvar *fuzz* 1d-6)
(defun make-interval (low high) (cons low high))
(defun low (interval) (car interval))
(defun high (interval) (cdr interval))
(defun size (interval) (- (high interval) (low interval)))
(defun small? (interval) (< (size interval) *fuzz*))
(defun disjoint? (interval1 interval2)
(or (< (high interval1) (low interval2))
(> (low interval1) (high interval2))))
(defun intersects? (interval1 interval2) (not (disjoint? interval1 interval2)))
(defun add-interval (interval1 interval2)
(make-interval (+ (low interval1) (low interval2))
(+ (high interval1) (high interval2))))
(defun multiply-interval (interval1 interval2)
(let ((l1 (low interval1))
(h1 (high interval1))
(l2 (low interval2))
(h2 (high interval2)))
(let ((b1 (* l1 l2))
(b2 (* l1 h2))
(b3 (* h1 l2))
(b4 (* h1 h2)))
(make-interval (min b1 b2 b3 b4) (max b1 b2 b3 b4)))))
(defun compute-interval (expression environment)
(cond ((numberp expression) (make-interval expression expression))
((symbolp expression) (cdr (assoc expression environment)))
((eq (car expression) '+)
(reduce #'add-interval
(mapcar #'(lambda (expression)
(compute-interval expression environment))
(rest expression))))
((eq (car expression) '*)
(reduce #'multiply-interval
(mapcar #'(lambda (expression)
(compute-interval expression environment))
(rest expression))))))
(defun impossible? (equation environment)
(disjoint? (compute-interval (second equation) environment)
(compute-interval (third equation) environment)))
(defun variables-of (expression)
(cond ((numberp expression) nil)
((member expression '(+ * =)) nil)
((null expression) nil)
((atom expression) (list expression))
(t (append (variables-of (car expression))
(variables-of (cdr expression))))))
(defun make-environment (variables interval)
(cond ((null variables) nil)
((member (first variables) (rest variables))
(make-environment (rest variables) interval))
(t (cons (cons (first variables) interval)
(make-environment (rest variables) interval)))))
(defun a-half-interval (interval)
(let ((midpoint (/ (+ (low interval) (high interval)) 2)))
(if (and (> midpoint (low interval))
(< midpoint (high interval)))
(either (make-interval (low interval) midpoint)
(make-interval midpoint (high interval)))
;; Due to floating-point inaccuracy, we can have a big interval
;; without any midpoints in it that can be represented.
(fail))))
(defun copy-all (tree)
(if (atom tree) tree (cons (copy-all (car tree)) (copy-all (cdr tree)))))
(defun biggest-cell (environment)
(if (null (rest environment))
(first environment)
(let ((big-cell (biggest-cell (rest environment))))
(let ((interval1 (cdr (first environment)))
(interval2 (cdr big-cell)))
(if (> (size interval1) (size interval2))
(first environment)
big-cell)))))
(defun refine-environment (environment equations)
(let ((cell (biggest-cell environment)))
(cond ((small? (cdr cell))
(copy-all environment))
(t
(local (setf (cdr cell) (a-half-interval (cdr cell))))
(if (some #'(lambda (equation) (impossible? equation environment))
equations)
(fail)
(refine-environment environment equations))))))
(defun solve (equations)
(refine-environment
(make-environment (variables-of equations)
(make-interval (- *infinity*) *infinity*))
equations))
(defun dam-test ()
(for-effects
(print (solve '((= (+ (* x x) (* y y)) (* z z))
(= x (+ 5 y))
(= z (+ x y 3)))))))
(defun dam-testa ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*))
(y (a-real-betweenv (- *infinity*) *infinity*))
(z (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (andv (=v (+v (*v x x) (*v y y)) (*v z z))
(=v x (+v 5 y))
(=v z (+v x y 3))))
(list x y z))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun eq4 () (for-effects (print (solve '((= (+ x x x x) 3))))))
(defun eq4a ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (=v (+v x x x x) 3.0))
x)
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun eq3 () (for-effects (print (solve '((= (* 4 x) 3))))))
(defun eq2 ()
(for-effects
(print (solve '((= (+ x x y y y) 17)
(= (+ x x x x x x x y y y y) 27))))))
(defun eq2a ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*))
(y (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (andv (=v (+v x x y y y) 17.0)
(=v (+v x x x x x x x y y y y) 27.0)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun eq1 ()
(for-effects
(print (solve '((= (+ (* 2 x) (* 3 y)) 17) (= (+ (* 7 x) (* 4 y)) 27))))))
(defun equation1 ()
(for-effects
(print
(solution
(let ((x (an-integer-betweenv -1000000 1000000))
(y (an-integer-betweenv -1000000 1000000)))
(assert! (andv (=v (+v (*v 2 x) (*v 3 y)) 17)
(=v (+v (*v 7 x) (*v 4 y)) 27)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
;;; FIXME: This blows stack on SBCL.
#+nil
(defun equation1a ()
;; note: This uses ratnums which eventually get turned into flonums.
(for-effects
(print
(solution
(let ((x (a-real-betweenv -1000000 1000000))
(y (a-real-betweenv -1000000 1000000)))
(assert! (andv (=v (+v (*v 2 x) (*v 3 y)) 17)
(=v (+v (*v 7 x) (*v 4 y)) 27)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation1b ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*))
(y (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (andv (=v (+v (*v 2 x) (*v 3 y)) 17)
(=v (+v (*v 7 x) (*v 4 y)) 27)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation1c ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*))
(y (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (andv (=v (+v (*v 2.0 x) (*v 3.0 y)) 17.0)
(=v (+v (*v 7.0 x) (*v 4.0 y)) 27.0)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation2 ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (=v (*v x x) 4.0))
x)
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation3 ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*))
(y (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (andv (=v (+v x x y y y) 17.0)
(=v (+v x x x x x x x y y y y) 27.0)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
;;; FIXME: This doesn't actually get (1 5), but ends up with
;;; variables very close to them -- even though the comment
;;; leads me to think it's supposed to get (1 5).
(defun equation3a ()
;; note: This uses ratnums which eventually get turned into flonums.
(for-effects
(print
(solution
(let ((x (a-real-betweenv -1000000 1000000))
(y (a-real-betweenv -1000000 1000000)))
(assert! (andv (=v (+v x x y y y) 17)
(=v (+v x x x x x x x y y y y) 27)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation3b ()
(for-effects
(print
(solution
(let ((x (an-integer-betweenv (- *infinity*) *infinity*))
(y (an-integer-betweenv (- *infinity*) *infinity*)))
(assert! (andv (=v (+v x x y y y) 17)
(=v (+v x x x x x x x y y y y) 27)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation3c ()
(for-effects
(print
(solution
(let ((x (an-integer-betweenv -1000000 1000000))
(y (an-integer-betweenv -1000000 1000000)))
(assert! (andv (=v (+v x x y y y) 17)
(=v (+v x x x x x x x y y y y) 27)))
(list x y))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation4 ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv -1000000 1000000)))
(assert! (=v (+v x x x x) 8.0))
x)
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation4a ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (=v (+v x x x x) 3))
x)
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation4b ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*)))
(assert! (=v (+v x x x x) 8.0))
x)
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
(defun equation4c ()
(for-effects
(print
(solution
(let ((x (an-integer-betweenv -1000000 1000000)))
(assert! (=v (+v x x x x) 8.0))
x)
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
;;; This doesn't appear to terminate very quickly, if at all...
(defun nonlinear1 ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*))
(y (a-real-betweenv (- *infinity*) *infinity*))
(z (a-real-betweenv (- *infinity*) *infinity*)))
(assert!
(andv (orv (=v (+v (*v 4 x x y) (*v 7 y z z) (*v 6 x x z z)) 2)
(=v (+v (*v 3 x y) (*v 2 y y) (*v 5 x y z)) -4))
(>=v (*v (+v x y) (+v y z)) -5)))
(list x y z))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
;;; FIXME: blows stack on SBCL.
#+nil
(defun nonlinear2 ()
(for-effects
(print
(solution
(let ((x (a-real-betweenv (- *infinity*) *infinity*))
(y (a-real-betweenv (- *infinity*) *infinity*))
(z (a-real-betweenv (- *infinity*) *infinity*)))
(assert!
(andv (=v (+v (*v 4 x x y) (*v 7 y z z) (*v 6 x x z z)) 2)
(=v (+v (*v 3 x y) (*v 2 y y) (*v 5 x y z)) -4)
(>=v (*v (+v x y) (+v y z)) -5)))
(list x y z))
(reorder #'range-size
#'(lambda (x) (< x *fuzz*))
#'>
#'divide-and-conquer-force)))))
;;; Tam V'Nishlam Shevah L'El Borei Olam
Something went wrong with that request. Please try again.