Files
Calabi-Yau-Period-Geometry/sage/toric_topdata.sage
T

341 lines
17 KiB
Python

### CLASS: ToricTopData ####
# Methods:
# - pointstopoly: converts list of points to CY-polytope
# - disc: computes discriminant for polytope
# - topdata: computes all sorts of topological data for polytope
# - CICYtopdata: computes topdata in the simple CICY format,
# e.g. [[3,3]] for two cubics in P5 or [[3,0,1],[0,3,1]] for the Tian--Yau manifold
import numpy as np
# For grep-commands
import re
class ToricTopData:
findzs = re.compile('z\d+')
findls = re.compile('l\d+')
findindexz = re.compile('\d+')
def pointstopoly(self, points):
# Takes a list of points and returns the points inside its convex hull while omitting points inside faces of co-dimension one
zeros = [0 for i in range(len(points[0]))] #origin
pc = LatticePolytope(points) # all points
polytope = [list(m) for m in pc.points()]
self.pcodim1 = [] #points in co-dimension one
for f in pc.facets():
for i in f.interior_point_indices():
try:
polytope.remove([list(m) for m in f.points(i)][0]) # remove those inside codim1 faces
self.pcodim1.append([list(m) for m in f.points(i)][0]) # save omitted points in pcodim1
except:
pass
# Moving zeros to the end
polytope.remove(zeros)
polytope.append(zeros)
return polytope
def disc(self, polytope, only_sc=False, no_triangulation=0):
# Takes a poltyope and returns a list of tuples [disc_i,codim_i] of discriminant factors disc_i coming from a relation inside a face of codimension codim_i
self.pc = LatticePolytope(polytope)
dimp = len(polytope[0])
zeros = [0 for i in range(dimp)]
# Compute all relations in all faces
self.globalks = []
if only_sc:
pcfaces = self.pc.faces()[:-1]
else:
pcfaces = self.pc.faces()
for facesofdim in pcfaces: # for each face-dimension
atdim = []
for face in facesofdim: # for each face of fixed dimension
atface = []
facepoints = [list(i) for i in matrix(face.points())]
try:
facepoints.remove(zeros)
except:
pass
for p in facepoints.copy():
try:
if p not in polytope:
facepoints.remove(p) # remove points that were omitted for polytope
except:
pass
rels = {str(i):polytope.index(i) for i in facepoints} # dictionary for getting positions in relations correct
kernel = matrix(facepoints).kernel().gens()
for k in kernel: # just formatting
globalk = [0 for i in range(len(polytope)-1)] # MAY HAVE TO BE ADAPTED FOR CASES WITH MORE VERTICES
for i in range(len(k)):
globalk[rels[str(facepoints[i])]] = k[i]
atdim.append(globalk)
self.globalks.append([[i[j] for j in range(len(i))] for i in np.unique(np.matrix(atdim),axis=0)])
for atdim in self.globalks: # inserting weight for innter point
for k in atdim:
if len(k)>0:
k.append(-sum(k))
# Computing discriminants
self.pc = PointConfiguration(polytope)
pc_star = self.pc.restrict_to_star_triangulations(zeros)
pc_star_and_fine=pc_star.restrict_to_fine_triangulations()
if(len(pc_star_and_fine.triangulations_list())==0):
print("No fine star triangulation!")
return 0
if(len(pc_star_and_fine.triangulations_list())>1):
print("More than one fine star triangulation! (",len(pc_star_and_fine.triangulations_list()),")")
triangulation=pc_star_and_fine.triangulations_list()[no_triangulation]
fan=triangulation.fan(zeros)
X=ToricVariety(fan)
ls = [l for l in X.Mori_cone().rays()]
if (len(matrix(ls).kernel().gens())>0):
print("Non-simplicical Mori-cone")
As = [var("a_{}".format(u), latex_name="a_{{}}".format(u)) for u in (1..len(ls))]
z = var('z', n=len(ls)+1, latex_name='z') # z[0] is superfluous
l = var('l', n=len(ls)+1, latex_name='l') # l[0] is superfluous
#self.globalks.reverse()
self.globalks = [k for k in self.globalks if k not in ([[]],)]
discs = []
for i,atdim in enumerate(self.globalks):
solsys = []
for kindex in range(len(atdim)):
sol=solve((matrix(As)*matrix(ls)-matrix(atdim[kindex])).list(),matrix(As).list())
zindex = [abs(i.subs(sol)) for i in As].index(1) +1
solsys.append(z[zindex]-prod([sum([atdim[j][i]*l[j] for j in range(len(atdim))])**(atdim[kindex][i]) for i in range(len(atdim[0]))]))
lambdas = list(set(self.findls.findall(str(solsys))))
lambdas.sort()
try:
pol = maxima.eliminate(solsys,[eval(i) for i in lambdas[1:]]).sage()
if pol[0]==0:
pol = maxima.eliminate(solsys,[eval(i) for i in lambdas[:-1]]).sage()
reversesolve = True
else:
reversesolve = False
except: # this is just for 1-parameter cases where there is nothing to solve
pol = solsys
pass
try:
if not reversesolve:
pol = [i/l0**(i.degree(l0)) for i in pol]
else:
pol = [i/eval(lambdas[-1])**(i.degree(eval(lambdas[-1]))) for i in pol]
except:
pass
for poli in pol:
if (poli not in [a[0] for a in discs] and poli!=0):
if only_sc:
discs.append([poli,len(self.globalks)-i]) # account for offset
else:
discs.append([poli,len(self.globalks)-1-i])
if only_sc: #insert an empty list for the codimension 0 face discriminant
discstmp = discs
discstmp.reverse()
discstmp.append([])
discstmp.reverse()
discs = discstmp
return discs
def topdata(self, polytope, nef_partition=0, lvec=0, no_triangulation=0, returnvars=False):
# Computes topological data for hypersurfaces and CICYs in toric ambient spaces.
# For hypersurfaces, nef_partition should remain untouched (=0).
# For CICYs with d polynomials, a nef_partition has to be supplied:
# its format should be a list of d lists as in the examples below
# which is a decomposition of the N points of $(polytope); the number
# should corresponds to an enumeration of the points of $(polytope) with
# the inner point omitted.
# Conditions for possible fibrations are:
# - elliptic: if t_i^n = 0 and t_i^{n-1} != 0
# - K3 (only 3-folds): if c2.t_i = 24
# input: - array of points in polytope
# - (for CICYs:) nef-partition
# - (optional:) set of l-vectors to be used
# - (optional:) index of triangulation
# output: - Mori-cone generators
# - Intersection rings of CY and ambient space
# - Topological data
# - Possible fibrations
# - instance attributes: triangulation and Mori-cone generators
pc = PointConfiguration(polytope)
dimp=len(polytope[0]) # dimension of polytope
if( nef_partition == 0):
no_polys = 1
nef_partition = [[i for i in range(len(polytope)-1)]] # trivial partition
else:
if( sorted(flatten(nef_partition)) != [ i for i in range(len(polytope)-1)] ):
print("Nef-partition not valid!")
no_polys = len(nef_partition)
zeros=zero_vector(dimp)
pc_star = pc.restrict_to_star_triangulations(zeros)
pc_star_and_fine=pc_star.restrict_to_fine_triangulations()
if(len(pc_star_and_fine.triangulations_list())==0):
print("No fine star triangulation!")
return 0
if(len(pc_star_and_fine.triangulations_list())>1):
print("More than one fine star triangulation! (",len(pc_star_and_fine.triangulations_list()),")")
self.triangulation=pc_star_and_fine.triangulations_list()[no_triangulation]
self.fan=self.triangulation.fan(zeros)
self.X=ToricVariety(self.fan)
# Change l-vectors if given:
if lvec==0:
self.lall=matrix(self.X.Mori_cone().rays())
else:
self.lall=lvec
nodivs = len(polytope)-dimp-1
self.MoriMatrix=matrix([[-sum([self.lall[j][i] for i in part]) for part in nef_partition] for j in range(nodivs)]), self.lall[:,:-1]
# Check whether Mori-cone is simplicial and continue with first $(nodivs) vectors
if ( self.lall.dimensions()[0] > nodivs ):
print("Non-simplicial Kähler-cone! (",self.lall.dimensions()[0]," > ",nodivs,")")
print("Picking ",nodivs," linearly independent vectors.")
l = transpose(transpose(self.lall)*(matrix(transpose(matrix(self.lall.kernel().gens())).kernel().gens()).transpose()))
self.MoriMatrix=matrix([[-sum([l[j][i] for i in part]) for part in nef_partition] for j in range(nodivs)]), l[:,:-1]
else:
l = self.lall
HH=self.X.cohomology_ring()
self.D = [HH(self.X.divisor(i)) for i in [0..len(polytope)-2]]
zs = list(set(self.findzs.findall(str(self.D)))) # gives the $(nodivs) z-variables present in divisor classes
# Find Kähler-cone generators as duals to l-vectors:
Bsinv=l.matrix_from_columns([eval(m) for m in sorted(list(self.findindexz.findall(str(zs))))])
if (Bsinv.det()==0):
print("l-vectors not independent in divisor basis. Pick them manually with ``lvec=...''")
return
Bs = transpose(Bsinv**(-1))
self.J = Bs*vector([self.D[i] for i in list(set([eval(i) for i in self.findindexz.findall(str(self.findzs.findall(str(self.D))))]))])
# Hyperplane divisor class:
self.Ydualform = product([sum([HH(self.D[p]) for p in part]) for part in nef_partition])
# Tuple lists for computations below
TupleListCY = UnorderedTuples(range(nodivs),dimp-no_polys)
TupleListCYordered = Tuples(range(nodivs),dimp-no_polys) # for intersection ring multiplicities must be included
TupleListAm = UnorderedTuples(range(nodivs),dimp)
# Find intersection numbers on CY and on ambient space:
intCY = [list((tl,self.X.integrate(self.Ydualform*product([self.J[tl[j]] for j in range(dimp-no_polys)])))) for tl in TupleListCY]
intAm = [list((tl,self.X.integrate(product([self.J[tl[j]] for j in range(dimp)])))) for tl in TupleListAm]
# Form intersection ring on CY:
t = var('t', n=nodivs, latex_name='t')
intring = sum([product([t[index] for index in tl])*self.X.integrate(self.Ydualform*product([self.J[tl[j]] for j in range(dimp-no_polys)])) for tl in TupleListCYordered])
intringnomults = sum([product([t[index] for index in tl])*self.X.integrate(self.Ydualform*product([self.J[tl[j]] for j in range(dimp-no_polys)])) for tl in set([tuple(sorted(i)) for i in TupleListCYordered])])
# Compute Chern-character for topological data:
var(self.findzs.findall(str([self.J[i] for i in range(nodivs)]))) # Introduces all z in Kähler-forms as variables...
zs = [eval(z) for z in self.findzs.findall(str([self.J[i] for i in range(nodivs)]))] # ... and puts them into a vector.
tsubs=solve([lift(self.J[i])==t[i] for i in range(nodivs)],zs) # Finds substitution rule for zs in terms Kähler-cone generators
cc = (product([1+eval(str((lift(d)))).subs(tsubs[0]) for d in self.D]))/(product([1+sum([eval(str(lift(HH(self.D[p])))) for p in part]) for part in nef_partition]).subs(tsubs[0])) # Adjunction-formula
print('--- Toric divisors (ambient space) -----------')
print(self.D)
print('\n--- Kähler cone generators (ambient space) --- ')
print(self.J)
print('\n--- Mori-cone-generators (ambient space) ------')
for j in range(nodivs):
print(self.MoriMatrix[0][j],self.MoriMatrix[1][j])
print('\n--- Intersection on CY ------------------------')
Ccy=[product([t[iden[0][i]] for i in range(dimp-no_polys)])==iden[1] for iden in intCY]
print(Ccy)
print('\nR = ',intring)
print('\nR (no multiplicities) = ',intringnomults)
print('\n--- Intersection in ambient space -------------')
Cam=[product([t[iden[0][i]] for i in range(dimp)])==iden[1] for iden in intAm]
print(Cam)
print('\n--- Topological data --------------------------')
#print((product([1+lift(d) for d in self.D]))/(product([1+sum([lift(HH(self.D[p])) for p in part]) for part in nef_partition])))
chi = self.X.integrate(self.Ydualform*(product([1+lift(d) for d in self.D]))/(product([1+sum([lift(HH(self.D[p])) for p in part]) for part in nef_partition])))
print('chi = ',chi)
c = var('c')
if( round((cc).subs({t:c*t for t in t}).taylor(c,0,1).coefficient(c,1).subs({t[i]:1 for i in range(nodivs)}),8) != 0):
return "First Chern-class not zero!"
# Print Chern classes:
cJ = [[0]]*(dimp-no_polys)
for i in range(2,dimp-no_polys):
print('\nc',i,' = ',(cc).subs({t:c*t for t in t}).taylor(c,0,i).coefficient(c,i))
TupleListi = UnorderedTuples(range(nodivs),i)
TupleListz = UnorderedTuples(range(nodivs),dimp-no_polys-i)
cJ[i] = [sum([((cc).subs({t:c*t for t in t}).taylor(c,0,i).coefficient(c,i)).coefficient(product([t[i] for i in tupl]))*((product([t[j] for j in tuplz])*product([t[i] for i in tupl])).subs(Ccy)) for tupl in TupleListi]) for tuplz in TupleListz]
print('integrated: '+str([cJ[i][j]*product([t[k] for k in TupleListz[j]]) for j in range(len(TupleListz))]))
print('\nc',dimp-no_polys,' = ',(cc).subs({t:c*t for t in t}).taylor(c,0,dimp-no_polys).coefficient(c,dimp-no_polys))
msg = ""
# elliptic:
## Indicates whether J_i^n==0 with J_i^(n-1)!=0 for n-folds and some J_i
for i in range(nodivs):
for intn in intCY:
if( intn[0] == [i]*(dimp-no_polys) and intn[1] == 0):
if not(all([round((t[i]**(dimp-no_polys-1)*t[j]).subs(Ccy),10) == 0 for j in range(nodivs) if j != i]) ):
msg += "Possible elliptic fibration in cycle dual to t"+str(i)+".\n"
# K3 for 3-folds
## Indicates whether c2.J_i =24 for some J_i
if( dimp+1-no_polys == 4):
TupleList2 = UnorderedTuples(range(nodivs),2)
for j in range(nodivs):
if( sum([((cc).subs({t:c*t for t in t}).taylor(c,0,2).coefficient(c,2)).coefficient(product([t[i] for i in tupl]))*((product([t[i] for i in tupl])*t[j]).subs(Ccy)) for tupl in TupleList2]) == 24 ):
msg += "Possible K3-fibration in divisor t"+str(j)+".\n"
if( msg != "" ):
print('\n--- Fibrations --------------------------------')
print(msg)
if returnvars:
return str([(self.MoriMatrix[0][j],self.MoriMatrix[1][j]) for j in range(nodivs)]).replace("(","{").replace(")","}").replace("[","{").replace("]","}") # might return more in the future if necessary
def CICYtopdata(self, CICY, justdata=False, justmori=False, lvec=0):
# Input: a list l with entries l[i,j] that give the weight in the
# ambient projective space i of polynomial j.
# E.g. (P^3| 3 1)
# (P^2| 2 0)
# corresponds to the list ((3,1),(2,0))
ambients=(np.array(CICY).transpose()).sum(axis=0)-1
if ambients in ZZ:
ambients=[ambients]
dimp=sum(ambients)
polytope=identity_matrix(int(dimp))
zeros=[0]*dimp
polytope=polytope.insert_row(0,zeros)
offset=0
for i in range(len(ambients)):
newrow=[0]*dimp
entry=[]
for j in range(ambients[i]):
newrow[offset]=-1
offset+=1
entry.append(offset)
polytope=polytope.insert_row(sum([ambients[n]+1 for n in [0..i]]),newrow)
partition=[]
eqs=matrix(CICY).transpose()
counters=[0]*len(ambients)
for eq in eqs:
part=[]
for j in range(len(eq)):
part.append([sum([ambients[m]+1 for m in [0..(j-1)]])+counters[j]+a for a in [0..(eq[j]-1)]])
counters[j]+=len(part[j])
partition.append(flatten(part))
if justdata:
return [list(polytope),partition]
elif justmori:
moricone((list(polytope)))
else:
self.topdata(list(polytope),nef_partition=partition,lvec=lvec)