;;; -*-  mode: lisp; package: maxima; syntax: common-lisp; base: 10 -*- 
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;                                                                              
;;;  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.	 
;;;										 
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

(defpackage "RING"
  (:use :cl)
  (:export "R-PARSE"
	   "UNIT-ELEMENT"
	   "R-ZEROP"
	   "R+"
	   "R-"
	   "R*"
	   "R+"
	   "R/"
	   "R-EXPT"
	   "R-LCM"
	   "R-EZGCD"
	   "R-GCD"
	   "R-TOTAL-DEGREE"
	   "R-DIMENSION"
	   "R-EXPONENTS"
	   "R-COEFF"
	   "R-SUGAR"
	   "R-DIVIDES-P"   
	   "R-DIVIDES-LCM-P"
	   "R-LCM-DIVIDES-LCM-P"
	   "R-LCM-EQUAL-LCM-P"
	   "R-REL-PRIME-P"    
	   "R-EQUALP"
	   "R-CLONE"
	   "R-ELT"
	   "R->LIST"
	   "R-DIVISIBLE-BY-P"	
	   "R-REL-PRIME-P"	
	   "R-DEPENDS-P"	
	   "LEFT-TENSOR-PRODUCT-BY"
	   "RIGHT-TENSOR-PRODUCT-BY"
	   "LEFT-CONTRACT"
	   "R-LENGTH"
	   "MULTIPLY-BY"
	   "DIVIDE-BY"
	   "ADD-TO"
	   "SUBTRACT-FROM"
	   "UNARY-MINUS"
	   "SCALAR"
	   "SCALAR-COEFF"
	   "INSERT-ITEM"
	   "APPEND-ITEM")
  (:shadowing-import-from 
   #+openmcl-native-threads #:ccl 
   #+cmu #:pcl
   #+sbcl #:sb-pcl 
   #+lispworks #:hcl 
   #+allegro #:mop 
   #+clisp #:clos
   #:class-slots #:slot-definition-name)
  (:documentation
   "Implements ring operations.  These are all operations that are
performed on the coefficients by the package, and thus the coefficient
ring can be changed by merely redefining these operations."))

(in-package :ring)

(defclass scalar () 
  ((coeff :initarg :coeff :accessor scalar-coeff))
  (:default-initargs :coeff nil)
  (:documentation "Wraps objects suitable as scalars/polynomial coefficients"))

(defmethod print-object ((self scalar) stream)
  (print-unreadable-object (self stream :type t :identity t)
    (format stream "COEFF=~A"
	    (slot-value self 'coeff))))

(defgeneric unit-element (class))

(defgeneric r-zerop (object)
  (:method ((self number)) (zerop self))
  (:documentation "Tests whether a ring element is 0."))

(defgeneric r+ (x y)
  (:method ((x number) (y number)) (+ x y))
  (:documentation "Adds ring elements.")) 

(defgeneric r- (minuend &rest subtrahends)
  (:method ((minuend number) &rest subtrahends) (apply #'- (cons minuend subtrahends)))
  (:documentation "Subtracts ring elements.")) 

(defgeneric r* (x y)
  (:method (x y) (* x y))
  (:documentation "Multiplies ring elements.")) 

(defgeneric left-tensor-product-by (self other)
  (:documentation "Takes a tensor product of SELF with OTHER, where
OTHER is the left factor."))

(defgeneric right-tensor-product-by (self other)
  (:documentation "Takes a tensor product of SELF with OTHER, where
OTHER is the right factor."))

(defgeneric r/ (numerator &rest denominators)
  (:method ((numerator number) &rest denominators) (apply #'/ (cons numerator denominators)))
  (:documentation "Divides ring elements.")) 

(defgeneric r-lcm (x y)
  (:method ((x integer) (y integer)) (lcm x y))
  (:documentation "Returns the least common multiple of ring elements.")) 

(defgeneric r-expt (x y)
  (:method ((x number) (y integer)) (expt x y))
  (:method ((x t) (y integer)) 
    (declare (type fixnum y))
    (cond
      ((minusp y) (error "r-expt: Negative exponent."))
      ((r-zerop x) (if (zerop y) 1))
      (t
       (do ((k 1 (ash k 1))
	    (q x (r* q q))		;keep squaring
	    (p 1 (if (not (zerop (logand k y))) (r* p q) p)))
	   ((> k y) p)
	 (declare (fixnum k))))))
  (:documentation "Raises X to power Y.")) 

(defgeneric r-ezgcd (x y)
  (:method ((x integer) (y integer)
	    &aux (c (gcd x y)))
    (values c (/ x c) (/ y c)))
  (:documentation "Solves the diophantine system: X=C*X1, Y=C*X2,
C=GCD(X,Y).  It returns C, X1 and Y1. The result may be obtained by
the Euclidean algorithm."))

(defgeneric r-gcd (x y)
  (:method ((x integer) (y integer))
    (gcd x y))
  (:documentation "Returns GCD(X,Y)."))

(defgeneric r-dimension (object))
(defgeneric r-exponents (object))

(defgeneric r-coeff (object))
(defgeneric (setf r-coeff) (new-value object))

(defgeneric r-total-degree (object &optional start end))

(defgeneric r-divides-p (object1 object2)
  (:method ((object1 integer) (object2 integer))
    (zerop (rem object2 object1)))
  (:documentation "Returns T if OBJECT1 divides OBJECT2"))

(defgeneric r-divides-lcm-p (object1 object2 object3)
  (:documentation "Returns T if OBJECT divides LCM(OBJECT2, OBJECT3), NIL otherwise."))

(defgeneric r-lcm-divides-lcm-p (object1 object2 object3 object4)
  (:documentation "Returns T if LCM(OBJECT1,OBJECT2) divides LCM(OBJECT3,OBJECT4), NIL otherwise."))

(defgeneric r-lcm-equal-lcm-p (object1 object2 object3 object4)
  (:documentation "Returns T if object LCM(OBJECT1,OBJECT2) equals LCM(OBJECT3,OBJECT4), NIL otherwise."))

(defgeneric r-equalp (object1 object2)
  (:method (object1 object2) (equalp object1 object2))
  (:method ((object1 list) (object2 list))
    (every #'r-equalp object1 object2))
  (:method ((object1 scalar) (object2 scalar))
    (r-equalp (scalar-coeff object1) (scalar-coeff object2)))
  (:documentation "Equality using deep comparison of object slots."))

(defgeneric r-elt (object index)
  (:documentation "Access a part of an object OBJECT with index INDEX."))

(defgeneric (setf r-elt) (new-value object index)
  (:documentation "A setter of a part of an object OBJECT with index INDEX."))

(defgeneric r-length (object))

(defgeneric r->list (object))
(defgeneric r-sugar (object))
(defgeneric r-rel-prime-p (object1 object2))
(defgeneric left-contract (object k))
(defgeneric r-divisible-by-p (object1 object2))
(defgeneric r-depends-p (object k))

(defgeneric multiply-by (self other)
  (:method (self other) (r* self other))
  (:documentation "Multiply object SELF and OTHER and store the result
into SELF.  It returns SELF. For instances of a class, this operation
may be destructive."))

(defgeneric divide-by (self other)
  (:method (self other) (r/ self other))
  (:documentation "Divided object SELF by OTHER and store the result
into SELF.  It returns SELF. For instances of a class, this operation
may be destructive."))

(defgeneric add-to (self other)
  (:documentation "Add to object SELF another object OTHER.  For
complex objects, it may destructively modify SELF and destructively
modify/invalidate object OTHER. For standard classes implementing this
method, the result should be an object which is EQ to SELF. For
built-in classes, such as NUMBER, the returned object may not be EQ to
the original, but it will be EQL to it.")
  (:method (self other) (r+ self other)))

(defgeneric subtract-from (self other)
  (:documentation "Subtract from an object SELF another object OTHER.
For complex objects, it may destructively modify SELF and
destructively modify/invalidate object OTHER. For standard classes
implementing this method, the result should be an object which is EQ
to SELF. For built-in classes, such as NUMBER, the returned object may
not be EQ to the original.")
  (:method (self other) (r- self other)))

(defgeneric unary-minus (self)
  (:method ((x number)) (- x)))

(defgeneric insert-item (self item))
(defgeneric append-item (self item))