Ask Your Question

Revision history [back]

click to hide/show revision 1
initial version

answered 2022-04-16 05:21:57 +0100

Max Alekseyev gravatar image

It appears to be an issue with symbolic calculation. However, there seems to be no point to deal with symbolic expressions here, while you in fact are working with polynomial/rational functions. It seems that it's more suitable to use LaurentPolynomialRing as a base class here. Also, it's better to leave dealing with matrices to Sage, not numpy. Here is an updated code:

import numpy as np
import sympy as sp
import itertools

R.<alpha_2, alpha_01, x_0, x_1, x_2> = LaurentPolynomialRing(QQ)
y=np.array([x_0, x_1, x_2])

f=np.zeros(shape=(3),dtype=object)
f[0]=alpha_01*(x_0 + x_1 + x_2 )*(x_1 + x_2 ) + alpha_2*x_1*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[1]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_2 ) + alpha_2*x_0*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[2]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_1 ) + alpha_2*x_0*x_1 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2

K=[4,1,1]
list2=[0,1,2]

F=np.zeros(shape=(3),dtype=object)     
for p in range(3):
    b=list2[p]
    g=(f[b])^K[b]/y[b]^(2*K[b])*(1/K[b])
    F[p]=(-g.subs({x_0:1}))    #F[i]=\hat{f}_i

A = Matrix(2,2,lambda i,j: diff(F[i+1],y[list2[j+1]]))

xd=det(A)

print("F[0]=",F[0])
print( " ---------------------------------- " )

print("xd=",xd)

dat=F[0]*xd
#dat=expand(dat)
print( " ---------------------------------- " )

print("DAT=F[0]*xd=",dat)

The coefficient of 1/x_1 (as a polynomial in other variables) in dat can be obtained as

sum(c*t for c,t in dat if t.degree(x_1)==-1)

It appears to be an issue with symbolic calculation. However, there seems to be no point to deal with symbolic expressions here, while you in fact are working with polynomial/rational functions. It seems that it's more suitable to use LaurentPolynomialRing as a base class here. Also, it's better to leave dealing with matrices to Sage, not numpy. Here is an updated code:

import numpy as np
import sympy as sp
import itertools

R.<alpha_2, alpha_01, x_0, x_1, x_2> = LaurentPolynomialRing(QQ)
y=np.array([x_0, x_1, x_2])

f=np.zeros(shape=(3),dtype=object)
#f=np.zeros(shape=(3),dtype=object)
f=[0]*3
f[0]=alpha_01*(x_0 + x_1 + x_2 )*(x_1 + x_2 ) + alpha_2*x_1*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[1]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_2 ) + alpha_2*x_0*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[2]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_1 ) + alpha_2*x_0*x_1 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2

K=[4,1,1]
list2=[0,1,2]

F=np.zeros(shape=(3),dtype=object)     
#F=np.zeros(shape=(3),dtype=object)     
F=[0]*3
for p in range(3):
    b=list2[p]
    g=(f[b])^K[b]/y[b]^(2*K[b])*(1/K[b])
    F[p]=(-g.subs({x_0:1}))    #F[i]=\hat{f}_i

A = Matrix(2,2,lambda i,j: diff(F[i+1],y[list2[j+1]]))

xd=det(A)

print("F[0]=",F[0])
print( " ---------------------------------- " )

print("xd=",xd)

dat=F[0]*xd
#dat=expand(dat)
print( " ---------------------------------- " )

print("DAT=F[0]*xd=",dat)

The coefficient of 1/x_1 (as a polynomial in other variables) in dat can be obtained as

sum(c*t for c,t in dat if t.degree(x_1)==-1)

It appears to be an issue with symbolic calculation. However, there seems to be no point to deal with symbolic expressions here, while you in fact are working with polynomial/rational functions. It seems that it's more suitable to use LaurentPolynomialRing as a base class here. Also, it's better to leave dealing with matrices to Sage, not numpy. Here is an updated code:

R.<alpha_2, alpha_01, x_0, x_1, x_2> = LaurentPolynomialRing(QQ)
y=np.array([x_0, y=[x_0, x_1, x_2])
x_2]

#f=np.zeros(shape=(3),dtype=object)
f=[0]*3
f[0]=alpha_01*(x_0 + x_1 + x_2 )*(x_1 + x_2 ) + alpha_2*x_1*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[1]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_2 ) + alpha_2*x_0*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[2]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_1 ) + alpha_2*x_0*x_1 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2

K=[4,1,1]
list2=[0,1,2]

#F=np.zeros(shape=(3),dtype=object)     
F=[0]*3
for p in range(3):
    b=list2[p]
    g=(f[b])^K[b]/y[b]^(2*K[b])*(1/K[b])
    F[p]=(-g.subs({x_0:1}))    #F[i]=\hat{f}_i

A = Matrix(2,2,lambda i,j: diff(F[i+1],y[list2[j+1]]))

xd=det(A)

print("F[0]=",F[0])
print( " ---------------------------------- " )

print("xd=",xd)

dat=F[0]*xd
#dat=expand(dat)
print( " ---------------------------------- " )

print("DAT=F[0]*xd=",dat)

The coefficient of 1/x_1 (as a polynomial in other variables) in dat can be obtained as

sum(c*t for c,t in dat if t.degree(x_1)==-1)

It appears to be an issue with symbolic calculation. However, there seems to be no point to deal with symbolic expressions here, while you in fact are working with polynomial/rational functions. It seems that it's more suitable to use LaurentPolynomialRing as a base class here. Also, it's better to leave dealing with matrices to Sage, not numpy. Here is an updated code:

R.<alpha_2, alpha_01, x_0, x_1, x_2> = LaurentPolynomialRing(QQ)
y=[x_0, x_1, x_2]

#f=np.zeros(shape=(3),dtype=object)
f=[0]*3
f[0]=alpha_01*(x_0 + x_1 + x_2 )*(x_1 + x_2 ) + alpha_2*x_1*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[1]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_2 ) + alpha_2*x_0*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[2]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_1 ) + alpha_2*x_0*x_1 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2

K=[4,1,1]
list2=[0,1,2]

#F=np.zeros(shape=(3),dtype=object)     
F=[0]*3
for p in range(3):
    b=list2[p]
    g=(f[b])^K[b]/y[b]^(2*K[b])*(1/K[b])
    F[p]=(-g.subs({x_0:1}))    #F[i]=\hat{f}_i

A = Matrix(2,2,lambda i,j: diff(F[i+1],y[list2[j+1]]))

xd=det(A)

print("F[0]=",F[0])
print( " ---------------------------------- " )

print("xd=",xd)

dat=F[0]*xd
#dat=expand(dat)
print( " ---------------------------------- " )

print("DAT=F[0]*xd=",dat)

The coefficient of 1/x_1 (as a polynomial in other variables) in dat can be obtained as

sum(c*t for c,t in dat if t.degree(x_1)==-1)
t.degree(x_1)==-1) * x_1

It appears to be an issue with symbolic calculation. However, there seems to be no point to deal with symbolic expressions here, while you in fact are working with polynomial/rational functions. It seems that it's more suitable to use LaurentPolynomialRing as a base class here. Also, it's better to leave dealing with matrices to Sage, not numpy. Here is an updated code:

R.<alpha_2, alpha_01, x_0, x_1, x_2> = LaurentPolynomialRing(QQ)
y=[x_0, x_1, x_2]

#f=np.zeros(shape=(3),dtype=object)
f=[0]*3
f[0]=alpha_01*(x_0 + x_1 + x_2 )*(x_1 + x_2 ) + alpha_2*x_1*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[1]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_2 ) + alpha_2*x_0*x_2 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2
f[2]=alpha_01*(x_0 + x_1 + x_2 )*(x_0 + x_1 ) + alpha_2*x_0*x_1 + (x_0*x_1 + x_0*x_2 + x_1*x_2 )*alpha_01 + 3*(x_0 + x_1 + x_2)^2

K=[4,1,1]
list2=[0,1,2]

#F=np.zeros(shape=(3),dtype=object)     
F=[0]*3
for p in range(3):
    b=list2[p]
    g=(f[b])^K[b]/y[b]^(2*K[b])*(1/K[b])
    F[p]=(-g.subs({x_0:1}))    #F[i]=\hat{f}_i

A = Matrix(2,2,lambda i,j: diff(F[i+1],y[list2[j+1]]))

xd=det(A)

print("F[0]=",F[0])
print( " ---------------------------------- " )

print("xd=",xd)

dat=F[0]*xd
#dat=expand(dat)
print( " ---------------------------------- " )

print("DAT=F[0]*xd=",dat)

The coefficient of 1/x_1 (as a polynomial in other variables) in dat can be obtained asas dat.coefficient(x_1^(-1)).

sum(c*t for c,t in dat if t.degree(x_1)==-1) * x_1
Copyright Sage, 2010. Some rights reserved under creative commons license. Content on this site is licensed under a Creative Commons Attribution Share Alike 3.0 license.
about | faq | help | privacy policy | terms of service
Powered by Askbot version 0.7.59