// copy and paste the first 210 lines as one block into a Magma, then // play with the worked examples from line 212 ffoll. function Pts(r,A) return [[a*k mod r : a in A] : k in [1..r-1]]; end function; function JunSuff(r,A) Points := Pts(r,A); Jun := [P : P in Points | &+P eq r]; Mid := [P : P in Points | &+P eq 2*r]; Sums := [[P1[i]+P2[i] : i in [1..#A]] : P1 in Jun, P2 in Jun]; if &and[P in Sums : P in Mid] then return true; else return false, [P : P in Mid | P notin Sums][1]; end if; end function; function Solutions(r); TempSolutions := [[0,0,0,0]]; // initiate list with nonsense for a1 in [i : i in [0..r div 4]], a2 in [i : i in [a1..(r-a1) div 3]], a3 in [i : i in [a2..(r-a1-a2) div 2]] do a4 := r-a1-a2-a3; A := [a1,a2,a3,a4]; if (a4 ne r) and (GreatestCommonDivisor([r] cat A) eq 1) and (&and[A ne Sort(B) : B in Pts(r,C), C in TempSolutions]) and (JunSuff(r,A)) then Append(~TempSolutions,A); end if; end for; return Remove(TempSolutions,1); // take out initial nonsense end function; /* r := 23; A := [1,2,4,16]; //try it for 1/23(1,2,4,16) Jun := [P : P in Pts(r,A) | &+P eq r]; Jun; [B : B in Jun | &and [JunSuff(B[i], [-r] cat Remove(B,i)) : i in [1..4]]]; */ KK := FiniteField(1009); RR := PolynomialRing(KK,4); function InvariantMonomials(r,A) Inv := [x^r,y^r,z^r,t^r,x*y*z*t]; for i in [0..r-1], j in [0..r-i-1], k in [0..r-i-j-1] do if i*j*k eq 0 then Append(~Inv,x^((-A[2]*i-A[3]*j-A[4]*k) mod r)*y^i*z^j*t^k); elif (A[2]*i+A[3]*j+A[4]*k) mod r eq 0 then Append(~Inv,y^i*z^j*t^k); end if; end for; Exclude(~Inv,1); // Omit 1 before testing divisibility. S := [[RR!1]]; // S is a list of lists Si of mons of deg i. for i in [1..r] do Si := [m : m in Inv | TotalDegree(m) eq i]; Append(~S,Si); // Omit those divisible by mons of deg i. Inv := [m : m in Inv | not &or[IsDivisibleBy(m,n) : n in Si]]; end for; if #Inv ne 0 then error("Invariant monomial has degree > r"); end if; return S; end function; function EigenSp(r,A) Inv := InvariantMonomials(r,A); Eig := [&cat[S : S in Inv]]; // Collate the list of lists. Rbar := quo; for s in [1..r-1] do temp := [x^((s-A[2]*i-A[3]*j-A[4]*k) mod r)*y^i*z^j*t^k : k in [0..r-i-j-1], j in [0..r-i-1], i in [0..r-1] | i*j*k eq 0] cat [y^i*z^j*t^k : k in [0..r-i-j-1], j in [0..r-i-1], i in [0..r-1] | (i*j*k ne 0) and ((A[2]*i+A[3]*j+A[4]*k) mod r eq s)]; temp := [m : m in temp | Rbar!m ne 0]; Append(~Eig,temp); end for; return Rotate(Eig,-1); // Put invariants in final rth place. end function; function NonZeroMonomials(Rbar,S) return [m : m in S | Rbar!m ne 0]; end function; function NotOne(R,E) return [[i] cat L : i in [1..#E] | #L ne 1 where L is NonZeroMonomials(R,E[i])]; end function; function ASets(r,A) Eig := EigenSp(r,A); // initialise dynamic constructions Bag := []; // To store the fruit. finished := false; I := []; // List of indexes i1, i2, ... where I make choices. C := []; // List of lists of candidates mi in Eig[i], i in I. // Cycle through monomials in last list then delete it. M := []; // List of number of current choices from C. // In other words, current choices are C[i,M[i]] for i in I. while not finished do // Make the current quotient ring by killing all monomials // in all Eig[i] except the current choices mi for i in I. Rbar := quo< RR | Exclude(Eig[r],1) cat &cat[Exclude(Eig[I[i]],C[i,M[i]]) : i in [1..#I]] >; // Analyse whether it is an A-cluster. tlist := NotOne(Rbar,Eig); if tlist eq [] then Include(~Bag,Basis(DivisorIdeal(Rbar))); end if; // If so, bag the fruit (intending to backtrack). // If each tlist[i] is [i, plus >= 1 monomial], so no // contradiction or result, go forward: add new index to I, // new list of candidates to C, and first choice to M. if (tlist ne []) and (&and[#l ne 1 : l in tlist]) then Append(~I,Integers()!tlist[1,1]); Append(~C,Remove(tlist[1],1)); Append(~M,1); end if; // I,M; // Uncomment this line for a diagnostic display. if (tlist eq []) or ((tlist ne []) and (&or[#l eq 1 : l in tlist])) then // Backtrack, include setting finished to true if can't. ss := #I; if &and[M[i] eq #(C[i]) : i in [1..ss]] then finished := true; else tt := Max([i : i in [1..ss] | M[i] ne #(C[i])]); I:=I[1..tt]; M:=M[1..tt]; M[tt]:=M[tt]+1; C:=C[1..tt]; end if; end if; end while; return Bag, Eig; end function; function EigenBasis(S,E) Rbar := quo; MonB := MonomialBasis(Rbar); return [[n : n in MonB | n in E[i]][1] : i in [1..#E]]; end function; function RelMatrix(rr,AA,SS,Eig) Rbar := quo; EB := EigenBasis(SS,Eig); M := Matrix(4,[Integers()|1,1,1,1]); // add invariant x*y*z*t for m in SS do if m eq x*y*z*t then; // don't add invariant x*y*z*t twice else B := [Degree(m,RR.i) : i in [1..4]]; c := &+[AA[i]*B[i] : i in [1..4]] mod rr; if c eq 0 then c := rr; end if; n := EB[c]; C := [Degree(n,RR.i) : i in [1..4]]; M := VerticalJoin(M,Matrix(4,[B[i]-C[i] : i in [1..4]])); end if; end for; return M; end function; // given matrix N, ask for index of positive entry in a row // with other entries negative, or vice versa. Return 0 if none. function convex_row(N) n := NumberOfRows(N); J := 0; for i in [1..n] do if J eq 0 then pos := [j : j in [1..n+4] | N[i,j] gt 0]; neg := [j : j in [1..n+4] | N[i,j] lt 0]; if (#neg eq 1) and (#pos gt 1) then found := true; J := neg[1]; elif (#pos eq 1) and (#neg gt 1) then found := true; J := pos[1]; end if; end if; end for; return J; end function; function IsNonSingularPiece(r,A,Si,Eig) M := RelMatrix(r,A,Si,Eig); M0 := M; good := true; while good and (NumberOfRows(M0) gt 4) do N := KernelMatrix(M0); j := convex_row(N); if j ne 0 then RemoveRow(~M0,j); else good := false; end if; end while; if good then N0 := Transpose(Adjoint(M0)); // dual cone if &and[N0[i,j] le 0 : j in [1..4], i in [1..4]] then N0 := -1*N0; end if; Discr := [1/r*&+[N0[i,j] : j in [1..4]]-1 : i in [1..4]]; // return true, M0, N0, Discr; return true, N0, Discr; else return false, M0, KernelMatrix(M0); end if; end function; // Hirzebruch-Jung continued fraction expansion. e.g. 20/3 = [7,3] = 7 - 1/3 function HJ(x); CF := []; while x ne Ceiling(x) do Append(~CF,Ceiling(x)); x := 1/(Ceiling(x)-x); end while; Append(~CF,Ceiling(x)); return CF; end function; // // typical way of using // Sol20:= [A : A in Solutions(20) | A[1] eq 1]; ASets20 := [ASets(20,A) : A in Sol20]; for i in [1..#Sol20] do A := Sol20[i]; Eig := EigenSp(20,A); i, #ASets20[i], &and[IsNonSingularPiece(20,A,ASets20[i][j],Eig) : j in [1..#ASets20[i]]]; end for; i:=8; r:=20; A:= Sol20[i]; Eig:= EigenSp(r,A); AS := ASets20[i]; for j in [1..5] do IsNonSingularPiece(r,A,AS[j],Eig); end for; for j in [15,17] do IsNonSingularPiece(r,A,AS[j],Eig); end for; Sol30:= [A : A in Solutions(30) | A[1] eq 1]; ASets30 := [ASets(30,A) : A in Sol30]; i:=27; A := Sol30[i]; Eig := EigenSp(30,A); for S in ASets30[i] do IsNonSingularPiece(30,A,S,Eig); end for; // 27th case 1/30(1,5,8,16), 56 affine pieces, only No.11 and No.21 singular r:=15; Solr:= [A : A in Solutions(r) | A[1] eq 1]; ASetsr := [ASets(r,A) : A in Solr]; for i in [1..#Solr] do A := Solr[i]; Eig := EigenSp(r,A); i, #ASetsr[i], &and[IsNonSingularPiece(r,A,ASetsr[i][j],Eig) : j in [1..#ASetsr[i]]]; end for; First interesting case is r = 12, with 7 cases, of which 2 have typical mild discrepancy. First singular case is 1/15(1,3,5,6), 28 affines r = 20: 3 singular case: 1/20(1,2,5,12) with 42 affines. 1/20(1,3,4,12) with 33. 1/20(1,4,5,10) with 30 r = 24: 23 case, of which 15 are nonsingular and 8 singular Do r = 24, Sol24[9] = 1/24(1,2,5,16). The 4 singular affine pieces are given by > [j : j in [1..#ASetsr[i]] | not IsNonSingularPiece(r,A,ASetsr[i][j],Eig)]; [ 19, 27, 28, 42 ] Case 19: > M0; [ 0 -1 2 1] m*c [ 0 2 -4 1] m*b [ 2 -1 0 0] m*l [ 1 2 -1 0] m*b*l [ 0 2 4 0] m*b*c*l [ 1 1 1 1] m^2*b*c*l disc m = 0 (presumably (m = d = 0) + (m = l = 0)) [-1 -1 -1 2] d t^2 = d*x*y*z [ 1 -1 -3 1] m x*t=m*y*z^3; [-1 3 -1 0] b y^3=b*x*z [-1 0 5 0] c z^5 = c*x [ 1 0 3 -1] l x*z^3 = l*t > KernelMatrix(M0); [ 1 0 0 0 0 0 0 -1 0 -1 0] z^2*t = m*c*y because y*z^3 [ 0 1 0 0 0 0 0 -1 -1 0 0] y^2*t = m*b*z^4 because y*z^4 [ 0 0 1 0 0 0 0 -1 0 0 -1] x^2 = m*l*y because y*z^3 [ 0 0 0 1 0 0 0 -1 -1 0 -1] x*y^2 = m*b*l*z because z*t [ 0 0 0 0 1 0 0 -1 -1 -1 -1] y^2*z^4 = m*b*c*l because t [ 0 0 0 0 0 1 0 -2 -1 -1 -1] x*y*z*t = m^2*b*c*l because y [ 0 0 0 0 0 0 1 -1 0 -1 1] d*l = m*c because x*y*z m*c*x*y*z = m*y*z^6 = x*z^3*t = l*t^2 = l*d*x*y*z Case 27: > M0; [ 0 5 -2 0] [ 1 -1 -3 1] [ 1 2 -1 0] [ 2 -1 0 0] [ 0 -1 2 1] [ 0 2 4 0] [ 1 1 1 1] [ 0 -4 0 2] d t^2 = d*y^4 [ 0 3 2 -1] mu y^3*z^2 = mu*t [ 0 2 -4 1] la y^2*t = la*z^4 [ 1 -3 1 0] l x*z = l*y^3 [-1 0 5 0] c z^5 = c*x > KernelMatrix(M0); [ 1 0 0 0 0 0 0 0 -1 -1 0 0] y^5 = mu*la*z^2 because z^4 [ 0 1 0 0 0 0 0 0 0 -1 -1 0] x*t = la*l*y*z^3 because y*z^4 [ 0 0 1 0 0 0 0 0 -1 -1 -1 0] x*y^2 = mu*la*l*z because z*t [ 0 0 0 1 0 0 0 0 -1 -1 -2 0] x^2 = mu*la*l^2*y because y^3 [ 0 0 0 0 1 0 0 0 0 -1 -1 -1] z^2*t = la*l*c*y because x*y [ 0 0 0 0 0 1 0 0 -1 -1 -1 -1] y^2*z^4 = mu*la*l*c*y because y [ 0 0 0 0 0 0 1 0 -1 -2 -2 -1] x*y*z*t = mu*la^2*l^2*c*y not basic [ 0 0 0 0 0 0 0 1 1 -1 -1 -1] d*mu = la*l*c because x*y z^2*t = d*mu*y because y^4: d*mu*y^4 = y^3*z^2*t Case 28 > M0; [ 0 3 2 -1] [ 1 -3 1 0] [ 1 2 -1 0] [ 2 -1 0 0] [ 0 -1 2 1] [ 0 4 0 1] [ 1 1 1 1] [ 0 -4 0 2] [ 0 -2 4 -1] [ 0 5 -2 0] [ 1 -1 -3 1] [-1 2 1 1] > KernelMatrix(M0); [ 1 0 0 0 0 0 0 0 -1 -1 0 0] [ 0 1 0 0 0 0 0 0 -1 0 -1 0] [ 0 0 1 0 0 0 0 0 -1 -1 -1 0] [ 0 0 0 1 0 0 0 0 -2 -1 -2 0] [ 0 0 0 0 1 0 0 0 -1 0 -1 -1] [ 0 0 0 0 0 1 0 0 -1 -1 -1 -1] [ 0 0 0 0 0 0 1 0 -2 -1 -2 -1] [ 0 0 0 0 0 0 0 1 0 1 -1 -1] Case 42; > M0; [ 0 -1 2 1] [ 0 5 -2 0] [ 0 2 4 0] [ 1 1 1 1] [ 0 -4 0 2] t^2 = d*y^4 [ 0 3 2 -1] y^3*z^2 = mu*t [ 1 0 -5 0] x = a*z^5 [ 0 2 -4 1] y^2*t = la*z^4 [ 0 -3 6 0] z^6 = c*y^3 > KernelMatrix(M0); [ 1 0 0 0 0 0 0 -1 -1] z^2*t = la*c*y because y^3 [ 0 1 0 0 0 -1 0 -1 0] y^5 = mu*la*z^2 because z^4 [ 0 0 1 0 0 -1 0 -1 -1] y^2*z^4 = mu*la*c because z^2 [ 0 0 0 1 0 -1 -1 -2 -2] x*y*z*t = a*mu*la^2*c^2 because 1 [ 0 0 0 0 1 1 0 -1 -1] d*mu = la*c because y^4 y^3*z^2*t = mu*t^2 = d*mu*y^4 Sol30:= [A : A in Solutions(30) | A[1] eq 1]; ASets30 := [ASets(30,A) : A in Sol30]; for i in [26..30] do A := Sol30[i]; Eig := EigenSp(30,A); i, #ASets30[i], &and[IsNonSingularPiece(30,A,ASets30[i][j],Eig) : j in [1..#ASets30[i]]]; end for; i:=27; A := Sol30[i]; Eig := EigenSp(30,A); for S in ASets30[i] do IsNonSingularPiece(30,A,S,Eig); end for; 27th case 1/30(1,5,8,16), 56 affine pieces of which No.11 and No.21 are singular // // Application 1: AHilb(1/23(1,2,4,16) is nonsingular with // a single exceptional divisor of discrepancy 1 // corresponding to P9 = 1/23(9,18,13,6) // appearing on affine pieces numbers 4,5,9,10,11,12 r := 23; A := [1,2,4,16]; Eig := EigenSp(r,A); AH23 := ASets(r,A); [#S : S in AH23]; // 26 [IsNonSingularPiece(r,A,Si,Eig) : Si in AH23]; // all true for i in [1..5] do IsNonSingularPiece(r,A,AH23[i],Eig); end for; // Application 2: r := 20; Sol20:= [A : A in Solutions(20) | A[1] eq 1]; ASets20 := [ASets(r,A) : A in Sol20]; [#S : S in ASets20]; // [ 20, 24, 20, 20, 20, 20, 20, 42, 27, 33, 20, 30 ] for A in Sol20 do Eig:=EigenSp(20,A); &and[IsNonSingularPiece(20,A,S,Eig) : S in ASets(r,A)]; end for; The three singular guys are: > Sol20[8]; [ 1, 2, 5, 12 ] > Sol20[10]; [ 1, 3, 4, 12 ] > Sol20[12]; [ 1, 4, 5, 10 ] ////////////////// A := Sol20[12]; AH12 := ASets20[12]; Eig := EigenSp(20,A); [i : i in [1..#AH12] | not IsNonSingularPiece(r,A,AH12[i],Eig)]; // [18,21] // i.e. only two affine pieces out of 30 are sing > IsNonSingularPiece(r,A,AH12[18],Eig); false [ 3 -2 -3 0] [ 2 -3 2 0] [-1 4 -3 0] [-2 3 2 0] [ 0 0 -2 1] [ 1 -1 -1 1 0] > IsNonSingularPiece(r,A,AH12[21],Eig); false [ 3 -2 -1 -1] [ 2 -3 0 1] [-1 4 -1 -1] [-2 3 0 1] [ 0 0 2 -1] [ 1 -1 -1 1 0] ////////////////// A := Sol20[10]; AH10 := ASets20[10]; Eig := EigenSp(20,A); [i : i in [1..#AH10] | not IsNonSingularPiece(r,A,AH10[i],Eig)]; // [24,28] // i.e. only two affine pieces out of 33 are sing > IsNonSingularPiece(r,A,AH10[24],Eig); false [ 2 -2 -4 0] [ 1 -3 2 0] [ 0 4 -3 0] [-1 3 3 0] [ 0 0 -3 1] [ 1 -1 -1 1 0] ////////////////// AH8 := ASets20[8]; [#S : S in AH8]; /* [ 4, 5, 7, 7, 4, 4, 9, 11, 9, 10, 8, 7, 5, 4, 12, 10, 11, 9, 10, 11, 10, 7, 9, 10, 9, 5, 12, 11, 11, 7, 9, 9, 7, 5, 4, 4, 7, 4, 4, 4, 7, 5 ] */ > [IsNonSingularPiece(r,A,S,Eig) : S in AH8]; // [ true, true, .. true, false, true, false, true, .. true ] > [i : i in [1..#AH8] | not IsNonSingularPiece(r,A,AH8[i],Eig)]; // [ 15, 17 ] ////////////////// // All the calc for r = 20 r := 20; Sol20:= [A : A in Solutions(20) | A[1] eq 1]; ASets20 := [ASets(r,A) : A in Sol20]; [#S : S in ASets20]; for A in Sol20 do Eig:=EigenSp(20,A); &and[IsNonSingularPiece(20,A,S,Eig) : S in ASets(r,A)]; end for; A := Sol20[12]; AH12 := ASets20[12]; Eig := EigenSp(20,A); [i : i in [1..#AH12] | not IsNonSingularPiece(r,A,AH12[i],Eig)]; IsNonSingularPiece(r,A,AH12[18],Eig); IsNonSingularPiece(r,A,AH12[21],Eig); A := Sol20[10]; AH10 := ASets20[10]; Eig := EigenSp(20,A); [i : i in [1..#AH10] | not IsNonSingularPiece(r,A,AH10[i],Eig)]; IsNonSingularPiece(r,A,AH10[24],Eig); IsNonSingularPiece(r,A,AH10[28],Eig); A := Sol20[8]; AH8 := ASets20[8]; Eig := EigenSp(20,A); [i : i in [1..#AH10] | not IsNonSingularPiece(r,A,AH8[i],Eig)]; IsNonSingularPiece(r,A,AH8[15],Eig); IsNonSingularPiece(r,A,AH8[17],Eig); Case 1/20(1,2,5,12). Has 42 affines, of which 2 are singular, making PP^1 x (3-fold node). for i in [1..5] do IsNonSingularPiece(r,A,AH8[i],Eig); end for; discrepant divisors correspond to the vectors: P6 = 1/20(6,12,10,12) disc = 1 on 3,4,9,10,11,18,19,23,24,25 P9 = 1/20(9,18,5,8) disc = 1 on 8,10,11,12,13,18 P14 = 1/20(14,8,10,8) disc = 1 on 16,19,20,21,24,25,27,28 P1 + e1 = 1/20(21,2,5,12) disc = 1 on 29,30,32,33,39,40,41,42 P2 + e1 = 1/20(22,4,10,4) disc = 1 on 20,21,22,27,28,29,31,32,37,38,41,42 P11 + e2 = 1/20(11,22,15,12) disc = 2 on 7,8,9,10,16,19 Disc := []; for i in [1..#AH8] do nonsing, M, De := IsNonSingularPiece(r,A,AH8[i],Eig); if nonsing then Disc := Disc cat [M[i] : i in [1..4] | De[i] ne 0]; end if; end for; > DiscDiv := SetToSequence(SequenceToSet(Disc)); > [&+[A[i] : i in [1..4]] div 20 : A in DiscDiv]; // [ 2, 2, 3, 2, 2, 2 ] > DiscDiv; [ ( 6 12 10 12), ( 9 18 5 8), (11 22 15 12), (14 8 10 8), (22 4 10 4), (21 2 5 12) ] i.e. The discrepant divisors on the nonsing pieces are P6, P9, P14, P2+e1, P1 + e1 and P11+e2 (disc 2). Disc := []; for i in [1..#AH8] do nonsing, M, De := IsNonSingularPiece(r,A,AH8[i],Eig); if nonsing then Disc := Disc cat [M[i] : i in [1..4] | De[i] ne 0]; end if; end for; P6 := DiscDiv[1]; P9 := DiscDiv[2]; Q11 := DiscDiv[3]; P14 := DiscDiv[4]; Q2 := DiscDiv[5]; Q1 := DiscDiv[6]; IB1302525034 for i in [17..20] do IsNonSingularPiece(r,A,AH8[i],Eig); end for; The affine pieces involving P6 are 3,4,5,10,11,18,19,23,24,25 and 17 (singular). The star of P6 in the fan of A-Hilb is e2 over the P6-centred pentagon P4 P1 P2 Q11 P9 i.e. the 5 simplexes 4. P6 e2 P4 P1 3. P6 e2 P1 P2 9. P6 e2 P2 Q11 10. P6 e2 Q11 P9 11. P6 e2 P9 P4 and 23. P6 P1 P4 P12 25. P6 P14 P1 P2 19. P6 P14 P2 Q11 24. P6 P14 P12 P1 17. P6 P14 Q11 P12 P9 (not simplex) 18. P6 P12 P9 P4