AplK/aplK.py
Miguel Marco 71cb7aef68
Make sure computing half is an integer, and caching method
In some cases taking half of an integer results in a fractional type
(even if it represents an integer number). This causes problems with indices,
so this change uses division with quotient and remainder instead.

Also, it caches the degeneration morphism, which is expensive to create.
2022-03-30 14:11:46 +02:00

1290 lines
47 KiB
Python

from sage.algebras.commutative_dga import *
from sage.structure.parent import Parent
from sage.structure.element import Element
from sage.structure.richcmp import rich_to_bool, op_EQ, op_NE, op_LT, op_LE, op_GT, op_GE
from sage.rings.polynomial.polydict import ETuple
import copy
class AplN(DifferentialGCAlgebra_multigraded, UniqueRepresentation):
"""
A differential multi-graded commutative algebra
This class extends the DifferentialGCAlgebra_multigraded class
to create the algebra of polynomial differential forms with coefficients in QQ
INPUT:
- ``A`` -- a multi-graded commutative algebra
- ``differential`` -- the differential rules
"""
def __init__(self, A, differential):
"""
Initialize ``self``.
INPUT:
- ``A`` -- a multi-graded commutative algebra
- ``differential`` -- a differential
"""
self.n = len(A._names) // 2 - 1
DifferentialGCAlgebra_multigraded.__init__(self, A, differential)
self.t_relation = sum(self.gens()[:self.n+1])-1
self.y_relation = sum(self.gens()[self.n+1:])
def _base_repr(self):
"""
Return the base string representation of ``self``.
"""
s = DifferentialGCAlgebra._base_repr(self)
old = '{}'.format(self._degrees)
new = '{}'.format([i[0] for i in self._degrees_multi])
return s.replace(old, new)
def integrand(self, value):
"""
Returns a polynomial whose definite integral is equal to 'value'
OUTPUT:
- A polynomial
"""
if (value in QQ) or ((value in A) and (value.degree().is_zero())):
return value*factorial(self.n)*prod(self.gens()[self.n+2:])
else:
raise ValueError("only field elements can be integrated.")
class Element(GCAlgebra_multigraded.Element, DifferentialGCAlgebra.Element):
"""
Element class of the commutative differential multi-graded algebra (Apl)n.
"""
def is_homogeneous(self):
"""
Return "True" if "self" is homogeneous and "False" otherwise.
"""
degrees = set([k[0] for k in self.homogeneous_parts()])
return len(degrees) == 1
def total_degree(self):
"""
Return the total degree of this element.
"""
degrees = set([k[0] for k in self.homogeneous_parts()])
return max(degrees)
def homogeneous_part(self,degree):
"""
Returns the homogeneous part of the element for a given degree.
"""
if (degree<0):
raise ValueError('`degree` must be a non-negative integer')
parts = [self.homogeneous_parts()[k] for k in self.homogeneous_parts() if list(k)[0]==degree]
return self.parent()(sum(parts))
def reduce(self, relations=[]):
""""
The method reduces the polynomial by the Ideal <sum(t_i)-1, sum(y_j)>
INPUT::
- p: polynomial in variables t_0,..t_n,y_0,...,y_n
- relations (optional): if relations are given, the polynomial is reduced following
this rules. If no relations are given, the polynomial is reduced taking into
account the relations of the parent algebra.
OUTPUT::
- The polynomial after applying the corresponding relations
"""
apl = self.parent()
if relations:
return DifferentialGCAlgebra.Element.reduce(self, relations)
else:
return DifferentialGCAlgebra.Element.reduce(self, [apl.t_relation, apl.y_relation])
def _integrable_part(self):
"""
The method returns the coefficient in Q[t_1,...,t_n] of the term
of degree `n`
OUTPUT:
- A polynomial in the ring Q<t_1,..t_n>.
(After reducing the variables t_{0},y_{0} by the corresponding relations,
the algorithm removes the variables y_{1},..,y_{n}).
"""
apl = self.parent()
pol = apl(self.reduce()) #needed in case reduce() returns an integer
relations = []
for i in range(apl.n+1):
relations.append(apl.gens()[apl.n+1+i]-1)
pol = apl(pol.homogeneous_part(degree=apl.n))
return apl(pol.reduce(relations=relations))
def integral(self):
"""
Returns the definite integral of the term of degree `n`
OUTPUT:
- The rational number obtained by the definite integral
of the term of degree `n`
"""
apl = self.parent()
pol = self._integrable_part()
result = 0
for k in pol.dict().keys(): # k is a monomial
numerator = 1
denominator = apl.n
for e in k: # e is an exponents
numerator *= factorial(e)
denominator += e
coef = pol.dict()[k]
result += coef*numerator/factorial(denominator)
#return ((-1) ** ((apl.n*(apl.n-1)) // 2)) * result
return result
def _get_term(self, coef, exp, gens):
""""
Returns the term specified by its coefficient, exponents and generators
INPUT::
- coef: The coefficient of the term
- exp: A tuple with the exponents of each variable
- gens: The generators
OUTPUT::
- A polynomial
"""
return coef*prod(list(map(lambda exp, gen: gen**exp, exp, gens)))
def _get_terms(self):
""""
Returns self seen as a collection of terms in the variables y_j with coefficients in QQ<t_1,..,t_n>
OUTPUT::
- A dictionary with all terms
"""
a = self.parent()
n = len(a.gens())//2
ts_only_gens = a.gens()[:n]+tuple(n*[1])
ys_only_gens = tuple(n*[1])+a.gens()[n:]
terms = {}
pol = self.dict()
for exp in pol.keys():
key = self._get_term(1, exp,ys_only_gens)
coef = self._get_term(pol[exp],exp,ts_only_gens)
if key not in terms.keys():
terms[key] = coef
else:
terms[key] += coef
return terms
def _get_leading_pol(self):
""""
Returns the leading polynomial
OUTPUT::
- A polynomial
"""
terms = self._get_terms()
keys = list(terms.keys())
keys.sort(reverse=True)
return terms[keys[0]]*keys[0]
def _get_idx(self, exp):
""""
Returns the idx 'i' of the first generator 'y_i' that appears in exp
INPUT::
- exp: A tuple with the exponents
OUTPUT::
- A non-negative number lower or equal to 'n'
"""
n = len(exp)//2
filter_ts = n*[0]+n*[1]
diff_forms = list(map(lambda v, f: v*f, exp, filter_ts))
return next((i for i, x in enumerate(list(diff_forms)) if x), None)-n
def _monomial_primitive(self, monomial, idx):
""""
Returns the primitive of the given monomial
INPUT::
- monomial: A monomial
- idx: the index of the variable y_j that will be used to calculate the pseudo-primitive
OUTPUT::
- A polynomial
"""
exp = list(monomial.dict().keys())[0]
coef = monomial.dict()[exp]
n = len(exp)//2
if (exp[n+idx]==0):
raise IndexError("It is not possible to get a primitive using the given idx.")
primitive = list(copy.deepcopy(exp))
primitive[n+idx] = 0
primitive[idx] += 1
new_coef = coef/primitive[idx]
a = monomial.parent()
return self._get_term(new_coef, primitive, a.gens())
def _leading_primitive(self):
""""
Returns the primitive of the leading term of self
OUTPUT::
- A polynomial
"""
polynomial = self._get_leading_pol()
# necessary condition: all monomials have apl.gens[n+idx] factor
exp = list(polynomial.dict().keys())[0]
idx = self._get_idx(exp)
aGens = polynomial.parent().gens()
primitive = 0
for exp in polynomial.dict().keys():
coef = polynomial.dict()[exp]
monomial = self._get_term(coef, exp, aGens)
primitive += self._monomial_primitive(monomial, idx)
terms = primitive._get_terms()
result = sum([terms[key]*key for key in terms.keys()])
return result
def primitive(self):
""""
Retuns a polynomial whose differential is self
OUTPUT::
- A polynomial
"""
if self.differential()!=0 :
raise ValueError("This polynomial is not closed.")
w = self._leading_primitive()
theta = self - w.differential()
if (theta == 0):
return w
else:
return w+theta.primitive()
def AlgebraPL(n):
"""
The graded commutative differential cochain algebra Apl of dimension ``n``
INPUT:
- ``n`` -- A non-negative integer corresponding to the number of generators of degree 0.
"""
if (n<0):
raise ValueError('`n` must be a non-negative integer')
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)])
gca = GradedCommutativeAlgebra(QQ, names=ts+ys, degrees=degrees)
differentials = {gca.gens_dict()[ts[i]]:gca.gens_dict()[ys[i]] for i in range(n+1)}
return AplN(gca,differentials)
#############################################
# TEST
# print("Test AlgebraPL")
#A = AlgebraPL(5)
#A.inject_variables()
#p = t0**2 + t1 + y3 + y1*y2
#q = p.differential()
#############################################
class SimplicialAplN(AplN, UniqueRepresentation):
"""
The simplicial commutative cochain algebra Apl
"""
def _base_repr(self):
"""
Return the base string representation of ``self``.
"""
s = DifferentialGCAlgebra._base_repr(self)
old = '{}'.format(self._degrees)
new = '{}'.format([i[0] for i in self._degrees_multi])
return "Simplicial " + s.replace(old, new)
class Element(AplN.Element):
"""
Element class of the simplicial commutative differential algebra (Apl)n.
"""
def face_morphism(self,i):
""""
Returns the morphism face 'i' of self
INPUT::
- i: A non-negative integer lower or equal to 'n'. Represents the face index
OUTPUT::
- A morphism of cochain algebras
"""
n = self.parent().n
if(n==0):
return 0
domain = self.parent()
codomain = SimplicialAlgebraplN(n-1)
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-1+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-2])
images = tImages+yImages
return domain.hom(codomain=codomain,im_gens=images)
def face(self,i):
""""
Returns the face 'i' of self
INPUT::
- i: A non-negative integer lower or equal to 'n'. Represents the face index
OUTPUT::
- A polynomial that is the face 'i' of self
"""
return self.face_morphism(i)(self)
@cached_method
def degeneration_morphism(self,j):
""""
Returns the morphism degeneration 'j' of self
INPUT::
- j: A non-negative integer lower or equal to 'n'. Represents the degeneration index
OUTPUT::
- A morphism of cochain algebras
"""
n = self.parent().n
domain = self.parent()
codomain = SimplicialAlgebraplN(n+1)
if(j>n):
print("the degeneration index must be a non-negative number lower or equal to "+str(n))
return False
elif(j<0):
print("the degeneration index must be non-negative")
return False
else:
tImages = []
yImages = []
n = len(domain.gens())//2
for k in range(n):
if(k<j):
tImages.append(codomain.gens()[k])
yImages.append(codomain.gens()[n+1+k])
elif (k==j):
tImages.append(codomain.gens()[k]+codomain.gens()[k+1])
yImages.append(codomain.gens()[n+1+k]+codomain.gens()[n+1+k+1])
elif (k>j):
tImages.append(codomain.gens()[k+1])
yImages.append(codomain.gens()[n+1+k+1])
images = tImages+yImages
return domain.hom(codomain=codomain,im_gens=images)
def degeneration(self,j):
""""
Returns the degeneration 'j' of self
INPUT::
- j: A non-negative integer lower or equal to 'n'. Represents the degeneration index
OUTPUT::
- A polynomial that is the degeneration 'j' of self
"""
return self.degeneration_morphism(j)(self)
def propagate(self, face_list, dest_dim):
"""
Returns a polynomial in dimension 'dest_dim' that has 'self' as a face
INPUT:
- face_list: indexes of the elements that form the face of the simplex we want to propagate
- dest_dim: Dimension of the codomain cochain algebra
- EXAMPLE: s = SimplicialAlgebraplN(2)
s.inject_variables()
p = (t1^2+t2^2)*y1*y2
p.propagate2([0,1,2],3)
t1^2*t2*y0*y1 + t2^3*y0*y1 - t1^3*y0*y2 - t1*t2^2*y0*y2 + t0*t1^2*y1*y2 + t0*t2^2*y1*y2
"""
if not (self.is_homogeneous() and self.total_degree()==self.parent().n):
raise ValueError("The polynomial is not homogeneous of maximum degree.")
self_dim = self.parent().n
indexes = list(range(self_dim+1))
apl = SimplicialAlgebraplN(dest_dim)
y_gens = [apl.gens()[dest_dim+1+i] for i in face_list]
t_gens = [apl.gens()[i] for i in face_list]
graded_part = apl(0)
for i in reversed(indexes):
t = t_gens[i]
y_prod = copy.deepcopy(y_gens)
del y_prod[i]
y_prod = prod(y_prod)
graded_part += ((-1)**i)*t*y_prod#((-1)**(simplex_dim+i))*t*y_prod
terms = self.dict()
coeff_in_t = apl(0)
for exp in terms.keys():
coef = terms[exp]
ts_exp = exp[:self_dim+1]
coeff_in_t += coef*prod(list(map(lambda t, e: t**e, t_gens, ts_exp)))
return coeff_in_t*graded_part
def _get_primitive_basis_of_a_monomial(self, monomial):
"""
Returns a basis of the primitive of a given monomial
Example:
"""
n_vars = len(monomial) // 2
combinations = []
for i in range(n_vars):
mono = list(copy.deepcopy(monomial))
if mono[n_vars+i]:
mono[i] += 1
mono[n_vars+i] = 0
combinations.append(tuple(mono))
return list(set(combinations))
def _add_kernel_terms_to_primitive_basis(self, basis):
new_elements = []
for e in basis:
monomial = list(copy.deepcopy(e))
n = len(e) // 2
monomial[:n] = n*[0]
new_elements.append(tuple(monomial))
return new_elements
def _get_primitive_basis(self):
"""
Returns a basis of the primitive of 'self'
"""
basis = []
p = self.reduce()
for monomial in p.dict().keys():
basis += p._get_primitive_basis_of_a_monomial(monomial)
basis = list(set(basis))
basis += p._add_kernel_terms_to_primitive_basis(basis)
basis = list(set(basis))
terms = []
apl = self.parent()
for b in basis:
dim = len(monomial)//2 - 1
terms.append(self._get_term(1, b, apl.gens()))
terms.sort(reverse=False)
return terms
def _get_differential_basis(self):
"""
Returns a basis of the differential of 'self'
"""
primitive_basis = self._get_primitive_basis()
differential_basis = []
for b in primitive_basis:
diff = b.differential().reduce()
if diff: # as we want a basis, we don't consider the coeffient of the differential
for monomial in diff.dict().keys():
term = self._get_term(1, monomial, self.parent().gens())
differential_basis.append(term)
differential_basis = list(set(differential_basis))
differential_basis.sort(reverse=False)
return differential_basis
def _get_coordinate_vector_differential(self, diff_basis):
"""
Returns the coordinates of 'self'
in terms of the differential basis
"""
p_reduced = self.reduce()
#dic = p.dict()
#terms = [p_reduced.parent()({t: dic[t]}) for t in dic.keys()]
monomials = list(p_reduced.dict().keys())
v = vector(QQ,len(diff_basis)*[0])
for i in range(len(diff_basis)):
t = diff_basis[i]
term = list(t.dict().keys())[0] # there is only one monomial
if term in monomials:
v[i] = p_reduced.dict()[term]
return v
def _differential_matrix(self, primitive_basis, differential_basis):
"""
Returns the matrix of the differential map
"""
#primitive_basis = get_primitive_basis_of_a_polynomial(p.reduce())
#differential_basis = get_differential_basis_of_a_polynomial(p.reduce())
nrows = len(differential_basis)
ncols = len(primitive_basis)
M = matrix(QQ, nrows, ncols, nrows*ncols*[0])
for j in range(ncols):
primitive = primitive_basis[j]
#diff_dict = self._tuple_differential(primitive)
diff_primitive = primitive.differential().reduce()
monomials = list(diff_primitive.dict().keys())
for i in range(nrows):
diff_basis_element = differential_basis[i]
term = list(diff_basis_element.dict().keys())[0] # there is only one monomial
if term in monomials:
M[i,j] = diff_primitive.dict()[term]
return M
def solve_primitive(self):
p = self.reduce()
primitive_basis = p._get_primitive_basis()
diff_basis = p._get_differential_basis()
v = p._get_coordinate_vector_differential(diff_basis)
Mdif = p._differential_matrix(primitive_basis, diff_basis)
sol = Mdif.solve_right(v)
return sum(map(lambda x,y: x*y, sol, primitive_basis))
def SimplicialAlgebraplN(n):
"""
Returns the graded commutative differential cochain algebra (Apl)_n as part of
the simplicial cochain algebra A_{PL}
INPUT:
- ``n`` -- A non-negative integer corresponding to the number of generators of degree 0.
"""
if (n<0):
raise ValueError('`n` must be a non-negative integer')
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)])
gca = GradedCommutativeAlgebra(QQ, names=ts+ys, degrees=degrees)
differentials = {gca.gens_dict()[ts[i]]:gca.gens_dict()[ys[i]] for i in range(n+1)}
return SimplicialAplN(gca,differentials)
#############################################
# TEST
#print("TEST SimplicialAlgebraN")
# s = SimplicialAlgebraplN(5)
# s.inject_variables()
# ## face and degeneration maps
# print("Test face and degeneration maps")
# p = t1+y2
# p.degeneration(0).face(1)==p
# m1 = p.degeneration(1).face_morphism(1)
# m2 = p.degeneration_morphism(1)
# print((m1*m2).is_identity()) # true
# ## Propagation
# print("Test progatation")
# s = SimplicialAlgebraplN(2)
# s.inject_variables()
# p = (t1**2+t2**2)*y1*y2
# q = p.propagate([1,5,7],10)
# print(q)
# print(q.face(10).face(9).face(8).face(6).face(4).face(3).face(2).face(0).reduce() == p)
#############################################
class AlgebraAplK(UniqueRepresentation, Parent):
"""
The graded commutative differential cochain algebra Apl(K)
"""
def __init__(self,simplicial_complex, category=None):
"""
The constructor of Apl(K)
INPUT:
- simplicial_complex: The simplicial complex used to build the algebra
"""
self.K = simplicial_complex
self.n = simplicial_complex.dimension()
Parent.__init__(self,category=category)
f = MorphismFromQQToAplK(self)
self.register_coercion(f)
def _repr_(self):
"""
Return the base string representation of ``self``.
"""
return "The graded commutative differential cochain algebra Apl(K) of %s" % self.K
def _element_constructor_(self,values):
"""
The constructor of an element of this class
"""
return self.element_class(self,values)
def dimension(self):
"""
The dimension of the algebra. It is equivalent to the maximum degree possible that a non-null element of this algebra can take.
"""
return self.n
def simplicialComplex(self):
"""
The simplicial complex of this algebra.
"""
return self.K
def _n_cells_sorted(self,n):
"""
Returns all simplices of dimension 'n', as a sorted list.
"""
return self.K._n_cells_sorted(n)
def _faces_of_length(self, simplex, length):
"""
Returns all faces of dimension 'length' for a given simplex
"""
return [tuple(e) for e in Set(simplex).subsets(length)]
def null_element(self):
"""
Returns the null element of the algebra respect to the addition operation
"""
values = {}
for i in range(self.K.dimension()+1):
s = SimplicialAlgebraplN(i)
for f in self.K._n_cells_sorted(i):
values[tuple(f)]=s(0)
return self._element_constructor_(values)
def lift_cochain(self, cochain):
"""
Returns an element of the algebra Apl(K) whose integral is the cochain given as input.
"""
# first we find polynomials compatible with the integration value
values = {}
max_dim = 0
for k in cochain.keys():
dim = len(k)-1
if dim > max_dim:
max_dim = dim
apl = SimplicialAlgebraplN(dim)
if dim>0:
values[k] = ((-1) ** ((dim * (dim-1)) // 2)) * apl.integrand(cochain[k])
else:
values[k]= apl(cochain[k])
# second, we propagate the simplicial information
simplicial_structure = {}
for n in reversed(range(max_dim+1)):
keys = [k for k in values.keys() if len(k) == n+1] # keys of the n-simplices
for k in keys:
simplicial_structure[k] = values[k] + self._propagate(cochain, k)
#print(simplicial_structure)
return AlgebraAplK.Element(self,simplicial_structure)
def _propagate(self, cochain, simplex):
"""
Returns the part of the polynomial assigned to 'simplex' related to the
propagation of the simplicial information of its faces
"""
dim = len(simplex)-1
apl = SimplicialAlgebraplN(dim)
propagation = apl(0)
for n in reversed(range(dim)):
faces = self._faces_of_length(simplex, n+1)
apl = SimplicialAlgebraplN(n)
pol = apl(prod(apl.gens()[n+2:])) # we propagate the polynomial (y_1...y_n)
for f in faces:
if f in cochain.keys():
value = cochain[f]
face_indexes = [simplex.index(i) for i in f]
#print("dim: ",n," pol: ",pol," face_idx: ",face_indexes," dim: ", n)
propagation += value*((-1) ** ((n * (n-1)) // 2)) * factorial(n) * pol.propagate(face_indexes, dim)
return propagation
def cohomology_gens(self,degree):
"""
Returns the cohomolosy generators of a given degree
"""
gens = {}
basis = self.K.cohomology_ring(QQ).basis()
for k in basis.keys():
if k[0] == degree:
gens[k] = basis[k]
return gens
def cohomology_group(self, n):
chainComplex = self.K.chain_complex(base_ring=QQ)
cochainComplex = chainComplex.dual()
Hn = cochainComplex.differential(n).right_kernel()/cochainComplex.differential(n-1).column_space()
return Hn
def cochain_to_vector(self, cochain, dim):
"""
Converts a cochain, given as a dictionary, to a vector
"""
basis = self.K.n_chains(dim, base_ring=QQ, cochains=True).basis()
v = vector(QQ,len(basis)*[0])
i=0
for simplex in basis.keys():
v[i] = cochain[tuple(simplex)]
i = i+1
return v
def vector_to_cochain(self, v, dim):
"""
Converts a cochain, given as a vector, to a dictionary
"""
cochain = {}
for d in range(self.K.dimension()+1):
basis = self.K.n_chains(d, base_ring=QQ, cochains=True).basis()
if d!=dim:
for simplex in basis.keys():
cochain[tuple(simplex)] = 0
else:
i=0
for simplex in basis.keys():
cochain[tuple(simplex)] = v[i]
i = i+1
return cochain
def lift_cochain_complex(self, cochain, dim):
"""
(REVISAR)
Returns an element of the preimage by the differential of a given cochain
INPUT:
- cochain: a cochain as a vector
- dim: dimension of the cochain
"""
# S = the Simplicial Complex
# cochain = the cochain you want to find a preimage by the differential
# dim = the dimension of the cochain
cochain_vector = self.cochain_to_vector(cochain, dim)
cochainComplex = self.K.chain_complex(base_ring=QQ).dual()
diff = cochainComplex.differential(dim-1)
show(diff.echelon_form())
result = diff.solve_right(cochain_vector)
return self.vector_to_cochain(result)
def lift_cohomology_generator(self, g):
"""
Returns a representative cochain associated to a given cohomology generator.
"""
cycle = g.to_cycle()
cochain = {}
for n in range(self.K.dimension()+1):
for simplex in self.K._n_cells_sorted(n):
cochain[simplex.tuple()] = 0
for i in cycle:
#print(type(i[0].tuple()))
cochain[i[0].tuple()] = i[1]
return cochain
def get_cohomology_aplk(self,degree):
"""
Returns a set of representative elements of Apl(K), one for each cohomology generator of a given degree
"""
representative_aplk_elements = {}
dic_gens = self.cohomology_gens(degree)
for key in dic_gens.keys():
gen = dic_gens[key]
cochain = self.lift_cohomology_generator(gen)
representative_aplk_elements[key] = self.lift_cochain(cochain)
return representative_aplk_elements
def first_cohomology_group(self):
"""
Returns the degree >= 1 of the first non-null cohomology generator.
"""
basis = self.K.cohomology_ring(base_ring=QQ).basis()
indexes = list(basis.keys())
indexes.sort(reverse=False)
# print(indexes)
for n in range(1,len(indexes)):
if indexes[n][0] > 0:
return indexes[n][0]
return -1
class Element(Element):
"""
Element class of the simplicial commutative differential algebra (Apl)n.
"""
def __init__(self, parent, values):
Element.__init__(self,parent)
self.elements = self._reduce(values)
self.check()
def _repr_(self):
"""
Return the base string representation of ``self``.
"""
return str(self.elements)
def dimension(self):
"""
Returns the dimension where 'self' lives
"""
return self.parent().dimension()
def _add_(self,other):
"""
Addition operation
"""
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):
"""
Substration operation
"""
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):
"""
Multiplication operation
"""
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 __pow__(self, n):
"""
Return ``self`` raised to the power ``n``.
"""
C = self.__class__
newValues = {}
for k in self.elements.keys():
newValues[k] = self.elements[k]**n
return C(self.parent(),newValues)
def _richcmp_(self, other, op):
"""
Comparison operations. Only implemented equality operations (EQ and NE)
"""
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 _reduce(self, values):
"""
Returns the polynomials given as input, after applying the corresponding relations
of the algebra where each polynomial belongs.
INPUT:
- values: a dictionary where the values are polynomials in (Apl)n
"""
reduced_values = {}
for k in values.keys():
n = len(k)-1
reduced_values[k] = values[k].reduce()
return reduced_values
def check(self):
"""
Checks if the initial assignation respects the simplicial structure. For every simplex S, the algorithm checks
if the polynomial assigned to each face of S is equal to the same face of the polynomial assigned to S.
"""
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.elements[simplex.tuple()].face(f).reduce()
#print("face: " + str(f),"simplex: " + str(simplex),"sFace: " + str(self.elements[sFace]), "aFace: " + str(aFace))
if(self.elements[sFace] != aFace):
raise ValueError("the assignation is not valid")
def differential(self):
"""
Returns the differential of 'self'
"""
C = self.__class__
newValues = {}
for k in self.elements.keys():
newValues[k] = self.elements[k].differential()
return C(self.parent(),newValues)
def is_closed(self):
"""
Returns 'True' if 'self' is closed, i.e., its differential is zero.
"""
diff = self.differential()
if(diff == self.parent().null_element()):
return True
else:
return False
def homogeneous_part(self,degree):
"""
Returns the homogeneous component of grade 'degree' of 'self'
"""
C = self.__class__
part = {}
for k in self.elements.keys():
s = SimplicialAlgebraplN(len(k)-1)
part[k] = s(self.elements[k].homogeneous_part(degree=degree))
return C(self.parent(),part)
def homogeneous_parts(self):
"""
Returns all graded components of 'self'
"""
parts = {}
for n in range(self.dimension()+1):
hom_part = self.homogeneous_part(degree=n)
if hom_part != self.parent().null_element():
parts[n] = hom_part
return parts
def is_homogeneous(self):
"""
Returns 'True' if self is composed by homogeneous polynomials of the same degree
"""
isHomogeneous = True
S = self.parent().simplicialComplex()
dim = S.dimension()
for n in range(dim+1):
for simplex in S._n_cells_sorted(n):
pol = self.elements[simplex.tuple()]
if pol:
isHomogeneous = isHomogeneous and pol.is_homogeneous()
if not isHomogeneous:
break
if not isHomogeneous:
break
return isHomogeneous
def degree(self):
if not self.is_homogeneous():
return -1
else:
degree = -1
for e in self.elements.values():
if e:
degree = e.degree()[0]
break
return degree
def integral(self):
"""
Returns the cochain obtained as integration of 'self'
"""
cochain = {}
#integration of each polynomial
for k in self.elements.keys():
dim = len(k)-1
# the sign is due to stoke's theorem
cochain[k] = ((-1) ** ((dim * (dim-1)) // 2)) * self.elements[k].integral()
return cochain
def _face_basis(self, basis, face_idx):
face_basis = []
for b in basis:
face_pol = b.face(face_idx).reduce()
if face_pol:
for monomial in face_pol.dict().keys():
#remove coefficients
term = face_pol._get_term(1, monomial, face_pol.parent().gens())
face_basis.append(term)
face_basis = list(set(face_basis))
face_basis.sort(reverse=False)
return face_basis
def _get_primitive_basis(self):
primitive_basis = {}
S = self.parent().simplicialComplex()
dim = S.dimension()
min_dim = self.degree()
for n in reversed(range(min_dim,dim+1)):
for simplex in S._n_cells_sorted(n):
if simplex.tuple() not in primitive_basis.keys():
primitive_basis[simplex.tuple()] = []
# add primitive elements related to the assigned polynomial
pol = self.elements[simplex.tuple()]
primitive_basis[simplex.tuple()] += pol._get_primitive_basis()
primitive_basis[simplex.tuple()] = list(set(primitive_basis[simplex.tuple()]))
primitive_basis[simplex.tuple()].sort(reverse=False)
# propagate the primitive basis through its faces
for f in range(n+1):
simplex_face = S.face(simplex,f).tuple()
if simplex_face not in primitive_basis.keys():
primitive_basis[simplex_face] = []
primitive_basis[simplex_face] += self._face_basis(primitive_basis[simplex.tuple()], f)
primitive_basis[simplex_face] = list(set(primitive_basis[simplex_face]))
primitive_basis[simplex_face].sort(reverse=False)
return primitive_basis
def _diff_basis(self, basis):
diff_basis = []
for b in basis:
b_diff = b.differential().reduce()
if b_diff:
for monomial in b_diff.dict().keys():
#remove coefficients
term = b_diff._get_term(1, monomial, b_diff.parent().gens())
diff_basis.append(term)
diff_basis = list(set(diff_basis))
diff_basis.sort(reverse=False)
return diff_basis
def _get_differential_basis(self, primitive_basis):
differential_basis = {}
S = self.parent().simplicialComplex()
dim = S.dimension()
min_dim = self.degree()
for n in reversed(range(min_dim,dim+1)):
for simplex in S._n_cells_sorted(n):
prim_diff_basis = primitive_basis[simplex.tuple()]
differential_basis[simplex.tuple()] = self._diff_basis(prim_diff_basis)
return differential_basis
def _get_differential_matrix(self, p_basis, d_basis):
S = self.parent().simplicialComplex()
dim = S.dimension()
min_dim = self.degree()
ncols = sum([len(l) for l in p_basis.values()])
nrows = sum([len(l) for l in d_basis.values()])
M = zero_matrix(QQ, nrows, ncols)
row_offset = 0
col_index = 0
for n in reversed(range(min_dim-1,dim+1)):
for simplex in S._n_cells_sorted(n):
# print(simplex.tuple(), row_offset, col_index)
if simplex.tuple() in p_basis.keys():
primitives = p_basis[simplex.tuple()]
if simplex.tuple() in d_basis.keys():
differentials = d_basis[simplex.tuple()]
if primitives and differentials:
for p in primitives:
dif_p = p.differential().reduce()
dic = dif_p.dict()
for monomial in dic.keys():
term = dif_p._get_term(1, monomial, dif_p.parent().gens())
if term in differentials:
row_index = differentials.index(term)
M[row_offset + row_index, col_index] = dic[monomial]
col_index += 1
row_offset += len(differentials)
elif primitives:
col_index += len(primitives)
elif primitives:
col_index += len(primitives)
return M
def _get_vector_diff_coordinates(self, d_basis):
len_v = sum([len(l) for l in d_basis.values()])
v = zero_vector(QQ, len_v)
offset = 0
S = self.parent().simplicialComplex()
dim = S.dimension()
min_dim = self.degree()
for n in reversed(range(min_dim,dim+1)):
for s in S._n_cells_sorted(n):
simplex = s.tuple()
if simplex in d_basis.keys():
differentials = d_basis[simplex]
pol = self.elements[simplex].reduce()
if pol and differentials:
dic = pol.dict()
for monomial in dic.keys():
term = pol._get_term(1, monomial, pol.parent().gens())
if term in differentials:
v_index = differentials.index(term)
v[offset + v_index] = dic[monomial]
offset += len(differentials)
return v
def _get_offsets_matrix(self, p_basis):
offsets = {}
S = self.parent().simplicialComplex()
dim = S.dimension()
min_dim = self.degree()
offset = 0
first_simplex = True
for n in reversed(range(min_dim-1,dim+1)):
for simplex in S._n_cells_sorted(n):
if simplex.tuple() in p_basis.keys():
primitives = p_basis[simplex.tuple()]
if first_simplex:
offsets[simplex.tuple()] = 0
first_simplex = False
else:
offsets[simplex.tuple()] = offset
offset += len(primitives)
return offsets
def _get_row_face_matrix(self, primitives, face_idx, face_basis,
offset_simplex, offset_face, len_row):
M = zero_matrix(QQ,1, len_row)
for fb in face_basis:
fb_idx = face_basis.index(fb)
row = zero_matrix(QQ,1, len_row)
row[0,offset_face+fb_idx] = -1
for term in primitives:
face_term = term.face(face_idx).reduce()
dic = face_term.dict()
i = primitives.index(term)
for monomial in dic:
t = face_term._get_term(1, monomial, face_term.parent().gens())
if t == fb:
coef = dic[monomial]
row[0,offset_simplex+i] = coef
M = M.stack(row)
return M[1:]
def _get_face_matrix(self, p_basis):
offsets = self._get_offsets_matrix(p_basis)
len_row = sum([len(l) for l in p_basis.values()])
M = zero_matrix(QQ,1, len_row)
min_dim = self.degree()
for n in reversed(range(min_dim,dim+1)):
for simplex in S._n_cells_sorted(n):
offset_simplex = offsets[simplex.tuple()]
primitives = p_basis[simplex.tuple()]
for f in range(n+1):
simplex_face = S.face(simplex,f).tuple()
face_primitives = p_basis[simplex_face]
offset_face = offsets[simplex_face]
rows = self._get_row_face_matrix(primitives, f, face_primitives,
offset_simplex, offset_face,
len_row)
M = M.stack(rows)
return M[1:]
def _get_solution(self, sol, p_basis):
S = self.parent().simplicialComplex()
dim = S.dimension()
min_dim = self.degree()-1
offset = 0
values = {}
for n in reversed(range(dim+1)):
apl = SimplicialAlgebraplN(n)
for s in S._n_cells_sorted(n):
simplex = s.tuple()
pol = apl(0)
if simplex in p_basis.keys():
primitives = p_basis[simplex]
if primitives:
for p in primitives:
i = primitives.index(p)
coef = sol[offset+i]
pol += coef*p
offset += len(primitives)
values[simplex] = pol
A = self.parent()
# return values
return A._element_constructor_(values)
def solve_primitive(self):
p_basis = self._get_primitive_basis()
d_basis = self._get_differential_basis(p_basis)
diff_matrix = self._get_differential_matrix(p_basis, d_basis)
v_dif = self._get_vector_diff_coordinates(d_basis)
face_matrix = self._get_face_matrix(p_basis)
nrows = face_matrix.nrows()
v = vector(QQ, list(v_dif)+nrows*[0])
M = diff_matrix.stack(face_matrix)
sol = M.solve_right(v)
return self._get_solution(sol, p_basis)
class MorphismFromQQToAplK(Morphism):
def __init__(self, codomain):
self.simplicialComplex = codomain.K
Morphism.__init__(self, QQ, codomain)
def _call_(self, x):
S = self.simplicialComplex
params = {}
for d in range(S.dimension()+1):
s = SimplicialAlgebraplN(d)
for c in S._n_cells_sorted(d):
params[tuple(c)]=s(x)
A = AlgebraAplK(S)
return AlgebraAplK.Element(A, params)
#############################################
# TEST
print("TEST AlgebraAplK")
S = SimplicialComplex([[0,1,2]],is_mutable=False)
A = AlgebraAplK(S)
print("Test check if the assignation is valid")
s0 = SimplicialAlgebraplN(0)
s0.inject_variables()
p0 = s0(2)
p1 = s0(1)
p2 = s0(3)
s1 = SimplicialAlgebraplN(1)
s1.inject_variables()
p01 = 2*t0+t1
p02 = 2*t0+3*t1
p12 = t0+3*t1
s2 = SimplicialAlgebraplN(2)
s2.inject_variables()
p012 = 2-t1+t2
values = {(0,): p0, (1,): p1, (2,): p2, (0,1): p01, (0,2): p02, (1,2): p12, (0,1,2): p012}
e1 = A._element_constructor_(values)
e1
#############################################