newPackage( "agcodes", Version => "0.1", Date => "October 5, 2008", Authors => {{Name => "Nathan Ilten", Email => "nilten@cs.uchicago.edu>"}}, HomePage => "http://people.cs.uchicago.edu/~nilten/", Headline => "A Macaulay2 package for computing with AG codes", DebuggingMode => true ) export{ Divisor, RPoint, ZeroCycle, divisor, normalForm, rpoints, pointideal, zerocycle, globalSections, sectionIdeal, linearlyEquivalent, effective, agcode, weights, localparameter, tochart, localorder, localremainder, iseffective, fieldelements, mytuples, projectiveorder } -- New types Divisor = new Type of BasicList RPoint = new Type of List ZeroCycle = new Type of List -- New methods rpoints=method() divisor = method() normalForm = method() globalSections = method() pointideal= method() zerocycle= method() effective= method() linearlyEquivalent = method() agcode = method() weights = method() localparameter = method() tochart = method() localorder = method() projectiveorder = method() localremainder = method() iseffective = method() --Purpose: Compute rational points of an affine or projective variety --Input: A Variety --Output: A list of RPoints rpoints(Variety):=(X) -> ( R:= ring(X); B:= ambient(R); I:= ideal(X); K:= coefficientRing R; n:= dim ambient R; L:=fieldelements K; points:={}; if n==0 then {} else if class(X)===AffineVariety then ( Aff:=mytuples(L,n); for f in Aff do ( F:=map(K,B,f); if F(I)==0 then points=points|{new RPoint from f}; ); points ) else if n==1 then ( -- if substitute(I,matrix{{1}})==0 then {{1}} if I==0 then {{1}} else {} ) else ( Cords:=mytuples(L,n-1); Chart:=apply(Cords,i->{1}|i); for f in Chart do ( F:=map(K,B,f); if F(I)==0 then points=points|{new RPoint from f}; ); S:=K[drop(generators(ambient(R)),1)]; edgepoints:=apply(rpoints(Proj(S/sub(I,S))),i-> new RPoint from {0}|i); points|edgepoints ) ) -- Purpose: Create a ZeroCycle from an RPoint or from a Divisor -- Input: An RPoint or a Divisor or a Divisor and List of rational points -- Output: A ZeroCycle zerocycle RPoint := P -> new ZeroCycle from {{1,P}} zerocycle(Divisor, List):=(D,L)-> ( Z:= new ZeroCycle from {}; for P in L do ( Z=Z+(projectiveorder(D#0,P)-projectiveorder(D#1,P))*P; ); Z ) zerocycle(Divisor):=D->zerocycle(D,rpoints Proj ring (D#0)) -- We define how to add and subtract ZeroCycles ZeroCycle + ZeroCycle := (W,Z) -> ( SUPP:=unique(flatten(subtable(toList(0..(#W-1)),{1},W))|flatten(subtable(toList(0..(#Z-1)),{1},Z))); newcycle:={}; for P in SUPP do ( coeffs:=positions(W|Z,i->i_1==P); coeff:=sum(apply(coeffs,i->((W|Z)_i)_0)); if coeff!=0 then newcycle=newcycle|{{coeff,P}}; ); new ZeroCycle from newcycle) - ZeroCycle := Z -> new ZeroCycle from apply(Z,i->{-i_0,i_1}) ZeroCycle - ZeroCycle := (W,Z) -> W+(-Z) ZZ ZeroCycle:= ZZ * ZeroCycle:= (n,Z) -> new ZeroCycle from apply(Z,i->{i_0*n,i_1}) -- We define how to add and subtract RPoints, resulting in ZeroCycles ZZ RPoint := ZZ * RPoint := (n,P) -> new ZeroCycle from {{n,P}} - RPoint := P -> new ZeroCycle from {{-1,P}} RPoint + RPoint:= (P,Q)->zerocycle(P)+zerocycle(Q) RPoint - RPoint:= (P,Q)->zerocycle(P)-zerocycle(Q) ZeroCycle + RPoint:= (P,Q)->P+zerocycle(Q) RPoint + ZeroCycle:= (P,Q)->zerocycle(P)+Q ZeroCycle - RPoint:= (P,Q)->P-zerocycle(Q) RPoint - ZeroCycle:= (P,Q)->zerocycle(P)-Q -- Purpose: Describe the ideal corresponding to a rational point -- Input: A projective variety X and a rational point P -- Output: The corresponding homogeneous maximal ideal pointideal(RPoint,ProjectiveVariety):= (P,X) -> ( p:=position(P,i->i!=0); ideal(apply(generators(ring(X)),i->i*P_p)-apply(P,i->i*((generators(ring(X)))_p))) ) --Purpose: Compute the degree of a ZeroCycle --Input: A ZeroCycle or RPoint --Ouput: Whatever kind of coefficients we had degree(ZeroCycle):=Z ->(sum apply(Z,L->L_0)) degree(RPoint):=P->degree(zerocycle(P)) -- Some of the following was adapted from the Divisors tutorial -- Purpose: Create a new Divisor on a projective variety Y -- Input: Either one or two ideals in the hom. coordinate ring of Y OR an RPoint or ZeroCycle and the ProjectiveVariety Y -- Output: The corresponding Divisor represented as a BasicList of two Ideals -- Note: We assume that Y is S2 and, for RPoint or ZeroCycle input, that Y is a curve. This isn't checked. divisor(Ideal,Ideal) := (I,J) -> new Divisor from {purify1S2 I,purify1S2 J} divisor Ideal := I -> divisor(I, ideal 1_(ring I)) divisor(RPoint,ProjectiveVariety) := (P,X) -> divisor(pointideal(P,X), ideal 1_(ring X)) divisor(ZeroCycle,ProjectiveVariety) := (Z,X) -> sum(apply(Z,i->((i_0)*divisor(i_1,X)))) -- Define a normal form, equality, and sums and differences of divisors normalForm Divisor := D -> new Divisor from {D#0 : D#1, D#1 : D#0} Divisor == Divisor := (D,E) -> toList normalForm D == toList normalForm E Divisor + Divisor := (D,E) -> divisor(D#0 * E#0, D#1 * E#1) - Divisor := (D) -> new Divisor from {D#1, D#0} Divisor - Divisor := (D,E) -> D + (-E) ZZ Divisor := ZZ * Divisor := (n,D) -> (if n>=0 then divisor((D#0)^n, (D#1)^n) else -((-n)*D)) -- Had to change this one to allow negative n --Purpose: Compute the global sections of a Divisor --Input: Either a Divisor OR a ZeroCycle or an RPoint and a ProjectiveVariety --Output: A list, the first element being a row vector containg the numerators of a basis, and the second element being a common denominator globalSections Divisor := (D) -> ( I := D#0; J := D#1; f := (gens I)_(0,0); LD := basis(degree f, purify1S2((f*J) : I)); LD = super (LD ** (ring target LD)); {LD, f}) globalSections (ZeroCycle,ProjectiveVariety):= (Z,Y) -> globalSections divisor (Z,Y) globalSections(RPoint,ProjectiveVariety):= globalSections(ZeroCycle,ProjectiveVariety):= (Z,X)->globalSections(divisor(Z,X)) --Purpose: Compute the homogeneous Ideal corresponding to a section --Input: Two homogeneous polynomials and a Divisor --Output: An Ideal in the homogeneous coordinate ring sectionIdeal = (f,g,D) -> ( I := D#0; J := D#1; purify1S2((f*I):g) : J ) --Purpose: Test if two Divisors are linearly equivalent --Input: Two Divisors OR Two ZeroCycles and a ProjectiveVariety OR Two RPoints and a ProjectiveVariety --OUtput: Boolean linearlyEquivalent(Divisor,Divisor):= (D,E) -> ( --we make this into a method F := normalForm(D-E); LB := globalSections F; L := LB#0; if numgens source L != 1 then false else ( R := ring L; V := sectionIdeal(L_(0,0), LB#1, F); if V == ideal(1_R) then L_(0,0)/LB#1 else false) ) linearlyEquivalent(ZeroCycle,ZeroCycle,ProjectiveVariety):=(W,Z,Y)->linearlyEquivalent(divisor(W,Y),divisor(Z,Y)) linearlyEquivalent(RPoint,RPoint,ProjectiveVariety):=linearlyEquivalent(ZeroCycle,ZeroCycle,ProjectiveVariety):=(W,Z,X)->linearlyEquivalent(divisor(W,X),divisor(Z,X)) --Purpose: Return an element of |D| --Input: A Divisor OR a RPoint and a ProjectiveVariety OR a ZeroCycle and a ProjectiveVariety --Output: An effective Divisor effective Divisor := (D) -> ( LB := globalSections D; L := LB#0; -- the matrix of numerators if numgens source L == 0 then error(toString D + " is not effective") else divisor sectionIdeal(L_(0,0), LB#1, D)) effective(RPoint,ProjectiveVariety):=effective(ZeroCycle,ProjectiveVariety):=(W,X)->effective(divisor(W,X)) --Purpose: Check if a divisor is effective --Input: A Divisor OR a RPoint and a ProjectiveVariety OR a ZeroCycle and a ProjectiveVariety --Output: A boolean iseffective Divisor := D -> ( if not D#1==1 then return false; true) iseffective(RPoint,ProjectiveVariety):=iseffective(ZeroCycle,ProjectiveVariety):= (Z,Y)-> iseffective(divisor(Z,Y)) --Purpose: Move to an affine chart in projective space --Input: A ProjectiveVariety and RPoint OR a homogeneous coordinate Ring and an RPoint --Output: A List consisting of a Map from the homogeneous coordinate ring to an affine coordinate ring, and a new RPoint in affine coordinates tochart(ProjectiveVariety,RPoint):=(Y,P)->( R:=ring(Y); S:= ambient R; K:=coefficientRing S; i:=position(P,p->(p!=0)); --first we find a nonzero coordinate S1:=K[drop(generators S, {i,i})]; -- then we switch to that affine chart F:=map(S1,S,(generators S1)_{0..(i-1)}|{P_i}|(generators S1)_{i..(#P-2)}); R1:=S1/F(ideal(R)); FF:=map(R1,R,(generators R1)_{0..(i-1)}|{P_i}|(generators R1)_{i..(#P-2)}); P1:=new RPoint from P_{0..(i-1)}|P_{(i+1)..(#P-1)}; {FF,P1} ) tochart(Ring,RPoint):=(R,P)->tochart(Proj(R),P) --Purpose: Compute a local parameter in the local ring at a point --Input: An affine coordinate Ring and and RPoint --Output: An element of the ring which locally generates the maximal ideal --Note: If the tanget space at the point isn't one-dimensional we return an error since the whole thing makes no sense. localparameter(Ring,RPoint):=(R1,P1)->( m:=ideal(generators(R1)-P1); --maximal ideal T:=m/(m*m); B:=basis(T); if numgens source B != 1 then error("Tangent space isn't one-dimensional"); (entries B_0)_0 ) --Purpose: Compute the order of a polynomial or ideal in a point --Input: A polynomial as RingElement or an Ideal, and an RPoint --Output: An integer (or infinity) localorder(Ideal,RPoint):=(I,P1)->( if I==0 then return infinity; R1=ring(I); m:=ideal(generators(R1)-P1); o:=0; while isSubset(I,m^(o+1)) do o=o+1; o ) localorder(RingElement,RPoint):=(f,P1)->localorder(ideal f,P1) --Purpose: Compute the order of a polynomial or ideal in a point on a projective variety --Input: A polynomial as RingElement or an Ideal, and an RPoint --Output: An integer (or infinity) projectiveorder(Ideal,RPoint):=(I,P)->( CHART:=tochart(ring I,P); localorder(CHART_0(I),CHART_1) ) projectiveorder(RingElement,RPoint):=(f,P)->projectiveorder(ideal f,P) --Purpose: Compute the first non-disappearing image of a polynomial in m/m^k --Input: The polynomial as RingElement, a local parameter as RingElement, an RPoint, and the local order as integer --Output: An element of the coefficient field localremainder(RingElement,RingElement,RPoint,ZZ):=(f,t,P1,o)->( R1:=ring(f); K:=coefficientRing(R1); m:=ideal(generators(R1)-P1); i:=0; L:=fieldelements coefficientRing ambient R1; while localorder(f-(L_i)*(t^o),P1)==o do ( i=i+1); (L_i)_K ) --Purpose: Evaluate a list of rational functions in a point --Input: An RPoint and a List of homogeneous rational functions --Output: A List of values in the coefficient field pointeval:=(P,L)->( CHART:=tochart(ring(numerator(L_0)),P); FF:=CHART_0; P1:=CHART_1; L1:=apply(L,s->{FF(numerator(s)),FF(denominator(s))});--switch to chart and split into num/denom t:=localparameter(ring(L1_0_0),P1); vals:=apply(L1,l->( o0:=localorder(l_0,P1); o1:=localorder(l_1,P1); if o0 ( SUPP:=flatten(subtable(toList(0..(#Z-1)),{1},Z)); G:=globalSections(D); H:=flatten(entries(G_0)); g:=G_1; sections:=apply(H,f->(f/g)); R:= ring(D#0); K:= coefficientRing R; output:={}; for P in SUPP do (output=output|{((pointeval(P,sections)))}); M:=transpose matrix output; M ) --Purpose: Calculate the weight vector of a linear code --Input: A generator Matrix --Output: A List corresponding to the weight vector --Note: This only works for very small codes weights Matrix := M -> ( K:=ring M; L:=fieldelements K; UR:=mytuples(L,numgens target M); wts:=apply(UR,v->(length positions(flatten(entries((matrix {v})*M)),i->i!=0))); apply(toList(0..(numgens source M)),i->#positions(wts,j->j==i)) ) --Some further helpful functions mytuples= (L,n) -> ( N:= pack(1,L); for h in 1..n-1 do ( N= flatten table(L,N,(i,j)->append(j,i)) ); N) fieldelements= K -> ( p:= char K; d:=1;e:=1; if class(K)===GaloisField then (d= degree ideal ambient K;e=K_0); L:=matrix(mytuples(toList(0..p-1),d))*matrix(pack(1,apply(toList(0..d-1),i->e^i))); flatten entries L) purify1S2 = I -> ( M := compress gens I; if numgens source M == 0 then error "purify1S2: expected nonzero ideal"; f := ideal M_(0,0); f:(f:I) ) beginDocumentation() document { Key => Points, "A Macaulay2 package for computing with AG codes", PARA{}, } end