RectasDef
system:sage


{{{id=1|
%time
def Delta(lista,B):
    res=B(1)
    l=copy(lista)[:-1]
    while len(l)>0:
        res=res*B(l)
        l=l[:-1]
    return res

def clean_conj(lista,F):
    n=len(lista)
    if n%2==0:
        print "Not conjugate to a generator"
        return None
    m=ZZ((n-1)/2)
    a=lista[m]
    b=lista[m+1:]
    if F(b)*F(lista[:m])!=F(1):
        print "Not conjugate to a generator"
        return None
    return [a,b]        
  
def list_of_relations_2(orden0,wiring0):
    wiring=[[B(_[0]),_[1]]for _ in wiring0]
    res1=[]
    trenza=B(1)
    orden=copy(orden0)
    for cruce in wiring:
        u,v=cruce
        orden=(u^-1).permutation()*orden
        trenza=trenza*u
        w=[orden.inverse()(_) for _ in v]
        res1.append([trenza,v,w,orden])
        trw=Delta(w,B)
        trenza=trenza*trw
        orden=trw.permutation().inverse()*orden
    final=[]
    for aa in res1:
        lis=[]
        for i in aa[2]:
            accion=[_.sign()*orden0(_.abs()) for _ in (F([i])*aa[0]^-1).Tietze()]
            listas=clean_conj(accion,F)
            if listas==None:
                print 'Problems'
                return None
            else:
                lis+=[listas]
        datos1=[[_[0] for _ in lis],[_[1] for _ in lis]]
        final+=[datos1]
    return final
    
 
def CommProd1(a,l1,F,LR,XX):
    ll1=list(F(l1).Tietze())
    if len(ll1)==0:
        return M(0)
    return XX[a,ll1[0]]+tt([ll1[0]],LR)*CommProd1(a,ll1[1:],F,LR,XX)


def CommCyclic1(L,F,LR,XX,dir='RL'):
    res=[]
    r=len(L[0])
    comb=L[0]
    mn=comb.index(min(comb))
    L1=[v[mn:]+v[:mn] for v in L]
    LF=[F([-u for u in reversed(L1[1][j])]+[L1[0][j]]+L1[1][j]) for j in range(r)]
    for j in range(1,r):
        u=F(L1[1][j])
        if dir=='RL':
            LFp=[_ for _ in reversed(LF[j+1:]+LF[:j])]
        elif dir=='LR':
            LFp=LF[j+1:]+LF[:j]
        w=u*prod(LFp)*u^-1
        res+=[CommProd1(L1[0][j],list(w.Tietze()),F,LR,XX)]
    return res
    
def tt(i,ring):
    return prod([LR.gen(ZZ(j).abs()-1)^(ZZ(j).sign()) for j in i])
    
def MatTrunc(mat):
    A=mat
    m=A.nrows()
    return Matrix(m,[v.truncate(2) for v in A.list()])
    
    
def EscalonarM2(matriz,n):
    U=matriz
    U0=Matrix(U.nrows(),[v.constant_coefficient() for v in U.list()])
    A,B=U0.echelon_form(transformation=True)
    U1=B*U
    Apivot=[]
    j=0
    for i in range(A.rank()):
        while A[i,j]==0:
            j=j+1
        Apivot.append(j)
    for i in range(A.rank()):
        j=Apivot[i]
        U1.rescale_row(i,U1[i,j]^-1)
        for k in range(i)+range(i+1,A.nrows()):
            U1.add_multiple_of_row(k,i,-U1[k,j])
    SustNeg=identity_matrix(T0,n*(n-1)/2)
    for i in range(A.rank()):
        for  j in range(n*(n-1)/2):
            SustNeg[j,Apivot[i]]=-U1[i,j].truncate(2).polynomial()
    SustNeg=SustNeg.delete_rows(Apivot)
    return SustNeg

def ImagenMorfismo(L,XX,YY,S,LRv,dicv,inc):
    i,j=L
    res=XX[L].change_ring(S)
    ai=tt([i],LRv).subs(dicv)
    aj=tt([j],LRv).subs(dicv)
    aij=(ai*aj).truncate(2)
    bij=(ai-aij)
    cij=(aj-aij)
    res+=sum([(bij*YY[i,u,v]-cij*YY[j,u,v])*XX[u,v].change_ring(S) for (u,v) in inc])
    return res
///
CPU time: 0.01 s,  Wall time: 0.01 s
}}}

{{{id=2|
%time
ordenPos=Permutation([1,3,7,5,8,11,4,10,9,6,2])
wiringPos=[[(), [10, 9]], [(), [5, 8]], [(), [5, 11, 4, 9]], [(), [5, 10]], [(),
[5, 6]], [(), [5, 2]], [(-7, -8), [7, 8]], [(), [7, 9]], [(), [7, 4]],
[(), [7, 10, 6]], [(), [7, 11, 2]], [(-8, -4, 7), [6, 11]], [(-7, -5,
-6, -4, 8), [3, 8, 11]], [(), [3, 4]], [(), [3, 10]], [(), [3, 9, 2,
6]], [(3,), [8, 10]], [(4, 5, 2, 3, 4, -2), [1, 8, 4, 6]], [(), [1, 11,
10]], [(), [1, 2]], [(), [1, 9]], [(-3, -4), [8, 2]]]
ordenNeg=Permutation([1,7,5,3,4,10,11,8,2,6,9])
wiringNeg=[[(), [3, 4]], [(), [3, 10]], [(), [3, 11, 8]], [(), [3, 2, 6, 9]],
[(-5, -6, -7, -4, -5), [5, 8]], [(), [5, 11, 4, 9]], [(), [5, 10]], [(),
[5, 6]], [(), [5, 2]], [(-4, 3, 6), [11, 6]], [(-6, -5), [9, 10]], [(4,
5, 4, 8, -7, 6), [2, 8]], [(7,), [7, 4]], [(), [7, 10, 6]], [(), [7,
8]], [(), [7, 9]], [(), [7, 2, 11]], [(-4, 5, -6, -3), [1, 4, 8, 6]],
[(), [1, 9]], [(), [1, 11, 10]], [(), [1, 2]], [(-2, -4), [8, 10]]]
n=11
B=BraidGroup(n)
F=FreeGroup(n)
///
CPU time: 0.01 s,  Wall time: 0.01 s
}}}

{{{id=3|
%time
GrupoListaPos=list_of_relations_2(ordenPos,wiringPos)
GrupoListaNeg=list_of_relations_2(ordenNeg,wiringNeg)
combPos=[sorted(v[0]) for v in GrupoListaPos]
combNeg=[sorted(v[0]) for v in GrupoListaNeg]
print "Do combinatorics coincide? ",sorted(combPos)==sorted(combNeg)
# Construct the list inc as a basis of commutators in G_1/G_2
comb=[]
for v in combPos:
    for j in v[1:]:
        comb.append([v[0],j])
tot=[tuple(sorted(_.list())) for _ in Subsets([1..n],2)] 
inc=[_ for _ in tot if list(_) not in comb]

# The variables of the ring S are the unknowns of the final integral equation
xx=var(['v%da%db%d'% (k,i,j) for k in [1..n] for (i,j) in inc])
S=PolynomialRing(ZZ,xx)

# The dictionnary YY relates the index of the unknowns with the morphism of groups.
YY={}
for k in [1..n]:
    for l in range(len(inc)):
        YY[k,inc[l][0],inc[l][1]]=S.gen((k-1)*len(inc)+l)

LR=LaurentPolynomialRing(ZZ,'t',n)
LRv=LR.change_ring(S)
R=LR.polynomial_ring()
T=PowerSeriesRing(S,'s',num_gens=n,default_prec=2)
T0=T.change_ring(ZZ)
HomLRT0=R.hom([1+v for v in T0.gens()],codomain=T0,check=True)
dic={v:HomLRT0(v) for v in LR.gens()}
dicv={v:HomLRT0(v).change_ring(S) for v in LRv.gens()}

#We work with several modules, with basis the set of commutators. The elements XX[i,j] gives the commutator (any sign and order) with respect to the ring of Laurent polynomials

M=FreeModule(LR,n*(n-1)/2)
#MT=M.change_ring(T0)
#Mv=M.change_ring(T)
#XX={(i,j):M.gen(n*(i-1)-i*(i+1)/2+j-1) for i in range(1,n) for j in range(i+1,n+1)}
XX={tot[i]:M.gen(i) for i in range(n*(n-1)/2)}
for i in range(1,n+1):
    XX[i,i]=M(0)
    XX[-i,-i]=M(0)
    XX[i,-i]=M(0)
    XX[-i,i]=M(0)
    for j in range(i+1,n+1):
        XX[-i,j]=-tt([i],LR)^-1*XX[i,j]
        XX[i,-j]=-tt([j],LR)^-1*XX[i,j]
        XX[-i,-j]=tt([i],LR)^-1*tt([j],LR)^-1*XX[i,j]
        XX[j,i]=-XX[i,j]
        XX[-j,i]=-XX[i,-j]
        XX[j,-i]=-XX[-i,j]
        XX[-j,-i]=-XX[-i,-j]

# Translate the relations in the groups as relations for the Alexander invariants
relAlexPos=[]
for v in GrupoListaPos:
    relAlexPos+=CommCyclic1(v,F,LR,XX,dir='LR')
print "Alexander invariant for first group"
relAlexNeg=[]
for v in GrupoListaNeg:
    relAlexNeg+=CommCyclic1(v,F,LR,XX,dir='LR')
nrels=len(relAlexNeg)
print "Alexander invariant for the second group"

# The above relations are written in the power series ring truncated at level 2
relAlexSeriePos=Matrix([vector(w.subs(dic) for w in v.list()) for v in relAlexPos])
relAlexSerieNeg=Matrix([vector(w.subs(dic) for w in v.list()) for v in relAlexNeg])

#This matrix codifies each commutator in terms of the basis inc
SustNeg=EscalonarM2(relAlexSerieNeg,n)
print "The relations are used to write every one in terms of the basis inc"

# Image of each commutator of the first group in the second group
TotImagenMorfismo=[]
for u in tot:
    v=ImagenMorfismo(u,XX,YY,S,LRv,dicv,inc)
    TotImagenMorfismo.append(v)
    if u[-1]==n:
        print "Finished pairs with line ", u[0]
        
relAlexVarPos=[]
for ser in relAlexSeriePos:
    vct=0
    for j in range(n*(n-1)/2):
        vct+=ser[j].change_ring(S)*TotImagenMorfismo[j]
    relAlexVarPos.append(vct)
relAlexVarPos=MatTrunc(Matrix(relAlexVarPos))
relAlexVarPos=MatTrunc(relAlexVarPos*SustNeg.change_ring(T).transpose())
print "Images of relations of first group in terms of the basis inc in the second group, done"

relZPos=[]
for prueba in relAlexVarPos:
    relZPos.append(vector(flatten([[pr1.derivative(v).constant_coefficient() for v in T.gens()] for pr1 in prueba])))
relZPos=Matrix(nrels,relZPos)
print "Integral matrix from relations"


#Adding Jacobi relations in the base inc and in the associated one in the next level
JCB=[]
for i in range(1,n-1):
    for j in range(i+1,n):
        for k in range(j+1,n+1):
            JCB+=[T0.gen(i-1)*XX[j,k].change_ring(T0)+T0.gen(j-1)*XX[k,i].change_ring(T0)+T0.gen(k-1)*XX[i,j].change_ring(T0)]
JCB=Matrix(JCB)
JCB1=MatTrunc(JCB*SustNeg.transpose())
JCB2=Matrix(JCB.nrows(),flatten([[v.derivative(w).constant_coefficient() for w in T0.gens()] for v in JCB1.list()]))
print "Simplified Jacobi relations as integral matrix done"

#Eliminating generators from Jacobi relations (if no torsion)
SF,U,V=JCB2.smith_form()
Jdiag=[SF[v,v] for v in range(SF.rank())]
if Set(Jdiag)!=Set([1]):
    print "Torsion at level 2"
print "Smith form of Jacobi relations, done", SF.ncols()-SF.rank()," generators left"

# Reduce the relations taking into account Jacobi relations (no torsion!)
eqs=(relZPos*V).matrix_from_columns([SF.rank()..JCB2.ncols()-1])
eqs1=[_ for _ in list(Set(eqs.list())) if _!=0]
print len(eqs1)," equations in ",len(inc)*n," unknowns."

# Put the system in matrix for, print the diagonal of the Smith form and a particular solution
Aeq=Matrix([vector(eq.monomial_coefficient(v) for v in S.gens()) for eq in eqs1])
Beq=vector(-eq.constant_coefficient() for eq in eqs1)
print "Ranks: ", Aeq.rank(),Aeq.augment(Beq).rank()
U1,U2,U3=Aeq.smith_form()
diageq=[U1[v,v] for v in range(U1.rank())]
diageq1=[_ for _ in diageq if _!=1]
print "Ones in the diagonal of the smith form: ",diageq.count(1)
print "Rest of the diagonal of the smith form: ",diageq1
print "Non integers in the Particular Solution:", [_ for _ in U1.solve_right(U2*Beq) if _ not in ZZ]
///
Do combinatorics coincide?  True
Alexander invariant for first group
Alexander invariant for the second group
The relations are used to write every one in terms of the basis inc
Finished pairs with line  1
Finished pairs with line  2
Finished pairs with line  3
Finished pairs with line  4
Finished pairs with line  5
Finished pairs with line  6
Finished pairs with line  7
Finished pairs with line  8
Finished pairs with line  9
Finished pairs with line  10
Images of relations of first group in terms of the basis inc in the second group, done
Integral matrix from relations
Simplified Jacobi relations as integral matrix done
Smith form of Jacobi relations, done 91  generators left
930  equations in  253  unknowns.
Ranks:  241 241
Ones in the diagonal of the smith form:  240
Rest of the diagonal of the smith form:  [5]
Non integers in the Particular Solution: [2/5]
CPU time: 314.02 s,  Wall time: 314.14 s
}}}

{{{id=4|
ordenPos1=copy(ordenPos)
wiringPos1=[[tuple(-i for i in _[0]),_[1]] for _ in wiringPos]
///
}}}

{{{id=5|
%time
GrupoListaPos1=list_of_relations_2(ordenPos1,wiringPos1)
combPos1=[sorted(v[0]) for v in GrupoListaPos1]
print "Do combinatorics coincide? ",sorted(combPos1)==sorted(combNeg)
# Construct the list inc as a basis of commutators in G_1/G_2

# Translate the relations in the groups as relations for the Alexander invariants
relAlexPos1=[]
for v in GrupoListaPos1:
    relAlexPos1+=CommCyclic1(v,F,LR,XX,dir='LR')
print "Alexander invariant for first group"

# The above relations are written in the power series ring truncated at level 2
relAlexSeriePos1=Matrix([vector(w.subs(dic) for w in v.list()) for v in relAlexPos1])


        
relAlexVarPos1=[]
for ser in relAlexSeriePos1:
    vct=0
    for j in range(n*(n-1)/2):
        vct+=ser[j].change_ring(S)*TotImagenMorfismo[j]
    relAlexVarPos1.append(vct)
relAlexVarPos1=MatTrunc(Matrix(relAlexVarPos1))
relAlexVarPos1=MatTrunc(relAlexVarPos1*SustNeg.change_ring(T).transpose())
print "Images of relations of first group in terms of the basis inc in the second group, done"

relZPos1=[]
for prueba in relAlexVarPos1:
    relZPos1.append(vector(flatten([[pr1.derivative(v).constant_coefficient() for v in T.gens()] for pr1 in prueba])))
relZPos1=Matrix(nrels,relZPos1)
print "Integral matrix from relations"


# Reduce the relations taking into account Jacobi relations (no torsion!)
eqs_1=(relZPos1*V).matrix_from_columns([SF.rank()..JCB2.ncols()-1])
eqs1_1=[_ for _ in list(Set(eqs_1.list())) if _!=0]
print len(eqs1_1)," equations in ",len(inc)*n," unknowns."

# Put the system in matrix for, print the diagonal of the Smith form and a particular solution
Aeq1=Matrix([vector(eq.monomial_coefficient(v) for v in S.gens()) for eq in eqs1_1])
Beq1=vector(-eq.constant_coefficient() for eq in eqs1_1)
print "Ranks: ", Aeq1.rank(),Aeq1.augment(Beq1).rank()
U1_1,U2_1,U3_1=Aeq1.smith_form()
diageq_1=[U1_1[v,v] for v in range(U1_1.rank())]
diageq1_1=[_ for _ in diageq_1 if _!=1]
print "Ones in the diagonal of the smith form: ",diageq_1.count(1)
print "Rest of the diagonal of the smith form: ",diageq1_1
print "Non integers in the Particular Solution:", [_ for _ in U1_1.solve_right(U2_1*Beq1) if _ not in ZZ]
///
Do combinatorics coincide?  True
Alexander invariant for first group
Images of relations of first group in terms of the basis inc in the second group, done
Integral matrix from relations
930  equations in  253  unknowns.
Ranks:  241 241
Ones in the diagonal of the smith form:  240
Rest of the diagonal of the smith form:  [5]
Non integers in the Particular Solution: [1/5]
CPU time: 44.18 s,  Wall time: 44.19 s
}}}

{{{id=6|

///
}}}