In a for loop I generate a list of functions, then I take their derivatives with respect to different variables. For example, in the follwoing, Y[0,0] is a 'MPolynomial_polydict', and is equal to
x_0_1*x_0_1 + x_0_0*x_1_1 + x_0_2*x_0_2 + x_1_2*x_1_2 + x_0_0*x_2_2 + x_1_1*x_2_2 - 1.00000000000000
and list_of_variables[0] = x_0_1. When I take the derivative of the mentioned polynomial with respect to the mentioned variable I should get 2*x_0_1, but instead I get:
sage: derivative(Y[0,0],list_of_variables[0])
x_0_1
When I try this manually, it gives me the correct answer though:
sage: derivative(x_0_1*x_0_1 + x_0_0*x_1_1 + x_0_2*x_0_2 + x_1_2*x_1_2 + x_0_0*x_2_2 + x_1_1*x_2_2 - 1.00000000000000,x_0_1)
2*x_0_1
So, for some reason it considers the two x_0_1's in the first term of my polynomial as independent variables when I construct them automatically. Any ideas on why this problem arises, and how to fix this?
Edit:
This is how I build my vaiables:
names = [ [[] for i in range(n)] for j in range(n) ]
for i in range(n):
for j in range(i+1):
names[i][j] = (SR('x_' + str(j) + '_' + str(i)))
names[j][i] = (SR('x_' + str(j) + '_' + str(i)))
R = PolynomialRing(CC,names[i][j] for i in range(n) for j in range(n)])
R.gens()
R.inject_variables()
More edit:
Here is a the complete code:
#########################################################################
# Build variables, and the matrix corresponding to it
#########################################################################
def build_variables(n):
names = [ [[] for i in range(n)] for j in range(n) ]
for i in range(n):
for j in range(i+1):
names[i][j] = (SR('x_' + str(j) + '_' + str(i)))
names[j][i] = (SR('x_' + str(j) + '_' + str(i)))
return(names)
#########################################################################
# Define the function f that maps a matrix to the coefficients of its characteristic polynomial
#########################################################################
def CharPoly(Mat):
X = matrix(Mat)
n = X.ncols()
C_X = X.characteristic_polynomial()
Y = []
for i in range(n):
Y.append(C_X[i])
return(Y)
############################################################################
# This solves that lambda problem
############################################################################
def lambda_siep(G,L,iter=100,epsilon = .1):
# G is any graph on n vertices
# L is the list of n desired distinct eigenvalues
# m is the number of itterations of the Newton's method
# epsilon: the off-diagonal entries will be equal to epsilon
n = G.order()
my_variables = build_variables(n)
R = PolynomialRing(CC,[my_variables[i][j] for i in range(n) for j in range(n)])
R.gens()
R.inject_variables()
X = [ [ R.gens()[n*i+j] for i in range(n) ] for j in range(n) ]
Y = matrix(CharPoly(X)) - matrix(CharPoly(diagonal_matrix(L)))
J = matrix(R,n)
for i in range(n):
for j in range(n):
J[i,j] = derivative(Y[0][i],my_variables[j][j])
B = diagonal_matrix(L) + epsilon * G.adjacency_matrix()
count = 0
while count < iter:
T = [ B[i,j] for i in range(n) for j in range(n)]
C = (J(T)).solve_right(Y(T).transpose())
LC = list(C)
B = B - diagonal_matrix([LC[i][0] for i in range(n)])
count = count + 1
return(B)
