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
------------------------------------------------------------------------------
The demand for IT networking professionals continues to grow, and the
demand for specialized networking skills is growing even more rapidly.
Take a complimentary Learning@Cisco Self-Assessment and learn
about Cisco certifications, training, and career opportunities.
http://p.sf.net/sfu/cisco-dev2dev
_______________________________________________
open-axiom-devel mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/open-axiom-devel
