diff options
Diffstat (limited to 'src/algebra/padic.spad.pamphlet')
| -rw-r--r-- | src/algebra/padic.spad.pamphlet | 624 |
1 files changed, 624 insertions, 0 deletions
diff --git a/src/algebra/padic.spad.pamphlet b/src/algebra/padic.spad.pamphlet new file mode 100644 index 00000000..a75faa77 --- /dev/null +++ b/src/algebra/padic.spad.pamphlet @@ -0,0 +1,624 @@ +\documentclass{article} +\usepackage{axiom} +\begin{document} +\title{\$SPAD/src/algebra padic.spad} +\author{Clifton J. Williamson} +\maketitle +\begin{abstract} +\end{abstract} +\eject +\tableofcontents +\eject +\section{category PADICCT PAdicIntegerCategory} +<<category PADICCT PAdicIntegerCategory>>= +)abbrev category PADICCT PAdicIntegerCategory +++ Author: Clifton J. Williamson +++ Date Created: 15 May 1990 +++ Date Last Updated: 15 May 1990 +++ Basic Operations: +++ Related Domains: +++ Also See: +++ AMS Classifications: +++ Keywords: p-adic, completion +++ Examples: +++ References: +++ Description: This is the catefory of stream-based representations of +++ the p-adic integers. +PAdicIntegerCategory(p): Category == Definition where + p : Integer + I ==> Integer + NNI ==> NonNegativeInteger + ST ==> Stream + SUP ==> SparseUnivariatePolynomial + + Definition ==> Join(EuclideanDomain,CharacteristicZero) with + digits: % -> ST I + ++ \spad{digits(x)} returns a stream of p-adic digits of x. + order: % -> NNI + ++ \spad{order(x)} returns the exponent of the highest power of p + ++ dividing x. + extend: (%,I) -> % + ++ \spad{extend(x,n)} forces the computation of digits up to order n. + complete: % -> % + ++ \spad{complete(x)} forces the computation of all digits. + modulus: () -> I + ++ \spad{modulus()} returns the value of p. + moduloP: % -> I + ++ \spad{modulo(x)} returns a, where \spad{x = a + b p}. + quotientByP: % -> % + ++ \spad{quotientByP(x)} returns b, where \spad{x = a + b p}. + approximate: (%,I) -> I + ++ \spad{approximate(x,n)} returns an integer y such that + ++ \spad{y = x (mod p^n)} + ++ when n is positive, and 0 otherwise. + sqrt: (%,I) -> % + ++ \spad{sqrt(b,a)} returns a square root of b. + ++ Argument \spad{a} is a square root of b \spad{(mod p)}. + root: (SUP I,I) -> % + ++ \spad{root(f,a)} returns a root of the polynomial \spad{f}. + ++ Argument \spad{a} must be a root of \spad{f} \spad{(mod p)}. + +@ +\section{domain IPADIC InnerPAdicInteger} +<<domain IPADIC InnerPAdicInteger>>= +)abbrev domain IPADIC InnerPAdicInteger +++ Author: Clifton J. Williamson +++ Date Created: 20 August 1989 +++ Date Last Updated: 15 May 1990 +++ Basic Operations: +++ Related Domains: +++ Also See: +++ AMS Classifications: +++ Keywords: +++ Keywords: p-adic, completion +++ Examples: +++ References: +++ Description: +++ This domain implements Zp, the p-adic completion of the integers. +++ This is an internal domain. +InnerPAdicInteger(p,unBalanced?): Exports == Implementation where + p : Integer + unBalanced? : Boolean + I ==> Integer + NNI ==> NonNegativeInteger + OUT ==> OutputForm + L ==> List + ST ==> Stream + SUP ==> SparseUnivariatePolynomial + + Exports ==> PAdicIntegerCategory p + + Implementation ==> add + + PEXPR := p :: OUT + + Rep := ST I + + characteristic() == 0 + euclideanSize(x) == order(x) + + stream(x:%):ST I == x pretend ST(I) + padic(x:ST I):% == x pretend % + digits x == stream x + + extend(x,n) == extend(x,n + 1)$Rep + complete x == complete(x)$Rep + +-- notBalanced?:() -> Boolean +-- notBalanced?() == unBalanced? + + modP:I -> I + modP n == + unBalanced? or (p = 2) => positiveRemainder(n,p) + symmetricRemainder(n,p) + + modPInfo:I -> Record(digit:I,carry:I) + modPInfo n == + dv := divide(n,p) + r0 := dv.remainder; q := dv.quotient + if (r := modP r0) ^= r0 then q := q + ((r0 - r) quo p) + [r,q] + + invModP: I -> I + invModP n == invmod(n,p) + + modulus() == p + moduloP x == (empty? x => 0; frst x) + quotientByP x == (empty? x => x; rst x) + + approximate(x,n) == + n <= 0 or empty? x => 0 + frst x + p * approximate(rst x,n - 1) + + x = y == + st : ST I := stream(x - y) + n : I := _$streamCount$Lisp + for i in 0..n repeat + empty? st => return true + frst st ^= 0 => return false + st := rst st + empty? st + + order x == + st := stream x + for i in 0..1000 repeat + empty? st => return 0 + frst st ^= 0 => return i + st := rst st + error "order: series has more than 1000 leading zero coefs" + + 0 == padic concat(0$I,empty()) + 1 == padic concat(1$I,empty()) + + intToPAdic: I -> ST I + intToPAdic n == delay + n = 0 => empty() + modp := modPInfo n + concat(modp.digit,intToPAdic modp.carry) + + intPlusPAdic: (I,ST I) -> ST I + intPlusPAdic(n,x) == delay + empty? x => intToPAdic n + modp := modPInfo(n + frst x) + concat(modp.digit,intPlusPAdic(modp.carry,rst x)) + + intMinusPAdic: (I,ST I) -> ST I + intMinusPAdic(n,x) == delay + empty? x => intToPAdic n + modp := modPInfo(n - frst x) + concat(modp.digit,intMinusPAdic(modp.carry,rst x)) + + plusAux: (I,ST I,ST I) -> ST I + plusAux(n,x,y) == delay + empty? x and empty? y => intToPAdic n + empty? x => intPlusPAdic(n,y) + empty? y => intPlusPAdic(n,x) + modp := modPInfo(n + frst x + frst y) + concat(modp.digit,plusAux(modp.carry,rst x,rst y)) + + minusAux: (I,ST I,ST I) -> ST I + minusAux(n,x,y) == delay + empty? x and empty? y => intToPAdic n + empty? x => intMinusPAdic(n,y) + empty? y => intPlusPAdic(n,x) + modp := modPInfo(n + frst x - frst y) + concat(modp.digit,minusAux(modp.carry,rst x,rst y)) + + x + y == padic plusAux(0,stream x,stream y) + x - y == padic minusAux(0,stream x,stream y) + - y == padic intMinusPAdic(0,stream y) + coerce(n:I) == padic intToPAdic n + + intMult:(I,ST I) -> ST I + intMult(n,x) == delay + empty? x => empty() + modp := modPInfo(n * frst x) + concat(modp.digit,intPlusPAdic(modp.carry,intMult(n,rst x))) + + (n:I) * (x:%) == + padic intMult(n,stream x) + + timesAux:(ST I,ST I) -> ST I + timesAux(x,y) == delay + empty? x or empty? y => empty() + modp := modPInfo(frst x * frst y) + car := modp.digit + cdr : ST I --!! + cdr := plusAux(modp.carry,intMult(frst x,rst y),timesAux(rst x,y)) + concat(car,cdr) + + (x:%) * (y:%) == padic timesAux(stream x,stream y) + + quotientAux:(ST I,ST I) -> ST I + quotientAux(x,y) == delay + empty? x => error "quotientAux: first argument" + empty? y => empty() + modP frst x = 0 => + modP frst y = 0 => quotientAux(rst x,rst y) + error "quotient: quotient not integral" + z0 := modP(invModP frst x * frst y) + yy : ST I --!! + yy := rest minusAux(0,y,intMult(z0,x)) + concat(z0,quotientAux(x,yy)) + + recip x == + empty? x or modP frst x = 0 => "failed" + padic quotientAux(stream x,concat(1,empty())) + + iExquo: (%,%,I) -> Union(%,"failed") + iExquo(xx,yy,n) == + n > 1000 => + error "exquo: quotient by series with many leading zero coefs" + empty? yy => "failed" + empty? xx => 0 + zero? frst yy => + zero? frst xx => iExquo(rst xx,rst yy,n + 1) + "failed" + (rec := recip yy) case "failed" => "failed" + xx * (rec :: %) + + x exquo y == iExquo(stream x,stream y,0) + + divide(x,y) == + (z:=x exquo y) case "failed" => [0,x] + [z, 0] + + iSqrt: (I,I,I,%) -> % + iSqrt(pn,an,bn,bSt) == delay + bn1 := (empty? bSt => bn; bn + pn * frst(bSt)) + c := (bn1 - an*an) quo pn + aa := modP(c * invmod(2*an,p)) + nSt := (empty? bSt => bSt; rst bSt) + concat(aa,iSqrt(pn*p,an + pn*aa,bn1,nSt)) + + sqrt(b,a) == + p = 2 => + error "sqrt: no square roots in Z2 yet" + not zero? modP(a*a - (bb := moduloP b)) => + error "sqrt: not a square root (mod p)" + b := (empty? b => b; rst b) + a := modP a + concat(a,iSqrt(p,a,bb,b)) + + iRoot: (SUP I,I,I,I) -> ST I + iRoot(f,alpha,invFpx0,pPow) == delay + num := -((f(alpha) exquo pPow) :: I) + digit := modP(num * invFpx0) + concat(digit,iRoot(f,alpha + digit * pPow,invFpx0,p * pPow)) + + root(f,x0) == + x0 := modP x0 + not zero? modP f(x0) => + error "root: not a root (mod p)" + fpx0 := modP (differentiate f)(x0) + zero? fpx0 => + error "root: approximate root must be a simple root (mod p)" + invFpx0 := modP invModP fpx0 + padic concat(x0,iRoot(f,x0,invFpx0,p)) + + termOutput:(I,I) -> OUT + termOutput(k,c) == + k = 0 => c :: OUT + mon := (k = 1 => PEXPR; PEXPR ** (k :: OUT)) + c = 1 => mon + c = -1 => -mon + (c :: OUT) * mon + + showAll?:() -> Boolean + -- check a global Lisp variable + showAll?() == true + + coerce(x:%):OUT == + empty?(st := stream x) => 0 :: OUT + n : NNI ; count : NNI := _$streamCount$Lisp + l : L OUT := empty() + for n in 0..count while not empty? st repeat + if frst(st) ^= 0 then + l := concat(termOutput(n :: I,frst st),l) + st := rst st + if showAll?() then + for n in (count + 1).. while explicitEntries? st and _ + not eq?(st,rst st) repeat + if frst(st) ^= 0 then + l := concat(termOutput(n pretend I,frst st),l) + st := rst st + l := + explicitlyEmpty? st => l + eq?(st,rst st) and frst st = 0 => l + concat(prefix("O" :: OUT,[PEXPR ** (n :: OUT)]),l) + empty? l => 0 :: OUT + reduce("+",reverse_! l) + +@ +\section{domain PADIC PAdicInteger} +<<domain PADIC PAdicInteger>>= +)abbrev domain PADIC PAdicInteger +++ Author: Clifton J. Williamson +++ Date Created: 20 August 1989 +++ Date Last Updated: 15 May 1990 +++ Basic Operations: +++ Related Domains: +++ Also See: +++ AMS Classifications: +++ Keywords: +++ Keywords: p-adic, completion +++ Examples: +++ References: +++ Description: +++ Stream-based implementation of Zp: p-adic numbers are represented as +++ sum(i = 0.., a[i] * p^i), where the a[i] lie in 0,1,...,(p - 1). +PAdicInteger(p:Integer) == InnerPAdicInteger(p,true$Boolean) + +@ +\section{domain BPADIC BalancedPAdicInteger} +<<domain BPADIC BalancedPAdicInteger>>= +)abbrev domain BPADIC BalancedPAdicInteger +++ Author: Clifton J. Williamson +++ Date Created: 15 May 1990 +++ Date Last Updated: 15 May 1990 +++ Basic Operations: +++ Related Domains: +++ Also See: +++ AMS Classifications: +++ Keywords: p-adic, complementation +++ Examples: +++ References: +++ Description: +++ Stream-based implementation of Zp: p-adic numbers are represented as +++ sum(i = 0.., a[i] * p^i), where the a[i] lie in -(p - 1)/2,...,(p - 1)/2. +BalancedPAdicInteger(p:Integer) == InnerPAdicInteger(p,false$Boolean) + +@ +\section{domain PADICRC PAdicRationalConstructor} +<<domain PADICRC PAdicRationalConstructor>>= +)abbrev domain PADICRC PAdicRationalConstructor +++ Author: Clifton J. Williamson +++ Date Created: 10 May 1990 +++ Date Last Updated: 10 May 1990 +++ Basic Operations: +++ Related Domains: +++ Also See: +++ AMS Classifications: +++ Keywords: +++ Keywords: p-adic, completion +++ Examples: +++ References: +++ Description: This is the category of stream-based representations of Qp. +PAdicRationalConstructor(p,PADIC): Exports == Implementation where + p : Integer + PADIC : PAdicIntegerCategory p + CF ==> ContinuedFraction + I ==> Integer + NNI ==> NonNegativeInteger + OUT ==> OutputForm + L ==> List + RN ==> Fraction Integer + ST ==> Stream + + Exports ==> QuotientFieldCategory(PADIC) with + approximate: (%,I) -> RN + ++ \spad{approximate(x,n)} returns a rational number y such that + ++ \spad{y = x (mod p^n)}. + continuedFraction: % -> CF RN + ++ \spad{continuedFraction(x)} converts the p-adic rational number x + ++ to a continued fraction. + removeZeroes: % -> % + ++ \spad{removeZeroes(x)} removes leading zeroes from the + ++ representation of the p-adic rational \spad{x}. + ++ A p-adic rational is represented by (1) an exponent and + ++ (2) a p-adic integer which may have leading zero digits. + ++ When the p-adic integer has a leading zero digit, a 'leading zero' + ++ is removed from the p-adic rational as follows: + ++ the number is rewritten by increasing the exponent by 1 and + ++ dividing the p-adic integer by p. + ++ Note: \spad{removeZeroes(f)} removes all leading zeroes from f. + removeZeroes: (I,%) -> % + ++ \spad{removeZeroes(n,x)} removes up to n leading zeroes from + ++ the p-adic rational \spad{x}. + + Implementation ==> add + + PEXPR := p :: OUT + +--% representation + + Rep := Record(expon:I,pint:PADIC) + + getExpon: % -> I + getZp : % -> PADIC + makeQp : (I,PADIC) -> % + + getExpon x == x.expon + getZp x == x.pint + makeQp(r,int) == [r,int] + +--% creation + + 0 == makeQp(0,0) + 1 == makeQp(0,1) + + coerce(x:I) == x :: PADIC :: % + coerce(r:RN) == (numer(r) :: %)/(denom(r) :: %) + coerce(x:PADIC) == makeQp(0,x) + +--% normalizations + + removeZeroes x == + empty? digits(xx := getZp x) => 0 + zero? moduloP xx => + removeZeroes makeQp(getExpon x + 1,quotientByP xx) + x + + removeZeroes(n,x) == + n <= 0 => x + empty? digits(xx := getZp x) => 0 + zero? moduloP xx => + removeZeroes(n - 1,makeQp(getExpon x + 1,quotientByP xx)) + x + +--% arithmetic + + x = y == + EQ(x,y)$Lisp => true + n := getExpon(x) - getExpon(y) + n >= 0 => + (p**(n :: NNI) * getZp(x)) = getZp(y) + (p**((- n) :: NNI) * getZp(y)) = getZp(x) + + x + y == + n := getExpon(x) - getExpon(y) + n >= 0 => + makeQp(getExpon y,getZp(y) + p**(n :: NNI) * getZp(x)) + makeQp(getExpon x,getZp(x) + p**((-n) :: NNI) * getZp(y)) + + -x == makeQp(getExpon x,-getZp(x)) + + x - y == + n := getExpon(x) - getExpon(y) + n >= 0 => + makeQp(getExpon y,p**(n :: NNI) * getZp(x) - getZp(y)) + makeQp(getExpon x,getZp(x) - p**((-n) :: NNI) * getZp(y)) + + n:I * x:% == makeQp(getExpon x,n * getZp x) + x:% * y:% == makeQp(getExpon x + getExpon y,getZp x * getZp y) + + x:% ** n:I == + zero? n => 1 + positive? n => expt(x,n :: PositiveInteger)$RepeatedSquaring(%) + inv expt(x,(-n) :: PositiveInteger)$RepeatedSquaring(%) + + recip x == + x := removeZeroes(1000,x) + zero? moduloP(xx := getZp x) => "failed" + (inv := recip xx) case "failed" => "failed" + makeQp(- getExpon x,inv :: PADIC) + + inv x == + (inv := recip x) case "failed" => error "inv: no inverse" + inv :: % + + x:% / y:% == x * inv y + x:PADIC / y:PADIC == (x :: %) / (y :: %) + x:PADIC * y:% == makeQp(getExpon y,x * getZp y) + + approximate(x,n) == + k := getExpon x + (p :: RN) ** k * approximate(getZp x,n - k) + + cfStream: % -> Stream RN + cfStream x == delay +-- zero? x => empty() + invx := inv x; x0 := approximate(invx,1) + concat(x0,cfStream(invx - (x0 :: %))) + + continuedFraction x == + x0 := approximate(x,1) + reducedContinuedFraction(x0,cfStream(x - (x0 :: %))) + + termOutput:(I,I) -> OUT + termOutput(k,c) == + k = 0 => c :: OUT + mon := (k = 1 => PEXPR; PEXPR ** (k :: OUT)) + c = 1 => mon + c = -1 => -mon + (c :: OUT) * mon + + showAll?:() -> Boolean + -- check a global Lisp variable + showAll?() == true + + coerce(x:%):OUT == + x := removeZeroes(_$streamCount$Lisp,x) + m := getExpon x; zp := getZp x + uu := digits zp + l : L OUT := empty() + empty? uu => 0 :: OUT + n : NNI ; count : NNI := _$streamCount$Lisp + for n in 0..count while not empty? uu repeat + if frst(uu) ^= 0 then + l := concat(termOutput((n :: I) + m,frst(uu)),l) + uu := rst uu + if showAll?() then + for n in (count + 1).. while explicitEntries? uu and _ + not eq?(uu,rst uu) repeat + if frst(uu) ^= 0 then + l := concat(termOutput((n::I) + m,frst(uu)),l) + uu := rst uu + l := + explicitlyEmpty? uu => l + eq?(uu,rst uu) and frst uu = 0 => l + concat(prefix("O" :: OUT,[PEXPR ** ((n :: I) + m) :: OUT]),l) + empty? l => 0 :: OUT + reduce("+",reverse_! l) + +@ +\section{domain PADICRAT PAdicRational} +<<domain PADICRAT PAdicRational>>= +)abbrev domain PADICRAT PAdicRational +++ Author: Clifton J. Williamson +++ Date Created: 15 May 1990 +++ Date Last Updated: 15 May 1990 +++ Keywords: p-adic, complementation +++ Basic Operations: +++ Related Domains: +++ Also See: +++ AMS Classifications: +++ Keywords: p-adic, completion +++ Examples: +++ References: +++ Description: +++ Stream-based implementation of Qp: numbers are represented as +++ sum(i = k.., a[i] * p^i) where the a[i] lie in 0,1,...,(p - 1). +PAdicRational(p:Integer) == PAdicRationalConstructor(p,PAdicInteger p) + +@ +\section{domain BPADICRT BalancedPAdicRational} +<<domain BPADICRT BalancedPAdicRational>>= +)abbrev domain BPADICRT BalancedPAdicRational +++ Author: Clifton J. Williamson +++ Date Created: 15 May 1990 +++ Date Last Updated: 15 May 1990 +++ Keywords: p-adic, complementation +++ Basic Operations: +++ Related Domains: +++ Also See: +++ AMS Classifications: +++ Keywords: p-adic, completion +++ Examples: +++ References: +++ Description: +++ Stream-based implementation of Qp: numbers are represented as +++ sum(i = k.., a[i] * p^i), where the a[i] lie in -(p - 1)/2,...,(p - 1)/2. +BalancedPAdicRational(p:Integer) == + PAdicRationalConstructor(p,BalancedPAdicInteger p) + +@ +\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>> + +<<category PADICCT PAdicIntegerCategory>> +<<domain IPADIC InnerPAdicInteger>> +<<domain PADIC PAdicInteger>> +<<domain BPADIC BalancedPAdicInteger>> +<<domain PADICRC PAdicRationalConstructor>> +<<domain PADICRAT PAdicRational>> +<<domain BPADICRT BalancedPAdicRational>> +@ +\eject +\begin{thebibliography}{99} +\bibitem{1} nothing +\end{thebibliography} +\end{document} |