Hi all, have added sqrt for prime fields, perhaps for next release.
(1) -> )r sqrtpf F := PrimeField(nextPrime 273845092384750923487509234587) (1) PrimeField(273845092384750923487509234619) Type: Type x : F := 78345093845709384750923478509238475 (2) 203675171223548534986558619527 Type: PrimeField(273845092384750923487509234619) y := sqrt x (3) 86011080573889143280363302745 Type: PrimeField(273845092384750923487509234619) y^2-x (4) 0 Type: PrimeField(273845092384750923487509234619) x : F := 8596703498567034985679679876 (5) 8596703498567034985679679876 Type: PrimeField(273845092384750923487509234619) y := sqrt x (6) 53436731848479301994740043535 Type: PrimeField(273845092384750923487509234619) y^2-x (7) 0 Type: PrimeField(273845092384750923487509234619) x : F := 8596703498567034985679 (8) 8596703498567034985679 Type: PrimeField(273845092384750923487509234619) y := sqrt x >> Error detected within library code: sqrt: argument does not have a square root by Jacobi symbol. -- Mit freundlichen Grüßen Johannes Grabmeier Prof. Dr. Johannes Grabmeier Köckstraße 1, D-94469 Deggendorf Tel. +49-(0)-991-2979584, Tel. +49-(0)-151-681-70756 Tel. +49-(0)-991-3615-141 (d), Fax: +49-(0)-32224-192688 -- You received this message because you are subscribed to the Google Groups "FriCAS - computer algebra system" group. To unsubscribe from this group and stop receiving emails from it, send an email to [email protected]. To post to this group, send email to [email protected]. Visit this group at https://groups.google.com/group/fricas-devel. For more options, visit https://groups.google.com/d/optout.
)abbrev domain IPF InnerPrimeField ++ Authors: N.N., J.Grabmeier, A.Scheerhorn ++ Date Created: ?, November 1990, 26.03.1991 ++ Revision 05.06.2018 by J. Grabmeier ++ Basic Operations: ++ Related Constructors: PrimeField ++ Also See: ++ AMS Classifications: ++ Keywords: prime characteristic, prime field, finite field ++ References: ++ R.Lidl, H.Niederreiter: Finite Field, Encycoldia of Mathematics and ++ Its Applications, Vol. 20, Cambridge Univ. Press, 1983, ISBN 0 521 30240 4 ++ J. Grabmeier, A. Scheerhorn: Finite Fields in Axiom. ++ Axiom Technical Report Series, ATR/5 NP2522. ++ K. Pommerening: lectures on Cryptology, http://www.staff.uni-mainz.de/pommeren/Cryptology ++ Description: ++ InnerPrimeField(p) implements the field with p elements by using IntegerMod p. ++ Note: argument p MUST be a prime (this domain does not check). ++ See \spadtype{PrimeField} for a domain that does check. ++ In addition to the herited operations of IntegerMod p, the domain ++ provides exploits the structure of the cyclic group of its invertible elements. ++ It stores a primitive element w, i.a. generator of this group and it stores a ++ logarithm table for w as soon as this is required. ++ sqrt was added in 2018. InnerPrimeField(p : PositiveInteger) : Exports == Implementation where I ==> Integer NNI ==> NonNegativeInteger PI ==> PositiveInteger TBL ==> Table(PI, NNI) R ==> Record(key : PI, entry : NNI) SUP ==> SparseUnivariatePolynomial OF ==> OutputForm Exports ==> Join(FiniteFieldCategory, FiniteAlgebraicExtensionField(%), _ ConvertibleTo(Integer)) with sqrt : % -> % ++ sqrt(x) computes one y such that y^2 = x, error if there is no ++ square root, i.e. jacobi(x,p) = -1. ++ Implementation according to ++ http://www.staff.uni-mainz.de/pommeren/Cryptology/Asymmetric/5_NTh/ quadraticNonResidue: () -> % ++ quadraticNonResidue() computes the smallest non negative integer, ++ which represents a quadratic non residue. Implementation ==> IntegerMod p add -- local variables ============================================================ primitiveElementNotPresent? : Boolean := true -- gets false after initialization of the primitive Element. primitiveElt : PI := 1 -- primitiveElt stores a primitive element as soon as it is computed by -- createPrimitiveElement(). cyclicGroupSize : NNI := qcoerce(p-1)@NNI -- the size of the cyclic group of invertible elements. factorsOfCyclicGroupSize: List Record(factor : Integer, exponent : Integer) := [] -- will hold the factorization of the cyclic group size of the group of -- invertible elements, as soon as the logarithm tables are required. logarithmTableNotPresent? : Boolean := true -- gets false after initialization of the logarithm tables. discLogTable : Table(PI, TBL) := table()$Table(PI, TBL) -- tables indexed by the prime factors of the size q of the cyclic group -- discLogTable.factor is a table with keys -- primitiveElement() ^ (i * (q quo factor)) and entries i for -- i in 0..n-1, n computed in initialize() in order to use -- the minimal size limit 'limit' optimal. -- local functions ============================================================ twoPowerOfCyclicGroupSize : NNI := cGS : NNI := cyclicGroupSize r : NNI := 0 while cGS rem 2 = 0 repeat cGS := cGS quo 2 r := r+1 r -- local functions ============================================================= twoPower: % -> NNI ++ twoPower(x: %) computes the highest power of two in the factorization of ++ the order of x ~= 0. twoPower(x: %): NNI == x = 0 => error "twoPower: argument must not be 0." ord : NNI := order x r : NNI := 0 while ord rem 2 = 0 repeat ord := ord quo 2 r := r+1 r initializePrimitiveElement : () -> Void ++ initializePrimitiveElement() initializes the computation of a ++ primitive Element to be stored. initializePrimitiveElement(): Void == factorsOfCyclicGroupSize := factors(factor(cyclicGroupSize)$I)$(Factored I) -- get a primitive element: primitiveElt := lookup(createPrimitiveElement()) -- set initialization flag: primitiveElementNotPresent? := false void initializeLogarithmTable : () -> Void ++ initializeLogarithmTable() initializes the computation of a ++ logarithm table for the ++ primitive Element to be stored. initializeLogarithmTable(): Void == if primitiveElementNotPresent? then initializePrimitiveElement() -- now set up tables for discrete logarithm computations limit : Integer := 30 -- the minimum size for the discrete logarithm table for rec in factorsOfCyclicGroupSize repeat primeDivisor := rec.factor base : % := primitiveElement() ^ (cyclicGroupSize quo primeDivisor) l : Integer := length primeDivisor n : Integer := if odd? l then n := shift(primeDivisor, -(l quo 2)) else n := shift(1, (l quo 2)) if n < limit then d := (primeDivisor - 1) quo limit + 1 n := (primeDivisor - 1) quo d + 1 tbl : TBL := table()$TBL a : % := 1 for i in 0..(n-1)::NNI repeat insert!([lookup a , i::NNI]$R, tbl)$TBL a := a * base insert!([primeDivisor :: PI, copy(tbl)$TBL]$Record(key: PI, entry: TBL), discLogTable) -- set initialization flag: logarithmTableNotPresent? := false void -- exported functions ========================================================== -- exported functions from % =================================================== quadraticNonResidue(): % == found? : Boolean := false q : I := 1 while not found? repeat q := q+1 found? := jacobi(q, p)$IntegerNumberTheoryFunctions = -1 q :: % sqrt(x: %): % == zero? x => x jacobi(convert(x)@I, p)$IntegerNumberTheoryFunctions = -1 => error "sqrt: argument does not have a square root by Jacobi symbol." 3 = (p rem 4) => y : % := x ^ ((p+1) quo 4) -- we choose the smaller one of y and -y of their canonical -- representatives in 0..p-1 if convert(y)@I < convert(-y)@I then y else -y -- now we must have 1 = (p rem 4) and jacobi(x,p) = 1 (provided p is -- really a prime). b : % := quadraticNonResidue() e : NNI := twoPowerOfCyclicGroupSize u : NNI := ((p-1) :: NNI) quo (2^e) z : % := x r : NNI := twoPower z lr : List NNI := [r] while r > 0 repeat z := z * b ^ (2 ^ ( (e-r) :: NNI)) r := twoPower z lr := cons(r, lr) y : % := z ^ ( (u+1) quo 2) for r in rest lr repeat y := y / b ^ (2 ^ (e-r-1) :: NNI) if convert(y)@I < convert(-y)@I then y else -y generator() == 1 -- This uses Fermat's little theorem x^(p-1)=1 (mod p), so x^(q(p-1)+r) -- = x^r (mod p) x : % ^ n : Integer == zero?(n) => 1 zero?(x) => 0 r := positiveRemainder(n, p-1)::NNI ((x pretend IntegerMod p) ^$IntegerMod(p) r) pretend % if p <= convert(max()$SingleInteger)@Integer then q := p::SingleInteger recip x == zero?(y := convert(x)@Integer :: SingleInteger) => "failed" invmod(y, q)::Integer::% else recip x == zero?(y := convert(x)@Integer) => "failed" invmod(y, p)::% convert(x : %) == x pretend I normalElement() == 1 createNormalElement() == 1 characteristic() == p factorsOfCyclicGroupSize() == -- this fixes an infinite loop of function calls, problem was that -- factors factor(1) is the empty list: p = 2 => factorsOfCyclicGroupSize if empty? factorsOfCyclicGroupSize then initializePrimitiveElement() factorsOfCyclicGroupSize representationType() == "prime" tableForDiscreteLogarithm(fac) == if logarithmTableNotPresent? then initializeLogarithmTable() tbl := search(fac::PI, discLogTable)$Table(PI, TBL) tbl case "failed" => error "tableForDiscreteLogarithm: argument must be prime divisor_ of the order of the multiplicative group" tbl pretend TBL primitiveElement() == if primitiveElementNotPresent? then initializePrimitiveElement() index(primitiveElt) degree(x) : PI == 1::PositiveInteger extensionDegree() : PI == 1::PositiveInteger -- sizeOfGroundField() == p::NonNegativeInteger inGroundField?(x) == true coordinates(x) == new(1, x)$(Vector %) represents(v) == v.1 retract(x) == x retractIfCan(x) == x basis() == new(1, 1::%)$(Vector %) basis(n : PI) == n = 1 => basis() error("basis: argument must divide extension degree") definingPolynomial() == monomial(1, 1)$(SUP %) - monomial(1, 0)$(SUP %) minimalPolynomial(x) == monomial(1, 1)$(SUP %) - monomial(x, 0)$(SUP %) charthRoot x == x )abbrev domain PF PrimeField ++ Authors: N.N., ++ Date Created: November 1990, 26.03.1991 ++ Basic Operations: ++ Related Constructors: ++ Also See: ++ AMS Classifications: ++ Keywords: prime characteristic, prime field, finite field ++ References: ++ R.Lidl, H.Niederreiter: Finite Field, Encycoldia of Mathematics and ++ Its Applications, Vol. 20, Cambridge Univ. Press, 1983, ISBN 0 521 30240 4 ++ Description: ++ PrimeField(p) implements the field with p elements if p is a ++ prime number. ++ Error: if p is not prime. --++ with new compiler, want to put the error check before the add PrimeField(p : PositiveInteger) : Exp == Impl where Exp ==> Join(FiniteFieldCategory, FiniteAlgebraicExtensionField(%), _ ConvertibleTo(Integer)) with sqrt : % -> % ++ sqrt(x) computes one y such that y^2 = x, error if there is no square ++ root, i.e. jacobi(x,p) = -1. Impl ==> InnerPrimeField(p) add if not prime?(p)$IntegerPrimesPackage(Integer) then error "Argument to prime field must be a prime" --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.
F := PrimeField(nextPrime 273845092384750923487509234587) x : F := 78345093845709384750923478509238475 y := sqrt x y^2-x x : F := 8596703498567034985679679876 y := sqrt x y^2-x x : F := 8596703498567034985679 y := sqrt x
