############################################################################### ## Characteristic cone given by equation ideal ## ## ## ## Author: Ayan Sengupta, Mousumi Mukherjee and Debasattam Pal ## ## Version: 1.0 ## ############################################################################### reset() import re import sys #code edited on 13 Jan 2017 # Written in SageMath version 7.4 and Python 2.7 ''' This program takes the generators of the cone and the partial difference equations of discrete autonomous nD system as inputs and gives a 'yes/no' certificate according to whether the given cone is a characteristic cone for the given system or not. First the dimension of the system, the number of partial difference equations and the number of cone generating vectors are taken as input. Then the difference equations and the vectors are given by the user. As output the program will give if the given cone is a characteristic cone or not. The dimension of the system determines the number of free variables of the system. According to the input given by the user, the free variables are defined in x, as x1,x2,x3... . So the variables that are valid during matrix input are x1, x1^-1, x2, x2^-1, x3, x3^-1 ... . To make the input more readable, and reduce typo errors, user can use y variables instead of x inverses i.e. user can give input in form of x1, y1, x2, y2, x3, y3 ..., where y1, y2, y3... are considered as x1^-1, x2^-1, x3^-1... Transposed column vectors are taken as input for generating the cone. For example -c1 = 0 1 0 is actually the vector [0 1 0]^T EXAMPLE 1: - Dimension of the system: 3 - Number of difference equations: 2 - Number of vectors for defining the cone: 3 - Enter the partial difference equations: - f1= x1^2+x2^2+x3^2-1 - f2= x2*x3-1 - please enter the transposed column vectors to define the cone: - v1= -1 0 0 - v2= 0 -1 0 - v3= 0 0 -1 - This is a characteristic cone EXAMPLE 1: -Dimension of the system: 3 -Number of partial difference equations: 2 -Number of vectors for generating the cone: 3 -Enter the partial difference equations: -f1= x1^5*x3 + x2^6 + x3^7*x2 -f2= x3*x1^-1 -please enter the transposed column vectors to define the cone: -v1= 1 0 0 -v2= 0 -1 0 -v3= 0 0 -1 -This is a characteristic cone Some predefined functions will give more information than just an 'yes/no' answer EXAMPLE: -remainder -[-t1*t2^2 - t1*t3^2 + t1, t3, t2, t1, t2, t3] -P_ring -Multivariate Polynomial Ring in x1, x2, x3, y1, y2, y3, t1, t2, t3 over Real Field with 53 bits of precision -difference_poly -[x1^2 + x2^2 + x3^2 - 1.00000000000000, x2*x3 - 1.00000000000000] -cone_poly -[-y1 + t1, -y2 + t2, -y3 + t3] -G -[x1 + t1*t2^2 + t1*t3^2 - t1, x2 - t3, x3 - t2, y1 - t1, y2 - t2, y3 - t3, t1^2*t2 + t1^2*t3^3 - t1^2*t3 + t3, t1^2*t3^4 - t1^2*t3^2 + t1^2 + t3^2, t2*t3 - 1.00000000000000] HOW TO RUN THE CODE IN SAGE: First open SAGE from the terminal. Then you can load a file in a Sage session with the command "load". example: if your file is in the current directory where you started sage then: load("your_file.sage") if your file is in another directory then : load("~/the_file_location/your_file.sage") ''' def cone(num,num_var,ring): '''This function takes the cone generating vectors as input. The vectors provided must be transposed column vectors. ''' print("\nplease enter the transposed column vectors to define the cone:") vector_poly = [0]*num for i in xrange(num): x_vector =1 y_vector =1 vector_monomial =1 print 'v' + str(i + 1) + '=', temp = raw_input().strip(" ").split(" ") vector = [int(num) for num in temp] #print vector if len(vector) != num_var: print("Wrong Vector length") break for p in xrange(len(vector)): if vector[p] > 0: x_vector = (ring.gens()[p])^vector[p] #print x_vector if vector[p] < 0: y_vector = (ring.gens()[num_var+p])^abs(vector[p]) #print y_vector vector_monomial = x_vector*y_vector*vector_monomial x_vector =1 y_vector =1 #print vector_monomial vector_poly[i] = ring.gens()[2*num_var + i] - vector_monomial #print vector_poly return vector_poly def polynomialRing(num, v_num): P = PolynomialRing( RR, ['x%s' %p for p in range(1, num + 1)]+['y%s' %p for p in range(1, num + 1)]+['t%s' %p for p in range(1, v_num + 1)], order='lex') P.inject_variables() return P def inputModification(num): '''It takes the partial differnece equations as inputs. The equations are given as a polynomial in x and x^-1. x^-1 terms are replaced by y variables. ''' poly = [0]*num for i in xrange(num): print 'f' + str(i + 1) + '=', temp = str(raw_input()) temp = ''.join( c for c in temp if c not in '{(})[]' ) replace_positions = [m.start() for m in re.finditer('\^-', temp)] a = list(temp) for p in xrange(len(replace_positions)): a[replace_positions[p]-2] = 'y' a[replace_positions[p]+1] = '' temp = ''.join(a) temp2 = re.sub('\^','**',temp) poly[i] = eval(temp2) return poly def variableSeparation(ring,num_var): '''some x, y, and t variables are the generators of the ring. This functions returns one list containing x, y variables and a set containing t variables ''' all_vars = ring.variable_names() laurent_constructors_temp= all_vars[0:num_var*2] laurent_constructors = [eval(laurent_constructors_temp[i]) for i in xrange(len(laurent_constructors_temp))] cone_constructors_temp = set(all_vars)-set(laurent_constructors_temp) cone_constructors = set([eval(list(cone_constructors_temp)[i]) for i in xrange(len(cone_constructors_temp))]) return laurent_constructors,cone_constructors def IdealandGroebnerBasis(difference_poly,laurent_constructors,cone_constructors,ring,num_var): ''' The ideal and groebner basis is calculated ''' laurent_poly = [0]*num_var for i in xrange(num_var): laurent_poly[i] = ring.gens()[i]*ring.gens()[i+num_var] - 1 I = ideal(laurent_poly)+ideal(cone_constructors) + ideal(difference_poly) groebner = I.groebner_basis() return groebner def remainderGrobnerBasis(laurent_constructors,groebner): ''' The remainder is calculated by dividing all the xi and yi by the Groebner bases ''' num = len(laurent_constructors) remainder_list = [laurent_constructors[i].reduce(groebner) for i in xrange(0,num)] return remainder_list def checkCharacteristicCone(remainder_list,cone_variables): '''This function checks whether the remainders are polynomials of t. ''' for i in xrange(len(remainder_list)): if not set(remainder_list[i].variables())<=cone_variables: return False return True ########################## main ############################# num_var = input("Dimension of the system: ") d_num = input("Number of partial difference equation(s): ") v_num = input("Number of vectors for generating the cone: ") P_ring = polynomialRing(num_var,v_num) print ("\nEnter the partial difference equations:") xy_variables , t_variables = variableSeparation(P_ring,num_var) difference_poly = inputModification(d_num) cone_poly = cone(v_num, num_var, P_ring) G = IdealandGroebnerBasis(difference_poly,xy_variables,cone_poly,P_ring,num_var) remainder = remainderGrobnerBasis(xy_variables,G) output = checkCharacteristicCone(remainder,t_variables) if output is True: print("\nThis is a characteristic cone") else: print("\nThis is NOT a characteristic cone")