;;; -*- Mode: Lisp -*- ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; ;;; Copyright (C) 1999, 2002, 2009, 2015 Marek Rychlik ;;; ;;; 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. ;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; ;; Operations in ideal theory ;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; (defpackage "IDEAL" (:use :cl :ring :monom :order :term :polynomial :division :grobner-wrap :ring-and-order) (:export "POLY-DEPENDS-P" "RING-INTERSECTION" "ELIMINATION-IDEAL" "COLON-IDEAL" "COLON-IDEAL-1" "IDEAL-INTERSECTION" "POLY-LCM" "GROBNER-EQUAL" "GROBNER-SUBSETP" "GROBNER-MEMBER" "IDEAL-SATURATION-1" "IDEAL-SATURATION" "IDEAL-POLYSATURATION-1" "IDEAL-POLYSATURATION" )) (in-package :ideal) ;; Does the term depend on variable K? (defun term-depends-p (term k) "Return T if the term TERM depends on variable number K." (monom-depends-p (term-monom term) k)) ;; Does the polynomial P depend on variable K? (defun poly-depends-p (p k) "Return T if the term polynomial P depends on variable number K." (some #'(lambda (term) (term-depends-p term k)) (poly-termlist p))) (defun ring-intersection (plist k) "This function assumes that polynomial list PLIST=(P[0],P[1],...,P[J-1]) is a Grobner basis and it calculates the intersection of Id({P[0],P[1],...,P[J-1]}) with the ring R[X[K],...,X[N-1]], i.e. it discards polynomials which depend on variables X[0], X[1], ..., X[K-1]." (dotimes (i k plist) (setf plist (remove-if #'(lambda (p) (poly-depends-p p i)) plist)))) (defun elimination-ideal (ring-and-order flist k &optional (top-reduction-only $poly_top_reduction_only) (start 0)) "Given a list of polynomials FLIST, and an integer K, tt finds and returns the Groebner basis the elimination ideal of Id({FLIST}) obtained by eliminating the first K variables. Optional argument TOP-REDUCTION-ONLY indicates whether to fully reduce or only top-reduce. Optional argument START, defaulting to 0, is used to indicate that the first START elements of F form a Groebner basis." (ring-intersection (reduced-grobner ring-and-order flist start top-reduction-only) k)) (defun colon-ideal (ring-and-order f g &optional (top-reduction-only $poly_top_reduction_only) &aux (ring (ro-ring ring-and-order))) "Returns the reduced Grobner basis of the colon ideal Id(F):Id(G), where F and G are two lists of polynomials. The colon ideal I:J is defined as the set of polynomials H such that for all polynomials W in J the polynomial W*H belongs to I." (declare (type ring-and-order ring-and-order)) (cond ((endp g) ;;Id(G) consists of 0 only so W*0=0 belongs to Id(F) (if (every #'poly-zerop f) (error "First ideal must be non-zero.") (list (make-poly-from-termlist (list (make-term (make-monom :dimension (monom-dimension (poly-lm (find-if-not #'poly-zerop f)))) (funcall (ring-unit ring)))))))) ((endp (cdr g)) (colon-ideal-1 ring-and-order f (car g) top-reduction-only)) (t (ideal-intersection ring-and-order (colon-ideal-1 ring-and-order f (car g) top-reduction-only) (colon-ideal ring-and-order f (rest g) top-reduction-only) top-reduction-only)))) (defun colon-ideal-1 (ring-and-order f g &optional (top-reduction-only $poly_top_reduction_only)) "Returns the reduced Grobner basis of the colon ideal Id(F):Id({G}), where F is a list of polynomials and G is a polynomial." (declare (type ring-and-order ring-and-order)) (mapcar #'(lambda (x) (poly-exact-divide ring-and-order x g)) (ideal-intersection ring-and-order f (list g) top-reduction-only))) (defun ideal-intersection (ring-and-order f g &optional (top-reduction-only $poly_top_reduction_only) (ring (ro-ring ring-and-order))) (declare (type ring-and-order ring-and-order)) (mapcar #'poly-contract (ring-intersection (reduced-grobner ring-and-order (append (mapcar #'(lambda (p) (poly-extend p (make-monom :dimension 1 :initial-exponent 1))) f) (mapcar #'(lambda (p) (poly-append (poly-extend (poly-uminus ring p) (make-monom :dimension 1 :initial-exponent 1)) (poly-extend p))) g)) 0 top-reduction-only) 1))) (defun poly-lcm (ring-and-order f g &aux (ring (ro-ring ring-and-order))) "Return LCM (least common multiple) of two polynomials F and G. The polynomials must be ordered according to monomial order PRED and their coefficients must be compatible with the RING structure defined in the COEFFICIENT-RING package." (cond ((poly-zerop f) f) ((poly-zerop g) g) ((and (endp (cdr (poly-termlist f))) (endp (cdr (poly-termlist g)))) (let ((m (monom-lcm (poly-lm f) (poly-lm g)))) (make-poly-from-termlist (list (make-term m (funcall (ring-lcm ring) (poly-lc f) (poly-lc g))))))) (t (multiple-value-bind (f f-cont) (poly-primitive-part ring f) (multiple-value-bind (g g-cont) (poly-primitive-part ring g) (scalar-times-poly ring (funcall (ring-lcm ring) f-cont g-cont) (poly-primitive-part ring (car (ideal-intersection ring-and-order (list f) (list g) nil))))))))) ;; Do two Grobner bases yield the same ideal? (defun grobner-equal (ring-and-order g1 g2) "Returns T if two lists of polynomials G1 and G2, assumed to be Grobner bases, generate the same ideal, and NIL otherwise." (declare (type ring-and-order ring-and-order)) (and (grobner-subsetp ring-and-order g1 g2) (grobner-subsetp ring-and-order g2 g1))) (defun grobner-subsetp (ring-and-order 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." (declare (type ring-and-order ring-and-order)) (every #'(lambda (p) (grobner-member ring-and-order p g2)) g1)) (defun grobner-member (ring-and-order 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." (declare (type ring-and-order ring-and-order)) (poly-zerop (normal-form ring-and-order p g nil))) ;; Calculate F : p^inf (defun ideal-saturation-1 (ring-and-order f p &optional (start 0) (top-reduction-only $poly_top_reduction_only) &aux (ring (ro-ring ring-and-order))) "Returns the reduced Grobner basis of the saturation of the ideal generated by a polynomial list F in the ideal generated by a single polynomial P. The saturation ideal is defined as the set of polynomials H such for some natural number n (* (EXPT P N) H) is in the ideal spanned by F. Geometrically, over an algebraically closed field, this is the set of polynomials in the ideal generated by F which do not identically vanish on the variety of P." (declare (type ring-and-order ring-and-order)) (mapcar #'poly-contract (ring-intersection (reduced-grobner ring-and-order (saturation-extension-1 ring f p) start top-reduction-only) 1))) ;; Calculate F : p1^inf : p2^inf : ... : ps^inf (defun ideal-polysaturation-1 (ring-and-order f plist &optional (start 0) (top-reduction-only $poly_top_reduction_only)) "Returns the reduced Grobner basis of the ideal obtained by a sequence of successive saturations in the polynomials of the polynomial list PLIST of the ideal generated by the polynomial list F." (cond ((endp plist) (reduced-grobner ring-and-order f start top-reduction-only)) (t (let ((g (ideal-saturation-1 ring-and-order f (car plist) start top-reduction-only))) (ideal-polysaturation-1 ring-and-order g (rest plist) (length g) top-reduction-only))))) (defun ideal-saturation (ring-and-order f g &optional (start 0) (top-reduction-only $poly_top_reduction_only) &aux (k (length g)) (ring (ro-ring ring-and-order))) "Returns the reduced Grobner basis of the saturation of the ideal generated by a polynomial list F in the ideal generated a polynomial list G. The saturation ideal is defined as the set of polynomials H such for some natural number n and some P in the ideal generated by G the polynomial P**N * H is in the ideal spanned by F. Geometrically, over an algebraically closed field, this is the set of polynomials in the ideal generated by F which do not identically vanish on the variety of G." (declare (type ring-and-order ring-and-order)) (mapcar #'(lambda (q) (poly-contract q k)) (ring-intersection (reduced-grobner ring-and-order (polysaturation-extension ring f g) start top-reduction-only) k))) (defun ideal-polysaturation (ring-and-order f ideal-list &optional (start 0) (top-reduction-only $poly_top_reduction_only)) "Returns the reduced Grobner basis of the ideal obtained by a successive applications of IDEAL-SATURATION to F and lists of polynomials in the list IDEAL-LIST." (declare (type ring-and-order ring-and-order)) (cond ((endp ideal-list) f) (t (let ((h (ideal-saturation ring-and-order f (car ideal-list) start top-reduction-only))) (ideal-polysaturation ring-and-order h (rest ideal-list) (length h) top-reduction-only)))))