source: branches/f4grobner/symbolic-polynomial.lisp@ 4353

;;; -*-  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.  
;;;                                                                              
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

(defpackage "SYMBOLIC-POLYNOMIAL"
  (:use :cl :utils :monom :polynomial :infix :infix-printer :ring)
  (:export "SYMBOLIC-POLY" "READ-INFIX-FORM" "STRING->POLY" "POLY->STRING" "->INFIX")
  (:documentation "Implements symbolic polynomials. A symbolic
polynomial is polynomial which uses symbolic variables for reading and
printing in standard human-readable (infix) form."))

(in-package :symbolic-polynomial)

(defclass symbolic-poly (poly)
  ((vars :initform nil 
         :initarg :vars 
         :accessor symbolic-poly-vars)
   )
  (:default-initargs :termlist nil :vars nil))

(defmethod print-object ((self symbolic-poly) stream)
  (print-unreadable-object (self stream :type t :identity t)
    (with-accessors ((dimension poly-dimension)
                     (termlist poly-termlist) 
                     (order poly-term-order)
                     (vars symbolic-poly-vars))
        self
      (format stream "DIMENSION=~A TERMLIST=~A ORDER=~A VARS=~A" 
              dimension termlist order vars))))


(defmethod universal-equalp ((self symbolic-poly) (other symbolic-poly))
  (when (universal-equalp (symbolic-poly-vars self) (symbolic-poly-vars other))
    (call-next-method)))

(defmethod universal-equalp ((self symbol) (other symbol))
  (eq self other))

(defmethod update-instance-for-different-class :after ((old poly) (new  symbolic-poly) &key)
  "After adding variables to NEW, we need to make sure that the number
of variables given by POLY-DIMENSION is consistent with VARS."
  (assert (= (length (symbolic-poly-vars new)) (poly-dimension new))))


#|
(defun poly-eval-scalar (expr 
                         &aux 
                           (order #'lex>))
  "Evaluate a scalar expression EXPR in ring RING."
  (declare (type ring ring))
  (poly-lc (poly-eval expr nil ring order)))
|#


(defun read-infix-form (&key (stream t))
  "Parser of infix expressions with integer/rational coefficients
The parser will recognize two kinds of polynomial expressions:

- polynomials in fully expanded forms with coefficients
  written in front of symbolic expressions; constants can be optionally
  enclosed in (); for example, the infix form
        X^2-Y^2+(-4/3)*U^2*W^3-5
  parses to
        (+ (- (EXPT X 2) (EXPT Y 2)) (* (- (/ 4 3)) (EXPT U 2) (EXPT W 3)) (- 5))

- lists of polynomials; for example
        [X-Y, X^2+3*Z]          
  parses to
          (:[ (- X Y) (+ (EXPT X 2) (* 3 Z)))
  where the first symbol [ marks a list of polynomials.

-other infix expressions, for example
        [(X-Y)*(X+Y)/Z,(X+1)^2]
parses to:
        (:[ (/ (* (- X Y) (+ X Y)) Z) (EXPT (+ X 1) 2))
Currently this function is implemented using M. Kantrowitz's INFIX package."
  (read-from-string
   (concatenate 'string
                "#I(" 
                (with-output-to-string (s)
                  (loop
                     (multiple-value-bind (line eof)
                         (read-line stream t)
                       (format s "~A" line)
                       (when eof (return)))))
                ")")))

(defun read-poly (vars &key
                         (stream t) 
                         (order #'lex>)
                         (coefficient-class 'rational-field))
  "Reads an expression in prefix form from a stream STREAM.
The expression read from the strem should represent a polynomial or a
list of polynomials in variables VARS, over the ring RING.  The
polynomial or list of polynomials is returned, with terms in each
polynomial ordered according to monomial order ORDER."
  (poly-eval (read-infix-form :stream stream) vars order coefficient-class))

(defun string->poly (str vars 
                     &optional
                       (order #'lex>)
                       (coefficient-class 'rational-field))
  "Converts a string STR to a polynomial in variables VARS."
  (with-input-from-string (s str)
    (let ((p-or-plist (read-poly vars :stream s :order order :coefficient-class coefficient-class)))
      (etypecase p-or-plist
        (poly (change-class p-or-plist 'symbolic-poly :vars vars))
        (cons 
         (setf (cdr p-or-plist) (mapcar #'(lambda (p) (change-class p 'symbolic-poly :vars vars)) (cdr p-or-plist)))
         p-or-plist)))))

(defun poly->alist (p)
  "Convert a polynomial P to an association list. Thus, the format of the
returned value is  ((MONOM[0] . COEFF[0]) (MONOM[1] . COEFF[1]) ...), where
MONOM[I] is a list of exponents in the monomial and COEFF[I] is the
corresponding coefficient in the ring."
  (cond
    ((poly-p p)
     (mapcar #'->list (poly-termlist p)))
    ((and (consp p) (eq (car p) :[))
     (cons :[ (mapcar #'poly->alist (cdr p))))))

(defun string->alist (str vars
                      &optional
                        (order #'lex>)
                        (coefficient-class 'rational-field))
  "Convert a string STR representing a polynomial or polynomial list to
an association list (... (MONOM . COEFF) ...)."
  (poly->alist (string->poly str vars order coefficient-class)))

(defun poly-equal-no-sugar-p (p q)
  "Compare polynomials for equality, ignoring sugar."
  (declare (type poly p q))
  (equalp (poly-termlist p) (poly-termlist q)))

(defun poly-set-equal-no-sugar-p (p q)
  "Compare polynomial sets P and Q for equality, ignoring sugar."
  (null (set-exclusive-or  p q :test #'poly-equal-no-sugar-p )))

(defun poly-list-equal-no-sugar-p (p q)
  "Compare polynomial lists P and Q for equality, ignoring sugar."
  (every #'poly-equal-no-sugar-p p q))

(defmethod ->sexp :around ((self symbolic-poly) &optional (vars (symbolic-poly-vars self)))
  "Convert a symbolic polynomial SELF to infix form, using variables VARS. The default
value of VARS is the corresponding slot value of SELF."
  (call-next-method self vars))

(defgeneric poly->string (self &optional vars)
  (:documentation "Render polynomial SELF as a string, using symbolic variables VARS.")
  (:method ((self list) &optional (vars nil vars-p))
    (assert (eql (car self) :[))
    (cond (vars-p
           (format nil "[~{~a~^, ~}]" (mapcar #'(lambda (p) (poly->string p vars)) (cdr self))))
          (t 
           (format nil "[~{~a~^, ~}]" (mapcar #'(lambda (p) (poly->string p)) (cdr self))))))
  (:method ((self poly) &optional (vars nil))
    ;; Ensure that the number of variables matches the dimension
    (assert (= (length vars) (poly-dimension self)))
    (infix-print-to-string (->sexp self vars)))
  (:method ((self symbolic-poly) &optional (vars (symbolic-poly-vars self)))
    (infix-print-to-string (->sexp self vars))))