source: branches/f4grobner/buchberger.lisp@ 4503

;;; -*-  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 "BUCHBERGER"
  (:use :cl :grobner-debug 
        :polynomial :division 
        :criterion :pair-queue :priority-queue
        :ring
        )
  (:export "BUCHBERGER" 
           "GROBNER-MEMBER" 
           "GROBNER-SUBSETP" 
           "GROBNER-EQUAL"
           )
  (:documentation "Buchberger Algorithm Implementation."))

(in-package "BUCHBERGER")


(defun buchberger (f 
                   &optional 
                     (start 0)
                     (top-reduction-only $poly_top_reduction_only))
  "An implementation of the Buchberger algorithm. Return Grobner basis
of the ideal generated by the polynomial list F.  Polynomials 0 to
START-1 are assumed to be a Grobner basis already, so that certain
critical pairs will not be examined. If TOP-REDUCTION-ONLY set, top
reduction will be preformed. This function assumes that all polynomials
in F are non-zero."
  (declare (type fixnum start))
  (when (endp f) (return-from buchberger f)) ;cut startup costs
  (debug-cgb "~&GROBNER BASIS - BUCHBERGER ALGORITHM")
  (when (plusp start) (debug-cgb "~&INCREMENTAL:~d done" start))
  #+grobner-check (when (plusp start) 
                    (grobner-test (subseq f 0 start) (subseq f 0 start)))
  ;;Initialize critical pairs
  (let ((b (make-critical-pair-queue *normal-strategy* f start))
        (b-done (make-hash-table :test #'equal)))
    (declare (type critical-pair-queue b) (type hash-table b-done))
    (dotimes (i (1- start))
      (do ((j (1+ i) (1+ j))) ((>= j start))
        (setf (gethash (list (elt f i) (elt f j)) b-done) t)))
    (do ()
        ((queue-empty-p b)
         #+grobner-check(grobner-test f f)
         (debug-cgb "~&GROBNER END")
         f)
      (let ((pair (dequeue b)))
        (declare (type critical-pair pair))
        (cond
          ((criterion-1 pair) nil)
          ((criterion-2 pair b-done f) nil)
          (t 
           (let ((sp (normal-form (s-polynomial
                                   (critical-pair-first pair) 
                                   (critical-pair-second pair))
                                  f top-reduction-only)))
             (declare (type poly sp))
             (cond
               ((universal-zerop sp)
                nil)
               (t
                (setf sp (poly-primitive-part sp)
                      f (nconc f (list sp)))
                ;; Add new critical pairs
                (dolist (h f)
                  (enqueue b (make-instance 'critical-pair :first h :second sp)))
                (debug-cgb "~&Polynomials: ~d; Pairs left: ~d; Pairs done: ~d;"
                           (length f) (queue-size b)
                           (hash-table-count b-done)))))))
        (setf (gethash (list (critical-pair-first pair) (critical-pair-second pair)) b-done)
              t)))))


(defun grobner-member (p g)
  "Returns T if a polynomial P belongs to the ideal generated by the
polynomial list G, which is assumed to be a Grobner basis. Returns NIL otherwise."
  (universal-zerop (normal-form p g nil)))


(defun grobner-subsetp (g1 g2)
  "Returns T if a list of polynomials G1 generates
an ideal contained in the ideal generated by a polynomial list G2,
both G1 and G2 assumed to be Grobner bases. Returns NIL otherwise."
  (every #'(lambda (p) (grobner-member p g2)) g1))


(defun grobner-equal (g1 g2)
  "Returns T if two lists of polynomials G1 and G2, assumed to be Grobner bases,
generate  the same ideal, and NIL otherwise."
  (and (grobner-subsetp g1 g2) (grobner-subsetp g2 g1)))