From cfb12ad960d73e9b7a1e3b890438a571d183cc1b Mon Sep 17 00:00:00 2001 From: Julian Piribauer Date: Fri, 17 Jul 2026 22:19:09 +0200 Subject: [PATCH] Initial commit of topdata --- playground/.gitignore | 4 + sage/toric_topdata.sage | 340 ++++++++++++++++++++++++++++++++++++ sage/toric_topdata.sage.py | 348 +++++++++++++++++++++++++++++++++++++ 3 files changed, 692 insertions(+) create mode 100644 playground/.gitignore create mode 100644 sage/toric_topdata.sage create mode 100644 sage/toric_topdata.sage.py diff --git a/playground/.gitignore b/playground/.gitignore new file mode 100644 index 0000000..15c3deb --- /dev/null +++ b/playground/.gitignore @@ -0,0 +1,4 @@ +# Exclude testing files + +* +!.gitignore \ No newline at end of file diff --git a/sage/toric_topdata.sage b/sage/toric_topdata.sage new file mode 100644 index 0000000..9ee0b96 --- /dev/null +++ b/sage/toric_topdata.sage @@ -0,0 +1,340 @@ +### 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) diff --git a/sage/toric_topdata.sage.py b/sage/toric_topdata.sage.py new file mode 100644 index 0000000..e39c79c --- /dev/null +++ b/sage/toric_topdata.sage.py @@ -0,0 +1,348 @@ + + +# This file was *autogenerated* from the file sage/toric_topdata.sage +from sage.all_cmdline import * # import sage library + +_sage_const_0 = Integer(0); _sage_const_1 = Integer(1); _sage_const_2 = Integer(2); _sage_const_8 = Integer(8); _sage_const_10 = Integer(10); _sage_const_4 = Integer(4); _sage_const_24 = Integer(24)### FUNCTIONS: #### +# - 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 +findzs = re.compile('z\d+') +findls = re.compile('l\d+') +findindexz = re.compile('\d+') + + +def pointstopoly(points): + # Takes a list of points and returns the points inside its convex hull while omitting points inside faces of co-dimension one + + global pcodim1 #points in co-dimension one + + zeros = [_sage_const_0 for i in range(len(points[_sage_const_0 ]))] #origin + pc = LatticePolytope(points) # all points + polytope = [list(m) for m in pc.points()] + pcodim1 = [] + for f in pc.facets(): + for i in f.interior_point_indices(): + try: + polytope.remove([list(m) for m in f.points(i)][_sage_const_0 ]) # remove those inside codim1 faces + pcodim1.append([list(m) for m in f.points(i)][_sage_const_0 ]) # save omitted points in pcodim1 + except: + pass + # Moving zeros to the end + polytope.remove(zeros) + polytope.append(zeros) + return polytope + +def disc(polytope, only_sc=False, no_triangulation=_sage_const_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 + + global globalks, pc + pc = LatticePolytope(polytope) + dimp = len(polytope[_sage_const_0 ]) + zeros = [_sage_const_0 for i in range(dimp)] + # Compute all relations in all faces + globalks = [] + if only_sc: + pcfaces = pc.faces()[:-_sage_const_1 ] + else: + pcfaces = 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 = [_sage_const_0 for i in range(len(polytope)-_sage_const_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) + globalks.append([[i[j] for j in range(len(i))] for i in np.unique(np.matrix(atdim),axis=_sage_const_0 )]) + for atdim in globalks: # inserting weight for innter point + for k in atdim: + if len(k)>_sage_const_0 : + k.append(-sum(k)) + # Computing discriminants + pc = PointConfiguration(polytope) + 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())==_sage_const_0 ): + print("No fine star triangulation!") + return _sage_const_0 + if(len(pc_star_and_fine.triangulations_list())>_sage_const_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())>_sage_const_0 ): + print("Non-simplicical Mori-cone") + + As = [var("a_{}".format(u), latex_name="a_{{}}".format(u)) for u in (ellipsis_iter(_sage_const_1 ,Ellipsis,len(ls)))] + z = var('z', n=len(ls)+_sage_const_1 , latex_name='z') # z[0] is superfluous + l = var('l', n=len(ls)+_sage_const_1 , latex_name='l') # l[0] is superfluous + #globalks.reverse() + globalks = [k for k in globalks if k not in ([[]],)] + discs = [] + for i,atdim in enumerate(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(_sage_const_1 ) +_sage_const_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[_sage_const_0 ]))])) + lambdas = list(set(findls.findall(str(solsys)))) + lambdas.sort() + try: + pol = maxima.eliminate(solsys,[eval(i) for i in lambdas[_sage_const_1 :]]).sage() + if pol[_sage_const_0 ]==_sage_const_0 : + pol = maxima.eliminate(solsys,[eval(i) for i in lambdas[:-_sage_const_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[-_sage_const_1 ])**(i.degree(eval(lambdas[-_sage_const_1 ]))) for i in pol] + except: + pass + for poli in pol: + if (poli not in [a[_sage_const_0 ] for a in discs] and poli!=_sage_const_0 ): + if only_sc: + discs.append([poli,len(globalks)-i]) # account for offset + else: + discs.append([poli,len(globalks)-_sage_const_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(polytope, nef_partition=_sage_const_0 , lvec=_sage_const_0 , no_triangulation=_sage_const_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 + # - GLOBALS: triangulation and Mori-cone generators + + global triangulation, lall, MoriMatrix, fan, X, Ydualform, J, D + + pc = PointConfiguration(polytope) + + dimp=len(polytope[_sage_const_0 ]) # dimension of polytope + if( nef_partition == _sage_const_0 ): + no_polys = _sage_const_1 + nef_partition = [[i for i in range(len(polytope)-_sage_const_1 )]] # trivial partition + else: + if( sorted(flatten(nef_partition)) != [ i for i in range(len(polytope)-_sage_const_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())==_sage_const_0 ): + print("No fine star triangulation!") + return _sage_const_0 + if(len(pc_star_and_fine.triangulations_list())>_sage_const_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) + + # Change l-vectors if given: + if lvec==_sage_const_0 : + lall=matrix(X.Mori_cone().rays()) + else: + lall=lvec + nodivs = len(polytope)-dimp-_sage_const_1 + MoriMatrix=matrix([[-sum([lall[j][i] for i in part]) for part in nef_partition] for j in range(nodivs)]), lall[:,:-_sage_const_1 ] + # Check whether Mori-cone is simplicial and continue with first $(nodivs) vectors + if ( lall.dimensions()[_sage_const_0 ] > nodivs ): + print("Non-simplicial Kähler-cone! (",lall.dimensions()[_sage_const_0 ]," > ",nodivs,")") + print("Picking ",nodivs," linearly independent vectors.") + l = transpose(transpose(lall)*(matrix(transpose(matrix(lall.kernel().gens())).kernel().gens()).transpose())) + MoriMatrix=matrix([[-sum([l[j][i] for i in part]) for part in nef_partition] for j in range(nodivs)]), l[:,:-_sage_const_1 ] + else: + l = lall + + HH=X.cohomology_ring() + D = [HH(X.divisor(i)) for i in (ellipsis_range(_sage_const_0 ,Ellipsis,len(polytope)-_sage_const_2 ))] + zs = list(set(findzs.findall(str(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(findindexz.findall(str(zs))))]) + if (Bsinv.det()==_sage_const_0 ): + print("l-vectors not independent in divisor basis. Pick them manually with ``lvec=...''") + return + Bs = transpose(Bsinv**(-_sage_const_1 )) + J = Bs*vector([D[i] for i in list(set([eval(i) for i in findindexz.findall(str(findzs.findall(str(D))))]))]) + # Hyperplane divisor class: + Ydualform = product([sum([HH(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,X.integrate(Ydualform*product([J[tl[j]] for j in range(dimp-no_polys)])))) for tl in TupleListCY] + intAm = [list((tl,X.integrate(product([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])*X.integrate(Ydualform*product([J[tl[j]] for j in range(dimp-no_polys)])) for tl in TupleListCYordered]) + intringnomults = sum([product([t[index] for index in tl])*X.integrate(Ydualform*product([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(findzs.findall(str([J[i] for i in range(nodivs)]))) # Introduces all z in Kähler-forms as variables... + zs = [eval(z) for z in findzs.findall(str([J[i] for i in range(nodivs)]))] # ... and puts them into a vector. + tsubs=solve([lift(J[i])==t[i] for i in range(nodivs)],zs) # Finds substitution rule for zs in terms Kähler-cone generators + cc = (product([_sage_const_1 +eval(str((lift(d)))).subs(tsubs[_sage_const_0 ]) for d in D]))/(product([_sage_const_1 +sum([eval(str(lift(HH(D[p])))) for p in part]) for part in nef_partition]).subs(tsubs[_sage_const_0 ])) # Adjunction-formula + + print('--- Toric divisors (ambient space) -----------') + print(D) + print('\n--- Kähler cone generators (ambient space) --- ') + print(J) + print('\n--- Mori-cone-generators (ambient space) ------') + for j in range(nodivs): + print(MoriMatrix[_sage_const_0 ][j],MoriMatrix[_sage_const_1 ][j]) + + print('\n--- Intersection on CY ------------------------') + Ccy=[product([t[iden[_sage_const_0 ][i]] for i in range(dimp-no_polys)])==iden[_sage_const_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[_sage_const_0 ][i]] for i in range(dimp)])==iden[_sage_const_1 ] for iden in intAm] + print(Cam) + + print('\n--- Topological data --------------------------') + #print((product([1+lift(d) for d in D]))/(product([1+sum([lift(HH(D[p])) for p in part]) for part in nef_partition]))) + chi = X.integrate(Ydualform*(product([_sage_const_1 +lift(d) for d in D]))/(product([_sage_const_1 +sum([lift(HH(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,_sage_const_0 ,_sage_const_1 ).coefficient(c,_sage_const_1 ).subs({t[i]:_sage_const_1 for i in range(nodivs)}),_sage_const_8 ) != _sage_const_0 ): + return "First Chern-class not zero!" + + # Print Chern classes: + cJ = [[_sage_const_0 ]]*(dimp-no_polys) + for i in range(_sage_const_2 ,dimp-no_polys): + print('\nc',i,' = ',(cc).subs({t:c*t for t in t}).taylor(c,_sage_const_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,_sage_const_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,_sage_const_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[_sage_const_0 ] == [i]*(dimp-no_polys) and intn[_sage_const_1 ] == _sage_const_0 ): + if not(all([round((t[i]**(dimp-no_polys-_sage_const_1 )*t[j]).subs(Ccy),_sage_const_10 ) == _sage_const_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+_sage_const_1 -no_polys == _sage_const_4 ): + TupleList2 = UnorderedTuples(range(nodivs),_sage_const_2 ) + for j in range(nodivs): + if( sum([((cc).subs({t:c*t for t in t}).taylor(c,_sage_const_0 ,_sage_const_2 ).coefficient(c,_sage_const_2 )).coefficient(product([t[i] for i in tupl]))*((product([t[i] for i in tupl])*t[j]).subs(Ccy)) for tupl in TupleList2]) == _sage_const_24 ): + msg += "Possible K3-fibration in divisor t"+str(j)+".\n" + + if( msg != "" ): + print('\n--- Fibrations --------------------------------') + print(msg) + if returnvars: + return str([(MoriMatrix[_sage_const_0 ][j],MoriMatrix[_sage_const_1 ][j]) for j in range(nodivs)]).replace("(","{").replace(")","}").replace("[","{").replace("]","}") # might return more in the future if necessary + + +def CICYtopdata(CICY,justdata=False,justmori=False, lvec=_sage_const_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=_sage_const_0 )-_sage_const_1 + if ambients in ZZ: + ambients=[ambients] + dimp=sum(ambients) + polytope=identity_matrix(int(dimp)) + zeros=[_sage_const_0 ]*dimp + polytope=polytope.insert_row(_sage_const_0 ,zeros) + + offset=_sage_const_0 + for i in range(len(ambients)): + newrow=[_sage_const_0 ]*dimp + entry=[] + for j in range(ambients[i]): + newrow[offset]=-_sage_const_1 + offset+=_sage_const_1 + entry.append(offset) + polytope=polytope.insert_row(sum([ambients[n]+_sage_const_1 for n in (ellipsis_range(_sage_const_0 ,Ellipsis,i))]),newrow) + + partition=[] + eqs=matrix(CICY).transpose() + counters=[_sage_const_0 ]*len(ambients) + for eq in eqs: + part=[] + for j in range(len(eq)): + part.append([sum([ambients[m]+_sage_const_1 for m in (ellipsis_range(_sage_const_0 ,Ellipsis,(j-_sage_const_1 )))])+counters[j]+a for a in (ellipsis_range(_sage_const_0 ,Ellipsis,(eq[j]-_sage_const_1 )))]) + counters[j]+=len(part[j]) + partition.append(flatten(part)) + + if justdata: + return [list(polytope),partition] + elif justmori: + moricone((list(polytope))) + else: + topdata(list(polytope),nef_partition=partition,lvec=lvec) + +topdata(pointstopoly([(_sage_const_1 ,_sage_const_0 ),(_sage_const_0 ,_sage_const_1 ),(-_sage_const_1 ,-_sage_const_1 )])) +