/usr/share/axiom-20170501/src/algebra/UPSQFREE.spad is in axiom-source 20170501-3.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 | )abbrev package UPSQFREE UnivariatePolynomialSquareFree
++ Author: Dave Barton, Barry Trager
++ Description:
++ This package provides for square-free decomposition of
++ univariate polynomials over arbitrary rings,
++ a partial factorization such that each factor is a product
++ of irreducibles with multiplicity one and the factors are
++ pairwise relatively prime. If the ring
++ has characteristic zero, the result is guaranteed to satisfy
++ this condition. If the ring is an infinite ring of
++ finite characteristic, then it may not be possible to decide when
++ polynomials contain factors which are pth powers. In this
++ case, the flag associated with that polynomial is set to "nil"
++ (meaning that that polynomials are not guaranteed to be square-free).
UnivariatePolynomialSquareFree(RC,P) : SIG == CODE where
RC : IntegralDomain
P : Join(UnivariatePolynomialCategory(RC),IntegralDomain) with
gcd : (%,%) -> %
++ gcd(p,q) computes the greatest-common-divisor of p and q.
fUnion ==> Union("nil", "sqfr", "irred", "prime")
FF ==> Record(flg:fUnion, fctr:P, xpnt:Integer)
SIG ==> with
squareFree : P -> Factored(P)
++ squareFree(p) computes the square-free factorization of the
++ univariate polynomial p. Each factor has no repeated roots, and the
++ factors are pairwise relatively prime.
squareFreePart : P -> P
++ squareFreePart(p) returns a polynomial which has the same
++ irreducible factors as the univariate polynomial p, but each
++ factor has multiplicity one.
BumInSepFFE : FF -> FF
++ BumInSepFFE(f) is a local function, exported only because
++ it has multiple conditional definitions.
CODE ==> add
if RC has CharacteristicZero then
squareFreePart(p:P) == (p exquo gcd(p, differentiate p))::P
else
squareFreePart(p:P) ==
unit(s := squareFree(p)$%) * */[f.factor for f in factors s]
if RC has FiniteFieldCategory then
BumInSepFFE(ffe:FF) ==
["sqfr", map(charthRoot,ffe.fctr), characteristic$P*ffe.xpnt]
else if RC has CharacteristicNonZero then
BumInSepFFE(ffe:FF) ==
np:=multiplyExponents(ffe.fctr,characteristic$P:NonNegativeInteger)
(nthrp := charthRoot(np)) case "failed" =>
["nil", np, ffe.xpnt]
["sqfr", nthrp, characteristic$P*ffe.xpnt]
else
BumInSepFFE(ffe:FF) ==
["nil",
multiplyExponents(ffe.fctr,characteristic$P:NonNegativeInteger),
ffe.xpnt]
if RC has CharacteristicZero then
squareFree(p:P) == --Yun's algorithm - see SYMSAC '76, p.27
--Note ci primitive is, so GCD's don't need to %do contents.
--Change gcd to return cofctrs also?
ci:=p; di:=differentiate(p); pi:=gcd(ci,di)
degree(pi)=0 =>
(u,c,a):=unitNormal(p)
makeFR(u,[["sqfr",c,1]])
i:NonNegativeInteger:=0; lffe:List FF:=[]
lcp := leadingCoefficient p
while degree(ci)^=0 repeat
ci:=(ci exquo pi)::P
di:=(di exquo pi)::P - differentiate(ci)
pi:=gcd(ci,di)
i:=i+1
degree(pi) > 0 =>
lcp:=(lcp exquo (leadingCoefficient(pi)**i))::RC
lffe:=[["sqfr",pi,i],:lffe]
makeFR(lcp::P,lffe)
else
squareFree(p:P) == --Musser's algorithm - see SYMSAC '76, p.27
--p MUST BE PRIMITIVE, Any characteristic.
--Note ci primitive, so GCD's don't need to %do contents.
--Change gcd to return cofctrs also?
ci := gcd(p,differentiate(p))
degree(ci)=0 =>
(u,c,a):=unitNormal(p)
makeFR(u,[["sqfr",c,1]])
di := (p exquo ci)::P
i:NonNegativeInteger:=0; lffe:List FF:=[]
dunit : P := 1
while degree(di)^=0 repeat
diprev := di
di := gcd(ci,di)
ci:=(ci exquo di)::P
i:=i+1
degree(diprev) = degree(di) =>
lc := (leadingCoefficient(diprev) exquo _
leadingCoefficient(di))::RC
dunit := lc**i * dunit
pi:=(diprev exquo di)::P
lffe:=[["sqfr",pi,i],:lffe]
dunit := dunit * di ** (i+1)
degree(ci)=0 => makeFR(dunit*ci,lffe)
redSqfr:=squareFree(divideExponents(ci,characteristic$P)::P)
lsnil:= [BumInSepFFE(ffe) for ffe in factorList redSqfr]
lffe:=append(lsnil,lffe)
makeFR(dunit*(unit redSqfr),lffe)
|