thanks for the honor - or not: I am not the author of this domain.
It is a "mutable domain"
)bo PUSH('ModMonic, $mutableDomains)
I used that, when we extended finite fields, along the lines of the existing
standard implemention at that time.
Am 24.10.2011 um 14:54 schrieb Gabriel Dos Reis:
>
> Consider the definition of ModMonoid:
> --------------8<---------------8<----------------8<--------------
> )abbrev domain MODMON ModMonic
> ++ Description:
> ++ This package \undocumented
>
> ModMonic(R,Rep): C == T
> where
> R: Ring
> Rep: UnivariatePolynomialCategory(R)
> C == Join(UnivariatePolynomialCategory(R),CoercibleFrom Rep) with
> --operations
> setPoly : Rep -> Rep
> ++ setPoly(x) \undocumented
> modulus : -> Rep
> ++ modulus() \undocumented
> reduce: Rep -> %
> ++ reduce(x) \undocumented
> lift: % -> Rep --reduce lift = identity
> ++ lift(x) \undocumented
> Vectorise: % -> Vector(R)
> ++ Vectorise(x) \undocumented
> UnVectorise: Vector(R) -> %
> ++ UnVectorise(v) \undocumented
> An: % -> Vector(R)
> ++ An(x) \undocumented
> pow : -> PrimitiveArray(%)
> ++ pow() \undocumented
> computePowers : -> PrimitiveArray(%)
> ++ computePowers() \undocumented
> if R has FiniteFieldCategory then
> frobenius: % -> %
> ++ frobenius(x) \undocumented
> --LinearTransf: (%,Vector(R)) -> SquareMatrix<deg> R
> --assertions
> if R has Finite then Finite
> T == add
> --constants
> m:Rep := monomial(1,1)$Rep --| degree(m) > 0 and LeadingCoef(m) =
> R$1
> d := degree(m)$Rep
> d1 := (d-1):NonNegativeInteger
> twod := 2*d1+1
> frobenius?:Boolean := R has FiniteFieldCategory
> --VectorRep:= DirectProduct(d:NonNegativeInteger,R)
> --declarations
> x,y: %
> p: Rep
> d,n: Integer
> e,k1,k2: NonNegativeInteger
> c: R
> --vect: Vector(R)
> power:PrimitiveArray(%)
> frobeniusPower:PrimitiveArray(%)
> computeFrobeniusPowers : () -> PrimitiveArray(%)
> --representations
> --mutable m --take this out??
> --define
> power := new(0,0)
> frobeniusPower := new(0,0)
> setPoly (mon : Rep) ==
> mon =$Rep m => mon
> oldm := m
> not one? leadingCoefficient mon => error "polynomial must be
> monic"
> -- following copy code needed since FFPOLY can modify mon
> copymon:Rep:= 0
> while not zero? mon repeat
> copymon := monomial(leadingCoefficient mon, degree mon)$Rep +
> copymon
> mon := reductum mon
> m := copymon
> d := degree(m)$Rep
> d1 := (d-1)::NonNegativeInteger
> twod := 2*d1+1
> power := computePowers()
> if frobenius? then
> degree(oldm)>1 and not((oldm exquo$Rep m) case "failed") =>
> for i in 1..d1 repeat
> frobeniusPower(i) := reduce lift frobeniusPower(i)
> frobeniusPower := computeFrobeniusPowers()
> m
> modulus == m
> if R has Finite then
> size == d * size()$R
> random == UnVectorise([random()$R for i in 0..d1])
> 0 == 0$Rep
> 1 == 1$Rep
> c * x == c *$Rep x
> n * x == (n::R) *$Rep x
> coerce(c:R):% == monomial(c,0)$Rep
> coerce(x:%):OutputForm == coerce(x)$Rep
> coefficient(x,e):R == coefficient(x,e)$Rep
> reductum(x) == reductum(x)$Rep
> leadingCoefficient x == (leadingCoefficient x)$Rep
> degree x == (degree x)$Rep
> lift(x) == x pretend Rep
> reduce(p) == (monicDivide(p,m)$Rep).remainder
> coerce(p) == reduce(p)
> x = y == x =$Rep y
> x + y == x +$Rep y
> - x == -$Rep x
> x * y ==
> p := x *$Rep y
> ans:=0$Rep
> while (n:=degree p)>d1 repeat
> ans:=ans + leadingCoefficient(p)*power.(n-d)
> p := reductum p
> ans+p
> Vectorise(x) == [coefficient(lift(x),i) for i in 0..d1]
> UnVectorise(vect) ==
> reduce(+/[monomial(vect.(i+1),i) for i in 0..d1])
> computePowers ==
> mat : PrimitiveArray(%):= new(d,0)
> mat.0:= reductum(-m)$Rep
> w: % := monomial$Rep (1,1)
> for i in 1..d1 repeat
> mat.i := w *$Rep mat.(i-1)
> if degree mat.i=d then
> mat.i:= reductum mat.i + leadingCoefficient mat.i * mat.0
> mat
> if frobenius? then
> computeFrobeniusPowers() ==
> mat : PrimitiveArray(%):= new(d,1)
> mat.1:= mult := monomial(1, size()$R)$%
> for i in 2..d1 repeat
> mat.i := mult * mat.(i-1)
> mat
>
> frobenius(a:%):% ==
> aq:% := 0
> while not zero? a repeat
> aq:= aq + leadingCoefficient(a)*frobeniusPower(degree a)
> a := reductum a
> aq
>
> pow == power
> monomial(c,e)==
> if e<d then monomial(c,e)$Rep
> else
> if e<=twod then
> c * power.(e-d)
> else
> k1:=e quo twod
> k2 := (e-k1*twod)::NonNegativeInteger
> reduce((power.d1 **k1)*monomial(c,k2))
> if R has Field then
>
> (x:% exquo y:%):Union(%, "failed") ==
> uv := extendedEuclidean(y, modulus(), x)$Rep
> uv case "failed" => "failed"
> return reduce(uv.coef1)
>
> recip(y:%):Union(%, "failed") == 1 exquo y
> divide(x:%, y:%) ==
> (q := (x exquo y)) case "failed" => error "not divisible"
> [q, 0]
>
> -- An(MM) == Vectorise(-(reduce(reductum(m))::MM))
> -- LinearTransf(vect,MM) ==
> -- ans:= 0::SquareMatrix<d>(R)
> -- for i in 1..d do setelt(ans,i,1,vect.i)
> -- for j in 2..d do
> -- setelt(ans,1,j, elt(ans,d,j-1) * An(MM).1)
> -- for i in 2..d do
> -- setelt(ans,i,j, elt(ans,i-1,j-1) + elt(ans,d,j-1) *
> An(MM).i)
> -- ans
>
> --------------8<---------------8<----------------8<--------------
>
> Its description does not have a claim for "author", but the intersection
> of authors of domains using this domain leaves me with Johannes Grabmeier :-)
>
>
> ModMonic is not a builtin domain, therefore one needs a definition of its
> representation, i.e. a definition for "Rep", is needed. Yet, none is
> provided. Rather, "Rep" here is specified as a parameter to the functor.
> (That is not the same as defining it locally to some value)
>
> Now, consider the capsule-level variable "power". It is declared as
>
> power:PrimitiveArray(%)
>
> then later assigned to as
>
> power := new(0,0)
>
> There are two lexical occurences of the literal '0' in there. Can you
> guess, which one is which? There are at least five interpretations for
> each constant:
>
> 0: R -- the parameter R satisfies Ring
> 0: Rep -- the parameter Rep satisfies UnivariatePolynomialCategory(R)
> 0: % -- the current domains satisfies same
> 0: NonNegativeInteger
> 0: Integer
>
>
> That this domain compiles at all is mysterious; it happens to depend on
> a very obscure implementation detail of AXIOM compilers, one that is not
> advertised at all -- and I don't think it should.
> It is voodoo. So, I am calling it Johannes's Rep voodoo :-).
>
> Brownie point to anyone who untangles it not to use the voodoo.
>
> NB: I ran into this domain while working on a completely different
> issue, one that related to fixing several shortcomings in the way
> constants are processed in many AXIOM systems.
>
> -- Gaby
>
> --
> You received this message because you are subscribed to the Google Groups
> "FriCAS - computer algebra system" group.
> To post to this group, send email to [email protected].
> To unsubscribe from this group, send email to
> [email protected].
> For more options, visit this group at
> http://groups.google.com/group/fricas-devel?hl=en.
>
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)-171-5503789
Tel. +49-(0)-991-3615-100 (d), Fax: +49-(0)-1803-5518-17745
--
You received this message because you are subscribed to the Google Groups
"FriCAS - computer algebra system" group.
To post to this group, send email to [email protected].
To unsubscribe from this group, send email to
[email protected].
For more options, visit this group at
http://groups.google.com/group/fricas-devel?hl=en.
