source: branches/f4grobner/5am-poly.lisp@ 4508

;;; -*-  Mode: Lisp -*- 
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;                                                                              
;;;  Copyright (C) 1999, 2002, 2009, 2015 Marek Rychlik <rychlik@u.arizona.edu>          
;;;                                                                              
;;;  This program is free software; you can redistribute it and/or modify        
;;;  it under the terms of the GNU General Public License as published by        
;;;  the Free Software Foundation; either version 2 of the License, or           
;;;  (at your option) any later version.                                         
;;;                                                                              
;;;  This program is distributed in the hope that it will be useful,             
;;;  but WITHOUT ANY WARRANTY; without even the implied warranty of              
;;;  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               
;;;  GNU General Public License for more details.                                
;;;                                                                              
;;;  You should have received a copy of the GNU General Public License           
;;;  along with this program; if not, write to the Free Software                 
;;;  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.  
;;;                                                                              
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;
;; Run tests using 5am unit testing framework
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

;; We assume that QuickLisp package manager is installed.
;; See :
;;      https://www.quicklisp.org/beta/
;;

;; The following is unnecessary after running:
;; * (ql:add-to-init-file)
;; at lisp prompt:
;;(load "~/quicklisp/setup")

(ql:quickload :fiveam)

(defpackage #:5am-poly
  (:use :cl :it.bese.fiveam :monom :polynomial :ring))

(in-package :5am-poly)

(def-suite poly-suite 
    :description "Polynomial package suite")

(in-suite poly-suite)

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;
;;        POLY class tests
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

(def-fixture poly-add-context ()
  (let ((p (make-instance 'poly))
        (q (make-instance 'poly :order nil))
        (p+q (make-instance 'poly))
        (p-q (make-instance 'poly))
        (p-uminus (make-instance 'poly)))
    ;; Populate the polynomials; the lists of (exponents . coefficient) pairs
    ;; must be in increasing order in Q, but Q is unordered (:ORDER NIL)
    ;; so it will be automatically sorted.
    (dolist (x '( ((2) . 22)  ((4) . 44) ((5) . 55) ((8) . 88) ((9) . 99) ))
      (poly-insert-term p (make-instance 'term :exponents (car x) :coeff (cdr x))))
    (dolist (x '( ((9) . 90) ((0) . 11)  ((2) . 20) ((3) . 33) ((4) . -44)  ((7) . 77) ((8) . 88) ))
      (poly-insert-term q (make-instance 'term :exponents (car x) :coeff (cdr x))))
    ;; P+Q
    (dolist (x '(((0) . 11) ((2) . 42)  ((3) . 33) ((5) . 55) ((7) . 77) ((8) . 176) ((9) . 189) ))
      (poly-insert-term p+q (make-instance 'term :exponents (car x) :coeff (cdr x))))
    ;; P-Q
    (dolist (x '(((0) . -11) ((2) . 2)  ((3) . -33) ((4) . 88) ((5) . 55) ((7) . -77) ((9) . 9)))
      (poly-insert-term p-q (make-instance 'term :exponents (car x) :coeff (cdr x))))
    ;; -P
    (dolist (x '( ((2) . -22)  ((4) . -44) ((5) . -55) ((8) . -88) ((9) . -99)))
      (poly-insert-term p-uminus (make-instance 'term :exponents (car x) :coeff (cdr x))))
    ;;(print p) (print q) (print p+q) (print p-q)
    (&body)))

(test poly-add
  "Polynomial addition"
  (with-fixture poly-add-context () (is (universal-equalp (add-to p q) p+q)))
  (with-fixture poly-add-context () (is (universal-equalp (add p q) p+q)))
  (with-fixture poly-add-context () (is (universal-equalp (subtract-from p q) p-q)))
  (with-fixture poly-add-context () (is (universal-equalp (subtract p q) p-q)))
  (with-fixture poly-add-context () (is (universal-equalp (unary-minus p) p-uminus)))
  )

(test poly-add-more
  "More polynomial addition"
  (is (equal (poly->alist (add-to 
                           (alist->poly '(((0) . 1)))
                           (alist->poly '(((1) . 2)))))
             '(((1) . 2) ((0) . 1))))
  (is (equal (poly->alist (add-to 
                           (alist->poly '(((0) . 1)))
                           (alist->poly '(((0) . 2)))))
           '(((0) . 3))))
  (is (equal (poly->alist (add-to 
                           (alist->poly '(((0) . 1) ((2) . 5)))
                           (alist->poly '(((1) . 2) ((3) . 7)))))
             '(((3) . 7) ((2) . 5) ((1) . 2) ((0) . 1)))))

(def-fixture poly-multiply-context ()
  (let ((p (make-instance 'poly))
        (q (make-instance 'poly :order nil))
        (p*q (make-instance 'poly)))
    ;; Populate the polynomials; the lists of (exponents . coefficient) pairs
    ;; must be in increasing order in Q, but Q is unordered (:ORDER NIL)
    ;; so it will be automatically sorted.
    (dolist (x '( ((0) . 1)  ((1) . 2) ))
      (poly-insert-term p (make-instance 'term :exponents (car x) :coeff (cdr x))))
    (dolist (x '( ((0) . 1)  ((1) . 3) ))
      (poly-insert-term q (make-instance 'term :exponents (car x) :coeff (cdr x))))
    ;; P*Q
    (dolist (x '( ((0) . 1) ((1) . 5) ((2) . 6)))
      (poly-insert-term p*q (make-instance 'term :exponents (car x) :coeff (cdr x))))
    (&body)))


(test poly-multiply
  "Polynomial multiplication"
  (with-fixture poly-multiply-context () (is (universal-equalp (multiply p q) p*q)))
  )

(test poly-standard-extension
  "Standard extension"
  (let* ((p (alist->poly '( ((0) . 1) ((1) . 2))))
         (q (alist->poly '( ((0) . 1) ((2) . 3))))
         (plist (list p q))
         (p-ext (alist->poly '( ((1 0 0) . 1) ((1 0 1) . 2))))
         (q-ext (alist->poly '( ((0 1 0) . 1) ((0 1 2) . 3))))
         (plist-st-ext (list p-ext q-ext)))
    (is (universal-equalp (standard-extension plist) plist-st-ext))))

(test poly-standard-extension-1
  "Standard extension 1"
  (let* ((p (alist->poly '( ((0) . 1) ((1) . 2))))
         (q (alist->poly '( ((0) . 1) ((2) . 3))))
         (plist (list p q))
         (p-ext (alist->poly '( ((0 0 0) . -1) ((1 0 0) . 1) ((1 0 1) . 2))))
         (q-ext (alist->poly '( ((0 0 0) . -1) ((0 1 0) . 1) ((0 1 2) . 3))))
         (plist-st-ext (list p-ext q-ext)))
    (is (universal-equalp (standard-extension-1 plist) plist-st-ext))))

(test poly-standard-sum
  "Standard sum"
  (let* ((p (alist->poly '( ((0) . 1) ((1) . 2))))
         (q (alist->poly '( ((0) . 1) ((2) . 3))))
         (plist (list p q))
         (std-sum (alist->poly '(((0 0 0) . -1) ((0 1 0) . 1) ((0 1 2) . 3)
                                 ((1 0 0) . 1) ((1 0 1) . 2)))))
    (is (universal-equalp (standard-sum plist) std-sum))))

(test poly-s-polynomial
  "S-Polynomial"
  (let* ((f (alist->poly '( ((1 0 0) . 1) ((1 0 1) . 2)))) ;x+2*x*z |*y*z
         (g (alist->poly '( ((0 1 0) . 1) ((0 1 2) . 2)))) ;y+2*y*z^2 |*x
         (s-poly (alist->poly '( ((1 1 0) . -1) ((1 1 1) . 1)))); x*y*z - x*y
         )
    (is (universal-equalp (s-polynomial f g) s-poly))))

(test poly-s-polynomial-alt
  "S-Polynomial with non-coprime coeffs"
  (let* ((f (alist->poly '( ((1 0 0) . 1) ((1 0 1) . 14)))) ;x+14*x*z |*3*y*z
         (g (alist->poly '( ((0 1 0) . 1) ((0 1 2) . 21)))) ;y+21*y*z^2 |*2*x
         (s-poly (alist->poly '( ((1 1 0) . -2) ((1 1 1) . 3)))) ;3*x*y*z - 2*x*y
         )
    (is (universal-equalp (s-polynomial f g) s-poly))))

(test poly-content
  "Poly-content"
  (let* ((p (alist->poly '( ((1 0 0) . 12) ((1 0 1) . 15)))) 
         (pc (make-instance 'integer-ring :value 3)))
    (is (universal-equalp (poly-content p) pc))))

(test poly-primitive-part
  "Poly-primitive-part"
  (let* ((p (alist->poly '( ((1 0 0) . 12) ((1 0 1) . 15)))) 
         (pp (alist->poly '( ((1 0 0) . 4) ((1 0 1) . 5)))))
    (is (universal-equalp (poly-primitive-part p) pp))))

(test saturation-extension
  "Saturation-extension"
  (let* ((f (list (alist->poly '( ((1 0 0) . 2) ((1 0 1) . 3))) ;2*x+3*x*z
                  (alist->poly '( ((1 0 0) . 5) ((0 1 2) . 7)))) ;5*x+7*y*z^2
           )
         (p (alist->poly '( ((1 1 1) . 11) ((2 3 4) . 13)))) ;11*x*y*z+13*x^2*y^3*z^4
         (sat-ext (list (alist->poly '( ((0 1 0 0) . 2) ((0 1 0 1) . 3)))
                        (alist->poly '( ((0 1 0 0) . 5) ((0 0 1 2) . 7)))
                        (alist->poly '( ((0 0 0 0) . -1) ((1 1 1 1) . 11) ((1 2 3 4) . 13) )))))
    (is (universal-equalp (saturation-extension f (list p)) sat-ext))))

(test saturation-extension-1
  "Saturation-extension-1"
  (let* ((f (list (alist->poly '( ((1 0 0) . 2) ((1 0 1) . 3))) ;2*x+3*x*z
                  (alist->poly '( ((1 0 0) . 5) ((0 1 2) . 7)))) ;5*x+7*y*z^2
           )
         (p (alist->poly '( ((1 1 1) . 11) ((2 3 4) . 13)))) ;11*x*y*z+13*x^2*y^3*z^4
         (sat-ext-1 (list (alist->poly '( ((0 1 0 0) . 2) ((0 1 0 1) . 3)))
                          (alist->poly '( ((0 1 0 0) . 5) ((0 0 1 2) . 7)))
                          (alist->poly '( ((0 0 0 0) . -1) ((1 1 1 1) . 11) ((1 2 3 4) . 13) )))))
    (is (universal-equalp (saturation-extension-1 f p) sat-ext-1))))

(test universal-expt
  "Universal-expt"
  (let ((f (alist->poly '( ((0) . 1) ((1) . 1))))
        (f2 (alist->poly '( ((0) . 1) ((1) . 2) ((2) . 1))))
        (f3 (alist->poly '( ((0) . 1) ((1) . 3) ((2) . 3) ((3) . 1)))))
    (is (universal-equalp (universal-expt f 2) f2))
    (is (universal-equalp (universal-expt f 3) f3))))

(run! 'poly-suite)
(format t "All tests done!~%")