Update 'aplK.py'

This commit is contained in:
calquezar 2021-04-26 20:14:59 +02:00
parent d014315b01
commit 21d24227b3

652
aplK.py
View File

@ -1,71 +1,105 @@
reset()
r"""
<Short one-line summary that ends with no period>
<Paragraph description>
EXAMPLES::
<Lots and lots of examples>
AUTHORS:
- YOUR NAME (2005-01-03): initial version
- person (date in ISO year-month-day format): short desc
"""
# ****************************************************************************
# Copyright (C) 2021 YOUR NAME <your email>
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 2 of the License, or
# (at your option) any later version.
# https://www.gnu.org/licenses/
# ****************************************************************************
from sage.structure.parent import Parent
from sage.structure.element import Element
from sage.structure.unique_representation import UniqueRepresentation
from sage.categories.rings import Rings
from sage.structure.richcmp import rich_to_bool, op_EQ, op_NE, op_LT, op_LE, op_GT, op_GE
from sage.categories.morphism import Morphism
from sage.rings.rational_field import QQ
from sage.algebras.commutative_dga import GradedCommutativeAlgebra
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
from sage.functions.other import factorial
#####################################################################
class Apln(UniqueRepresentation):
def __init__(self,n):
self.n = n
ts = ['t'+str(t) for t in range(n+1)]
ys = ['y'+str(t) for t in range(n+1)]
degrees = tuple((n+1)*[(0,2)])+tuple((n+1)*[(1,2)])
self.gca = GradedCommutativeAlgebra(QQ, names=ts+ys, degrees=degrees)
t_relation = -1
y_relation = 0
#t_relation = sum(self.gca.gens()[:n])-1
#y_relation = sum(self.gca.gens()[n:])
for i in range(n+1):
t_relation += self.gca.gens()[i]
y_relation += self.gca.gens()[n+1+i]
self.gca_with_relations = self.gca.quotient(self.gca.ideal(t_relation,y_relation),check=False)
dfs = {self.gca.gens_dict()[ts[i]]:self.gca.gens_dict()[ys[i]] for i in range(n+1)}
self.cdga = self.gca.cdg_algebra(dfs)
def _repr_(self):
return "Algebra Apln of dimension %s" % self.n
def inject_variables(self):
self.cdga.inject_variables()
def generators(self):
return self.cdga.gens()
def reduce_with_relations(self,pol):
return self.gca_with_relations(pol)
def homogeneousPart(self,pol,degree):
parts = [pol.homogeneous_parts()[k] for k in pol.homogeneous_parts() if list(k)[0]==degree]
hompart = 0
for p in parts:
hompart += p
return hompart
def face(self,i):
if(self.n==0):
return 0
domain = self.cdga
codomain = self.cdga
tImages = []
yImages = []
n = len(domain.gens())//2
for k in range(n):
if(k<i):
tImages.append(codomain.gens()[k])
yImages.append(codomain.gens()[n+k])
elif(k==i):
tImages.append(0)
yImages.append(0)
elif(k>i):
tImages.append(codomain.gens()[k-1])
yImages.append(codomain.gens()[n+k-1])
images = tImages+yImages
return domain.hom(codomain=codomain,im_gens=images)
def __init__(self,n):
self.n = n
ts = ['t'+str(t) for t in range(n+1)]
ys = ['y'+str(t) for t in range(n+1)]
degrees = tuple((n+1)*[(0,2)])+tuple((n+1)*[(1,2)])
self.gca = GradedCommutativeAlgebra(QQ, names=ts+ys, degrees=degrees)
t_relation = -1
y_relation = 0
#t_relation = sum(self.gca.gens()[:n])-1
#y_relation = sum(self.gca.gens()[n:])
for i in range(n+1):
t_relation += self.gca.gens()[i]
y_relation += self.gca.gens()[n+1+i]
self.gca_with_relations = self.gca.quotient(self.gca.ideal(t_relation,y_relation),check=False)
dfs = {self.gca.gens_dict()[ts[i]]:self.gca.gens_dict()[ys[i]] for i in range(n+1)}
self.cdga = self.gca.cdg_algebra(dfs)
#self.inject_variables()
def _repr_(self):
return "Algebra Apln of dimension %s" % self.n
def inject_variables(self):
self.cdga.inject_variables()
def _reduce_polynomial(self,p,n):
def generators(self):
return self.cdga.gens()
def reduce_with_relations(self,pol):
return self.gca_with_relations(pol)
def homogeneousPart(self,pol,degree):
parts = [pol.homogeneous_parts()[k] for k in pol.homogeneous_parts() if list(k)[0]==degree]
hompart = 0
for p in parts:
hompart += p
return hompart
def face(self,i):
if(self.n==0):
return 0
domain = self.cdga
codomain = self.cdga
tImages = []
yImages = []
n = len(domain.gens())//2
for k in range(n):
if(k<i):
tImages.append(codomain.gens()[k])
yImages.append(codomain.gens()[n+k])
elif(k==i):
tImages.append(0)
yImages.append(0)
elif(k>i):
tImages.append(codomain.gens()[k-1])
yImages.append(codomain.gens()[n+k-1])
images = tImages+yImages
return domain.hom(codomain=codomain,im_gens=images)
def _reduce_polynomial(self,p,n):
## Input:
# - p: polynomial in variables t_0,..t_n,y_0,...,y_n
# - n: dimension of the algebra (apl)_n where p is defined
@ -73,25 +107,25 @@ class Apln(UniqueRepresentation):
# - The polynomial after applying the corresponding relations to reduce the variables t_n, y_n
# - For practical reasons, we change the nomenclature in this way: t_i ==> xs_i, y_i ==> zs_i
xs = ['x'+str(x) for x in range(n+1)]
zs = ['z'+str(z) for z in range(n+1)]
R = PolynomialRing(QQ,names=xs+zs)
r_gens = R.gens()
gca_gens = self.gca.gens()
pol = R(0)
for monomial_exponents in list(p.dict().keys()):
monomial = R(1)
coef = p.dict()[monomial_exponents]
for i in range(n+1):
monomial *= r_gens[i]**monomial_exponents[i]
monomial *= r_gens[n+1+i]**monomial_exponents[self.n+1+i]
pol += coef*monomial
xn_relation = 1-sum(r_gens[:n])
zn_relation = -sum(r_gens[n+1:-1])
pol = pol.subs({r_gens[n]:xn_relation,r_gens[2*n+1]:zn_relation})
return pol
xs = ['x'+str(x) for x in range(n+1)]
zs = ['z'+str(z) for z in range(n+1)]
R = PolynomialRing(QQ,names=xs+zs)
r_gens = R.gens()
gca_gens = self.gca.gens()
pol = R(0)
for monomial_exponents in list(p.dict().keys()):
monomial = R(1)
coef = p.dict()[monomial_exponents]
for i in range(n+1):
monomial *= r_gens[i]**monomial_exponents[i]
monomial *= r_gens[n+1+i]**monomial_exponents[self.n+1+i]
pol += coef*monomial
xn_relation = 1-sum(r_gens[:n])
zn_relation = -sum(r_gens[n+1:-1])
pol = pol.subs({r_gens[n]:xn_relation,r_gens[2*n+1]:zn_relation})
return pol
def _convert_integrable_polynomial(self,p,n):
def _convert_integrable_polynomial(self,p,n):
## Input:
# - p: polynomial in variables t_0,..t_n,y_0,...,y_n
# - n: dimension of the algebra (apl)_n where p is defined
@ -100,284 +134,260 @@ class Apln(UniqueRepresentation):
# (After reducing the variables t_{n},y_{n} by the corresponding relations,
# the algorithm removes the variables y_{0},..,y_{n-1}).
xs = ['x'+str(x) for x in range(n)]
zs = ['z'+str(z) for z in range(n)]
R = PolynomialRing(QQ,names=xs+zs)
pol = R(self._reduce_polynomial(p,n))
r_gens = R.gens()
relations = {}
for i in range(n):
relations[r_gens[n+i]]=1
return pol.subs(relations)
xs = ['x'+str(x) for x in range(n)]
zs = ['z'+str(z) for z in range(n)]
R = PolynomialRing(QQ,names=xs+zs)
pol = R(self._reduce_polynomial(p,n))
r_gens = R.gens()
relations = {}
for i in range(n):
relations[r_gens[n+i]]=1
return pol.subs(relations)
def integral(self,p,n):
def integral(self,p,n):
## Input:
# - p: polynomial in variables t_0,..t_n,y_0,...,y_n
# - n: dimension of the algebra (apl)_n where p is defined
## Output:
# - the integral of the polynomial after reducing the variables t_n, y_n by the corresponding relations
#
xs = ['x'+str(x) for x in range(n)]
R = PolynomialRing(QQ,names=xs)
r_gens = R.gens()
pol = R(self._convert_integrable_polynomial(p,n))
result = pol
for i in range(n):
result = result.integral(r_gens[n-1-i])
result = result.subs({r_gens[n-1-i]:1-sum(r_gens[:n-1-i])})
return R(((-1)**(n*(n-1)/2)))*result
xs = ['x'+str(x) for x in range(n)]
R = PolynomialRing(QQ,names=xs)
r_gens = R.gens()
pol = R(self._convert_integrable_polynomial(p,n))
result = pol
for i in range(n):
result = result.integral(r_gens[n-1-i])
result = result.subs({r_gens[n-1-i]:1-sum(r_gens[:n-1-i])})
return R(((-1)**(n*(n-1)/2)))*result
def integral2(self,p,n):
## Input:
# - p: polynomial in variables t_0,..t_n,y_0,...,y_n
# - n: dimension of the algebra (apl)_n where p is defined
## Output:
# - the integral of the polynomial after reducing the variables t_n, y_n by the corresponding relations
#
xs = ['x'+str(x) for x in range(n)]
R = PolynomialRing(QQ,names=xs)
r_gens = R.gens()
pol = R(self._convert_integrable_polynomial(p,n))
result = 0
for k in pol.dict().keys(): # k is a monomial
numerator = 1
denominator = n
for e in k: # e is an exponents
numerator *= factorial(e)
denominator += e
coef = pol.dict()[k]
result += coef*numerator/factorial(denominator)
return R(((-1)**(n*(n-1)/2)))*result
#####################################################################
class MorphismFromQQToAplK(Morphism):
def __init__(self, codomain):
self.simplicialComplex= codomain.K
Morphism.__init__(self, QQ, codomain)
def _call_(self, x):
S = self.simplicialComplex
params1 = {}
for d in range(S.dimension()+1):
for c in S._n_cells_sorted(d):
params1[tuple(c)]=x
A = AplK(S)
return AplKElement(A,params1)
#####################################################################
class AplK(UniqueRepresentation, Parent):
def __init__(self,K, category=None):
self.K = K
self.n = K.dimension()
Parent.__init__(self,category=category )
def __init__(self,K, category=None):
self.K = K
self.n = K.dimension()
Parent.__init__(self,category=category )
f = MorphismFromQQToAplK(self)
self.register_coercion(f)
def _repr_(self):
return "Algebra AplK over %s" % self.K
def _repr_(self):
return "Algebra AplK over %s" % self.K
def _element_constructor_(self,values):
return self.element_class(self,values)
def _element_constructor_(self,values):
return self.element_class(self,values)
def dimension(self):
return self.K.dimension()
def dimension(self):
return self.K.dimension()
def simplicialComplex(self):
return self.K
def simplicialComplex(self):
return self.K
def _n_cells_sorted(self,n):
return self.K._n_cells_sorted(n)
def _n_cells_sorted(self,n):
return self.K._n_cells_sorted(n)
#####################################################################
class AplKElement(Element):
def __init__(self, parent,values):
self.elements = values
self.apl = Apln(parent.dimension())
Element.__init__(self,parent)
self.check()
def __init__(self, parent,values):
self.elements = values
self.apl = Apln(parent.dimension())
Element.__init__(self,parent)
self.check()
def _repr_(self):
return str(self.elements)
def _repr_(self):
return str(self.elements)
def dimension(self):
return self.parent().dimension()
def dimension(self):
return self.parent().dimension()
def nullElement(self):
C = self.__class__
values = {}
for k in self.elements.keys():
values[k]=0
return C(self.parent(),values)
def nullElement(self):
C = self.__class__
values = {}
for k in self.elements.keys():
values[k]=0
return C(self.parent(),values)
def _add_(self,other):
C = self.__class__
if(self.elements.keys()!=other.elements.keys()):
raise ValueError("error: keys do not match")
newValues = {}
for k in self.elements.keys():
newValues[k] = self.elements[k]+other.elements[k]
return C(self.parent(),newValues)
def _add_(self,other):
r"""
Return sum of the elements
INPUT:
- ``other`` -- An AplKEelement
OUTPUT: An AplKElement equal to the sum of self and the element passed as argument
EXAMPLES:
This example illustrates ... ::
sage: S = SimplicialComplex([[0,1],[0,2],[1,2]],is_mutable=False)
sage: apl = Apln(S.dimension())
sage: (t0, t1, y0, y1) = apl.generators()
sage: params1 = {(0,): t0, (1,): t0, (2,):3*t0, (0,1): t0+t1,(0,2): t0+3*t1, (1,2): t0+3*t1}
sage: params2 = {(0,): 3*t0, (1,): 4*t0, (2,):5*t0, (0,1): 3*t0+4*t1,(0,2): 3*t0+5*t1, (1,2): 4*t0+5*t1}
sage: A = AplK(S)
sage: e1 = AplKElement(A,params1)
sage: e2 = AplKElement(A,params2)
sage: e1+e2
{(0,): 4*t0, (1,): 5*t0, (2,): 8*t0, (0, 1): 4*t0 + 5*t1, (0, 2): 4*t0 + 8*t1, (1, 2): 5*t0 + 8*t1}
"""
C = self.__class__
if(self.elements.keys()!=other.elements.keys()):
raise ValueError("error: keys do not match")
newValues = {}
for k in self.elements.keys():
newValues[k] = self.elements[k]+other.elements[k]
return C(self.parent(),newValues)
def _sub_(self,other):
C = self.__class__
if(self.elements.keys()!=other.elements.keys()):
raise ValueError("error: keys do not match")
newValues = {}
for k in self.elements.keys():
newValues[k] = self.elements[k]-other.elements[k]
return C(self.parent(),newValues)
def _sub_(self,other):
C = self.__class__
if(self.elements.keys()!=other.elements.keys()):
raise ValueError("error: keys do not match")
newValues = {}
for k in self.elements.keys():
newValues[k] = self.elements[k]-other.elements[k]
return C(self.parent(),newValues)
def _mul_(self,other):
C = self.__class__
if(self.elements.keys()!=other.elements.keys()):
raise ValueError("error: keys do not match")
newValues = {}
for k in self.elements.keys():
newValues[k] = self.elements[k]*other.elements[k]
return C(self.parent(),newValues)
def _mul_(self,other):
C = self.__class__
if(self.elements.keys()!=other.elements.keys()):
raise ValueError("error: keys do not match")
newValues = {}
for k in self.elements.keys():
newValues[k] = self.elements[k]*other.elements[k]
return C(self.parent(),newValues)
def _richcmp_(self, other, op):
if(op==op_EQ):
if(self.elements.keys()!=other.elements.keys()):
return False
for k in self.elements.keys():
if(self.elements[k]!=other.elements[k]):
return False
return True
elif(op==op_NE):
return not self._richcmp_(other,op_EQ)
else:
raise NotImplementedError("operation not implemented for this element class")
def _richcmp_(self, other, op):
if(op==op_EQ):
if(self.elements.keys()!=other.elements.keys()):
return False
for k in self.elements.keys():
if(self.elements[k]!=other.elements[k]):
return False
return True
elif(op==op_NE):
return not self._richcmp_(other,op_EQ)
else:
raise NotImplementedError("operation not implemented for this element class")
def differential(self):
C = self.__class__
newValues = {}
for k in self.elements.keys():
newValues[k] = self.apl.cdga(self.elements[k]).differential()
return C(self.parent(),newValues)
def differential(self):
C = self.__class__
newValues = {}
for k in self.elements.keys():
newValues[k] = self.apl.cdga(self.elements[k]).differential()
return C(self.parent(),newValues)
def check(self):
S = self.parent().simplicialComplex()
dim = S.dimension()
for n in range(1,dim+1):
print("Checking dimension " + str(n))
for simplex in S._n_cells_sorted(n):
for f in range(n+1):
sFace = tuple(S.face(simplex,f))
aFace = self.apl.face(f)(self.elements[tuple(simplex)])
if(self.elements[sFace] != aFace):
raise ValueError("the assignation is not valid")
def check(self):
S = self.parent().simplicialComplex()
dim = S.dimension()
for n in range(1,dim+1):
#print("Checking dimension " + str(n))
for simplex in S._n_cells_sorted(n):
for f in range(n+1):
sFace = tuple(S.face(simplex,f))
aFace = self.apl.face(f)(self.elements[tuple(simplex)])
if(self.elements[sFace] != aFace):
raise ValueError("the assignation is not valid")
def homogeneous_component(self,degree):
C = self.__class__
component = {}
for k in self.elements.keys():
component[k] = self.apl.homogeneousPart(self.elements[k],degree)
return C(self.parent(),component)
def homogeneous_component(self,degree):
C = self.__class__
component = {}
for k in self.elements.keys():
component[k] = self.apl.homogeneousPart(self.elements[k],degree)
return C(self.parent(),component)
def homogeneous_components(self):
components = {}
for n in range(self.dimension()+1):
components[n] = self.homogeneous_component(n)
return components
def homogeneous_components(self):
components = {}
for n in range(self.dimension()+1):
components[n] = self.homogeneous_component(n)
return components
def is_homogeneous(self):
isHomogeneous = True
for k in self.elements.keys():
isHomogeneous = isHomogeneous and self.apl.cdga(self.elements[k]).is_homogeneous()
return isHomogeneous
def is_homogeneous(self):
isHomogeneous = True
for k in self.elements.keys():
isHomogeneous = isHomogeneous and self.apl.cdga(self.elements[k]).is_homogeneous()
return isHomogeneous
def is_in_kernel(self):
diff = self.differential()
if(diff == self.nullElement()):
return True
else:
return False
def is_in_kernel(self):
diff = self.differential()
if(diff == self.nullElement()):
return True
else:
return False
def _integrable_component(self):
C = self.__class__
component = {}
for k in self.elements.keys():
degree = len(k)-1
component[k] = self.apl.homogeneousPart(self.elements[k],degree)
#component[k] = self.homogeneous_component(dim).elements[k]
return component
def _integrable_component(self):
C = self.__class__
component = {}
for k in self.elements.keys():
degree = len(k)-1
component[k] = self.apl.homogeneousPart(self.elements[k],degree)
#component[k] = self.homogeneous_component(dim).elements[k]
return component
def integral(self):
component = self._integrable_component()
cochain = {}
#integration of each polynomial
for k in component.keys():
print("Integrating simplex " + str(k))
n = len(k)
pol = self.apl.cdga(component[k])
cochain[k] = self.apl.integral(pol,n)
return cochain
def integral(self):
component = self._integrable_component()
cochain = {}
#integration of each polynomial
for k in component.keys():
print("Integrating simplex " + str(k))
n = len(k)-1
pol = self.apl.cdga(component[k])
cochain[k] = self.apl.integral(pol,n)
return cochain
#####################################################################
#S = SimplicialComplex([[0,1,2]],is_mutable=False)
S = SimplicialComplex([[0,1],[0,2],[1,2]],is_mutable=False)
apl = Apln(S.dimension())
apl.inject_variables()
params1 = {}
params1[(0,)] = t0
params1[(1,)] = t0
params1[(2,)] = 3*t0
params1[(0,1)] = t0+t1
params1[(0,2)] = t0+3*t1
params1[(1,2)] = t0+3*t1
params2 = {}
params2[(0,)] = 3*t0
params2[(1,)] = 4*t0
params2[(2,)] = 5*t0
params2[(0,1)] = 3*t0+4*t1
params2[(0,2)] = 3*t0+5*t1
params2[(1,2)] = 4*t0+5*t1
params3 = {}
for d in range(S.dimension()+1):
for c in S._n_cells_sorted(d):
params3[tuple(c)]=randint(1,100)*t0
params4 = {}
params4[(0,)] = 3*t0+y0
params4[(1,)] = 4*t0
params4[(2,)] = 5*t0
params4[(0,1)] = 3*t0+4*t1+y0
params4[(0,2)] = 3*t0+5*t1+y0
params4[(1,2)] = 4*t0+5*t1
A = AplK(S)
e1 = AplKElement(A,params1)
e2 = AplKElement(A,params2)
#e3 = AplKElement(A,params3)
e4 = AplKElement(A,params4)
class morphism_from_QQ_to_AplK(Morphism):
def __init__(self, codomain):
self.simplicialComplex= codomain.K
Morphism.__init__(self, QQ, codomain)
def _call_(self, x):
S = self.simplicialComplex
params1 = {}
for d in range(S.dimension()+1):
for c in S._n_cells_sorted(d):
params1[tuple(c)]=x
A = AplK(S)
return AplKElement(A,params1)
f = morphism_from_QQ_to_AplK(A)
A.register_coercion(f)
#p = t0*t1**2*y0+t0*y1
S = SimplicialComplex([[0,1,2,3,4,5,6]],is_mutable=False)
apl = Apln(S.dimension())
apl.inject_variables()
p = y0*y1*y2*y3*y4*y5+y0*y1*y2*y3*y6+t4*t5
initial_assignation = {}
initial_assignation[(0,1,2,3,4,5,6)] = p
def generateElement(S,initial_assignation):
apl = Apln(S.dimension())
aplk_structure = initial_assignation.copy()
for n in reversed(range(1,S.dimension()+1)):
print("Generating dimension "+ str(n))
for simplex in S._n_cells_sorted(n):
p = aplk_structure[tuple(simplex)]
for f in range(n+1):
sFace = tuple(S.face(simplex,f))
aFace = apl.face(f)(p)
if(sFace not in aplk_structure.keys()):
aplk_structure[tuple(sFace)] = aFace
elif(aplk_structure[tuple(sFace)] != aFace):
raise ValueError("Error during the checking. The assignation is not valid")
A = AplK(S)
return AplKElement(A,aplk_structure)
g = generateElement(S,initial_assignation)
#g.check()
#cochain = g.integral()
#n = 1
#a = Apln(n)
#a.inject_variables()
#p = t0*y1
##print(p)
#print(a._reduce_polynomial(p,n))
#print(a._convert_integrable_polynomial(p,n))
#print(a.integral(p,n))