[1201] | 1 | ;;; -*- Mode: Lisp -*-
|
---|
[147] | 2 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
|
---|
| 3 | ;;;
|
---|
| 4 | ;;; Copyright (C) 1999, 2002, 2009, 2015 Marek Rychlik <rychlik@u.arizona.edu>
|
---|
| 5 | ;;;
|
---|
| 6 | ;;; This program is free software; you can redistribute it and/or modify
|
---|
| 7 | ;;; it under the terms of the GNU General Public License as published by
|
---|
| 8 | ;;; the Free Software Foundation; either version 2 of the License, or
|
---|
| 9 | ;;; (at your option) any later version.
|
---|
| 10 | ;;;
|
---|
| 11 | ;;; This program is distributed in the hope that it will be useful,
|
---|
| 12 | ;;; but WITHOUT ANY WARRANTY; without even the implied warranty of
|
---|
| 13 | ;;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
---|
| 14 | ;;; GNU General Public License for more details.
|
---|
| 15 | ;;;
|
---|
| 16 | ;;; You should have received a copy of the GNU General Public License
|
---|
| 17 | ;;; along with this program; if not, write to the Free Software
|
---|
| 18 | ;;; Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
|
---|
| 19 | ;;;
|
---|
| 20 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
|
---|
| 21 |
|
---|
[63] | 22 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
|
---|
| 23 | ;;
|
---|
| 24 | ;; An implementation of the algorithm of Gebauer and Moeller, as
|
---|
| 25 | ;; described in the book of Becker-Weispfenning, p. 232
|
---|
| 26 | ;;
|
---|
| 27 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
|
---|
| 28 |
|
---|
[494] | 29 | (defpackage "GEBAUER-MOELLER"
|
---|
[581] | 30 | (:use :cl :grobner-debug
|
---|
[1614] | 31 | :division :ring :monom :polynomial :order :ring-and-order
|
---|
[581] | 32 | :pair-queue :priority-queue
|
---|
| 33 | )
|
---|
[494] | 34 | (:export "GEBAUER-MOELLER"))
|
---|
| 35 |
|
---|
| 36 | (in-package :gebauer-moeller)
|
---|
| 37 |
|
---|
[1314] | 38 | (defun gebauer-moeller (ring-and-order f
|
---|
| 39 | &optional
|
---|
| 40 | (start 0)
|
---|
| 41 | (top-reduction-only $poly_top_reduction_only)
|
---|
| 42 | &aux
|
---|
[1317] | 43 | (ring (ro-ring ring-and-order)))
|
---|
[63] | 44 | "Compute Grobner basis by using the algorithm of Gebauer and
|
---|
| 45 | Moeller. This algorithm is described as BUCHBERGERNEW2 in the book by
|
---|
| 46 | Becker-Weispfenning entitled ``Grobner Bases''. This function assumes
|
---|
| 47 | that all polynomials in F are non-zero."
|
---|
| 48 | (declare (ignore top-reduction-only)
|
---|
[1313] | 49 | (type fixnum start)
|
---|
| 50 | (type ring-and-order ring-and-order))
|
---|
[63] | 51 | (cond
|
---|
| 52 | ((endp f) (return-from gebauer-moeller nil))
|
---|
| 53 | ((endp (cdr f))
|
---|
| 54 | (return-from gebauer-moeller (list (poly-primitive-part ring (car f))))))
|
---|
| 55 | (debug-cgb "~&GROBNER BASIS - GEBAUER MOELLER ALGORITHM")
|
---|
| 56 | (when (plusp start) (debug-cgb "~&INCREMENTAL:~d done" start))
|
---|
| 57 | #+grobner-check (when (plusp start)
|
---|
[1316] | 58 | (grobner-test ring-and-order (subseq f 0 start) (subseq f 0 start)))
|
---|
[63] | 59 | (let ((b (make-pair-queue))
|
---|
| 60 | (g (subseq f 0 start))
|
---|
| 61 | (f1 (subseq f start)))
|
---|
| 62 | (do () ((endp f1))
|
---|
| 63 | (multiple-value-setq (g b)
|
---|
| 64 | (gebauer-moeller-update g b (poly-primitive-part ring (pop f1)))))
|
---|
| 65 | (do () ((pair-queue-empty-p b))
|
---|
| 66 | (let* ((pair (pair-queue-remove b))
|
---|
| 67 | (g1 (pair-first pair))
|
---|
| 68 | (g2 (pair-second pair))
|
---|
[1320] | 69 | (h (normal-form ring-and-order (spoly ring-and-order g1 g2)
|
---|
[63] | 70 | g
|
---|
| 71 | nil #| Always fully reduce! |#
|
---|
| 72 | )))
|
---|
| 73 | (unless (poly-zerop h)
|
---|
| 74 | (setf h (poly-primitive-part ring h))
|
---|
| 75 | (multiple-value-setq (g b)
|
---|
| 76 | (gebauer-moeller-update g b h))
|
---|
| 77 | (debug-cgb "~&Sugar: ~d Polynomials: ~d; Pairs left: ~d~%"
|
---|
| 78 | (pair-sugar pair) (length g) (pair-queue-size b))
|
---|
| 79 | )))
|
---|
[1316] | 80 | #+grobner-check(grobner-test ring-and-order g f)
|
---|
[63] | 81 | (debug-cgb "~&GROBNER END")
|
---|
| 82 | g))
|
---|
| 83 |
|
---|
| 84 | (defun gebauer-moeller-update (g b h
|
---|
| 85 | &aux
|
---|
| 86 | c d e
|
---|
| 87 | (b-new (make-pair-queue))
|
---|
| 88 | g-new)
|
---|
| 89 | "An implementation of the auxillary UPDATE algorithm used by the
|
---|
| 90 | Gebauer-Moeller algorithm. G is a list of polynomials, B is a list of
|
---|
| 91 | critical pairs and H is a new polynomial which possibly will be added
|
---|
| 92 | to G. The naming conventions used are very close to the one used in
|
---|
| 93 | the book of Becker-Weispfenning."
|
---|
| 94 | (declare
|
---|
| 95 | #+allegro (dynamic-extent b)
|
---|
| 96 | (type poly h)
|
---|
| 97 | (type priority-queue b))
|
---|
| 98 | (setf c g d nil)
|
---|
| 99 | (do () ((endp c))
|
---|
| 100 | (let ((g1 (pop c)))
|
---|
| 101 | (declare (type poly g1))
|
---|
| 102 | (when (or (monom-rel-prime-p (poly-lm h) (poly-lm g1))
|
---|
| 103 | (and
|
---|
| 104 | (notany #'(lambda (g2) (monom-lcm-divides-monom-lcm-p
|
---|
| 105 | (poly-lm h) (poly-lm g2)
|
---|
| 106 | (poly-lm h) (poly-lm g1)))
|
---|
| 107 | c)
|
---|
| 108 | (notany #'(lambda (g2) (monom-lcm-divides-monom-lcm-p
|
---|
| 109 | (poly-lm h) (poly-lm g2)
|
---|
| 110 | (poly-lm h) (poly-lm g1)))
|
---|
| 111 | d)))
|
---|
| 112 | (push g1 d))))
|
---|
| 113 | (setf e nil)
|
---|
| 114 | (do () ((endp d))
|
---|
| 115 | (let ((g1 (pop d)))
|
---|
| 116 | (declare (type poly g1))
|
---|
| 117 | (unless (monom-rel-prime-p (poly-lm h) (poly-lm g1))
|
---|
| 118 | (push g1 e))))
|
---|
| 119 | (do () ((pair-queue-empty-p b))
|
---|
| 120 | (let* ((pair (pair-queue-remove b))
|
---|
| 121 | (g1 (pair-first pair))
|
---|
| 122 | (g2 (pair-second pair)))
|
---|
| 123 | (declare (type pair pair)
|
---|
| 124 | (type poly g1 g2))
|
---|
| 125 | (when (or (not (monom-divides-monom-lcm-p
|
---|
| 126 | (poly-lm h)
|
---|
| 127 | (poly-lm g1) (poly-lm g2)))
|
---|
| 128 | (monom-lcm-equal-monom-lcm-p
|
---|
| 129 | (poly-lm g1) (poly-lm h)
|
---|
| 130 | (poly-lm g1) (poly-lm g2))
|
---|
| 131 | (monom-lcm-equal-monom-lcm-p
|
---|
| 132 | (poly-lm h) (poly-lm g2)
|
---|
| 133 | (poly-lm g1) (poly-lm g2)))
|
---|
| 134 | (pair-queue-insert b-new (make-pair g1 g2)))))
|
---|
| 135 | (dolist (g3 e)
|
---|
| 136 | (pair-queue-insert b-new (make-pair h g3)))
|
---|
| 137 | (setf g-new nil)
|
---|
| 138 | (do () ((endp g))
|
---|
| 139 | (let ((g1 (pop g)))
|
---|
| 140 | (declare (type poly g1))
|
---|
| 141 | (unless (monom-divides-p (poly-lm h) (poly-lm g1))
|
---|
| 142 | (push g1 g-new))))
|
---|
| 143 | (push h g-new)
|
---|
| 144 | (values g-new b-new))
|
---|