diff options
Diffstat (limited to 'src/algebra/bezout.spad.pamphlet')
| -rw-r--r-- | src/algebra/bezout.spad.pamphlet | 206 |
1 files changed, 206 insertions, 0 deletions
diff --git a/src/algebra/bezout.spad.pamphlet b/src/algebra/bezout.spad.pamphlet new file mode 100644 index 00000000..7f6b5400 --- /dev/null +++ b/src/algebra/bezout.spad.pamphlet @@ -0,0 +1,206 @@ +\documentclass{article} +\usepackage{axiom} +\begin{document} +\title{\$SPAD/src/algebra bezout.spad} +\author{Clifton J. Williamson} +\maketitle +\begin{abstract} +\end{abstract} +\eject +\tableofcontents +\eject +\section{package BEZOUT BezoutMatrix} +<<package BEZOUT BezoutMatrix>>= +)abbrev package BEZOUT BezoutMatrix +++ Author: Clifton J. Williamson +++ Date Created: 2 August 1988 +++ Date Last Updated: 3 November 1993 +++ Basic Operations: bezoutMatrix, bezoutResultant, bezoutDiscriminant +++ Related Domains +++ Also See: +++ AMS Classifiactions: +++ Keywords: Bezout matrix, resultant, discriminant +++ Examples: +++ Reference: Knuth, The Art of Computer Programming, 2nd edition, +++ Vol. 2, p. 619, problem 12. +++ Description: +++ \spadtype{BezoutMatrix} contains functions for computing resultants and +++ discriminants using Bezout matrices. + +BezoutMatrix(R,UP,M,Row,Col): Exports == Implementation where + R : Ring + UP : UnivariatePolynomialCategory R + Row : FiniteLinearAggregate R + Col : FiniteLinearAggregate R + M : MatrixCategory(R,Row,Col) + I ==> Integer + lc ==> leadingCoefficient + + Exports ==> with + sylvesterMatrix: (UP,UP) -> M + ++ sylvesterMatrix(p,q) returns the Sylvester matrix for the two + ++ polynomials p and q. + bezoutMatrix: (UP,UP) -> M + ++ bezoutMatrix(p,q) returns the Bezout matrix for the two + ++ polynomials p and q. + + if R has commutative("*") then + bezoutResultant: (UP,UP) -> R + ++ bezoutResultant(p,q) computes the resultant of the two + ++ polynomials p and q by computing the determinant of a Bezout matrix. + + bezoutDiscriminant: UP -> R + ++ bezoutDiscriminant(p) computes the discriminant of a polynomial p + ++ by computing the determinant of a Bezout matrix. + + Implementation ==> add + + sylvesterMatrix(p,q) == + n1 := degree p; n2 := degree q; n := n1 + n2 + sylmat : M := new(n,n,0) + minR := minRowIndex sylmat; minC := minColIndex sylmat + maxR := maxRowIndex sylmat; maxC := maxColIndex sylmat + p0 := p + -- fill in coefficients of 'p' + while not zero? p0 repeat + coef := lc p0; deg := degree p0; p0 := reductum p0 + -- put bk = coef(p,k) in sylmat(minR + i,minC + i + (n1 - k)) + for i in 0..n2 - 1 repeat + qsetelt_!(sylmat,minR + i,minC + n1 - deg + i,coef) + q0 := q + -- fill in coefficients of 'q' + while not zero? q0 repeat + coef := lc q0; deg := degree q0; q0 := reductum q0 + for i in 0..n1-1 repeat + qsetelt_!(sylmat,minR + n2 + i,minC + n2 - deg + i,coef) + sylmat + + bezoutMatrix(p,q) == + -- This function computes the Bezout matrix for 'p' and 'q'. + -- See Knuth, The Art of Computer Programming, Vol. 2, p. 619, # 12. + -- One must have deg(p) >= deg(q), so the arguments are reversed + -- if this is not the case. + n1 := degree p; n2 := degree q; n := n1 + n2 + n1 < n2 => bezoutMatrix(q,p) + m1 : I := n1 - 1; m2 : I := n2 - 1; m : I := n - 1 + -- 'sylmat' will be a matrix consisting of the first n1 columns + -- of the standard Sylvester matrix for 'p' and 'q' + sylmat : M := new(n,n1,0) + minR := minRowIndex sylmat; minC := minColIndex sylmat + maxR := maxRowIndex sylmat; maxC := maxColIndex sylmat + p0 := p + -- fill in coefficients of 'p' + while not ground? p0 repeat + coef := lc p0; deg := degree p0; p0 := reductum p0 + -- put bk = coef(p,k) in sylmat(minR + i,minC + i + (n1 - k)) + -- for i = 0... + -- quit when i > m2 or when i + (n1 - k) > m1, whichever happens first + for i in 0..min(m2,deg - 1) repeat + qsetelt_!(sylmat,minR + i,minC + n1 - deg + i,coef) + q0 := q + -- fill in coefficients of 'q' + while not zero? q0 repeat + coef := lc q0; deg := degree q0; q0 := reductum q0 + -- put ak = coef(q,k) in sylmat(minR + n1 + i,minC + i + (n2 - k)) + -- for i = 0... + -- quit when i > m1 or when i + (n2 - k) > m1, whichever happens first + -- since n2 - k >= 0, we quit when i + (n2 - k) > m1 + for i in 0..(deg + n1 - n2 - 1) repeat + qsetelt_!(sylmat,minR + n2 + i,minC + n2 - deg + i,coef) + -- 'bezmat' will be the 'Bezout matrix' as described in Knuth + bezmat : M := new(n1,n1,0) + for i in 0..m2 repeat + -- replace A_i by (b_0 A_i + ... + b_{n_2-1-i} A_{n_2 - 1}) - + -- (a_0 B_i + ... + a_{n_2-1-i} B_{n_2-1}), as in Knuth + bound : I := n2 - i; q0 := q + while not zero? q0 repeat + deg := degree q0 + if (deg < bound) then + -- add b_deg A_{n_2 - deg} to the new A_i + coef := lc q0 + for k in minC..maxC repeat + c := coef * qelt(sylmat,minR + m2 - i - deg,k) + + qelt(bezmat,minR + m2 - i,k) + qsetelt_!(bezmat,minR + m2 - i,k,c) + q0 := reductum q0 + p0 := p + while not zero? p0 repeat + deg := degree p0 + if deg < bound then + coef := lc p0 + -- subtract a_deg B_{n_2 - deg} from the new A_i + for k in minC..maxC repeat + c := -coef * qelt(sylmat,minR + m - i - deg,k) + + qelt(bezmat,minR + m2 - i,k) + qsetelt_!(bezmat,minR + m2 - i,k,c) + p0 := reductum p0 + for i in n2..m1 repeat for k in minC..maxC repeat + qsetelt_!(bezmat,minR + i,k,qelt(sylmat,minR + i,k)) + bezmat + + if R has commutative("*") then + + bezoutResultant(f,g) == determinant bezoutMatrix(f,g) + + if R has IntegralDomain then + + bezoutDiscriminant f == + degMod4 := (degree f) rem 4 + (degMod4 = 0) or (degMod4 = 1) => + (bezoutResultant(f,differentiate f) exquo (lc f)) :: R + -((bezoutResultant(f,differentiate f) exquo (lc f)) :: R) + + else + + bezoutDiscriminant f == + lc f = 1 => + degMod4 := (degree f) rem 4 + (degMod4 = 0) or (degMod4 = 1) => + bezoutResultant(f,differentiate f) + -bezoutResultant(f,differentiate f) + error "bezoutDiscriminant: leading coefficient must be 1" + +@ +\section{License} +<<license>>= +--Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd. +--All rights reserved. +-- +--Redistribution and use in source and binary forms, with or without +--modification, are permitted provided that the following conditions are +--met: +-- +-- - Redistributions of source code must retain the above copyright +-- notice, this list of conditions and the following disclaimer. +-- +-- - Redistributions in binary form must reproduce the above copyright +-- notice, this list of conditions and the following disclaimer in +-- the documentation and/or other materials provided with the +-- distribution. +-- +-- - Neither the name of The Numerical ALgorithms Group Ltd. nor the +-- names of its contributors may be used to endorse or promote products +-- derived from this software without specific prior written permission. +-- +--THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS +--IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +--TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A +--PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER +--OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, +--EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, +--PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR +--PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF +--LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING +--NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS +--SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +@ +<<*>>= +<<license>> + +<<package BEZOUT BezoutMatrix>> +@ +\eject +\begin{thebibliography}{99} +\bibitem{1} nothing +\end{thebibliography} +\end{document} |