2019-03-19 16:47:52 -0500 answered a question Converting real numbers to rational format Write R.ideal( instead of Ideal( which tries to guess the ring that you want (and gets it wrong here). 2019-03-19 10:25:50 -0500 commented question Major index of skew-SYT Can you give a reference to a (mathematical) definition? 2019-03-18 13:16:58 -0500 answered a question Set a REAL function in SAGE A simple alternative to 'setting it to be real' is to make a substitution such as this one: whatever.subs({conjugate(H) : H, real_part(H) : H, imag_part(H) : 0})  It depends on what kind of expressions you want to get rid of. 2019-03-17 17:29:43 -0500 answered a question Factoring modulo prime ideal Nonzero prime ideals in a Dedekind domain – such as $\mathfrak{p} = (1-\sqrt{-2})$ in $O_K$ – are maximal, so the quotient is a field, called the residue field in Sage: K. = QuadraticField(-2) OK = K.ring_of_integers() P = OK.ideal(1-a) OK_mod_P = P.residue_field()  Then you can factor polynomials over $O_K/\mathfrak{p}$, e.g. as follows: sage: R. = PolynomialRing(ZZ) sage: f = x^3 + 5 sage: f_mod_P = f.change_ring(OK_mod_P); f_mod_P x^3 + 2 sage: f_mod_P.factor() (x + 2)^3  You can also lift the factorization and double-check the result (e.g. for fun): sage: S. = PolynomialRing(OK) sage: factorization_lift = prod(S(g)^e for (g,e) in f_mod_P.factor()); factorization_lift y^3 + 6*y^2 + 12*y + 8 sage: S(f) - factorization_lift -6*y^2 - 12*y - 3 sage: S(f) - factorization_lift in S.ideal(P) True  2019-03-16 04:11:56 -0500 commented answer transformation matrix for variable matrix of given jordan type You're welcome! :) 2019-03-14 09:25:34 -0500 marked best answer Free nilpotent Lie algebra of step 3 with 11 generators This works: S = LieAlgebra(QQ, 10, step=3, names='X', naming='linear')  and this works: S = LieAlgebra(QQ, 11, step=2, names='X', naming='linear')  but this doesn't work: S = LieAlgebra(QQ, 11, step=3, names='X', naming='linear')  though it should? Traceback: --------------------------------------------------------------------------- AttributeError Traceback (most recent call last) in () ----> 1 S = LieAlgebra(QQ, Integer(11), step=Integer(3), names='X', naming='linear') /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/misc/classcall_metaclass.pyx in sage.misc.classcall_metaclass.ClasscallMetaclass.__call__ (build/cythonized/sage/misc/classcall_metaclass.c:1640)() 327 """ 328 if cls.classcall is not None: --> 329 return cls.classcall(cls, *args, **kwds) 330 else: 331 # Fast version of type.__call__(cls, *args, **kwds) /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/algebras/lie_algebras/lie_algebra.pyc in __classcall_private__(cls, R, arg0, arg1, names, index_set, abelian, nilpotent, category, **kwds) 439 from sage.algebras.lie_algebras.nilpotent_lie_algebra import FreeNilpotentLieAlgebra 440 del kwds["step"] --> 441 return FreeNilpotentLieAlgebra(R, arg1, step, names=names, **kwds) 442 elif nilpotent: 443 raise ValueError("free nilpotent Lie algebras must have a" /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/misc/classcall_metaclass.pyx in sage.misc.classcall_metaclass.ClasscallMetaclass.__call__ (build/cythonized/sage/misc/classcall_metaclass.c:1640)() 327 """ 328 if cls.classcall is not None: --> 329 return cls.classcall(cls, *args, **kwds) 330 else: 331 # Fast version of type.__call__(cls, *args, **kwds) /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/algebras/lie_algebras/nilpotent_lie_algebra.pyc in __classcall_private__(cls, R, r, s, names, naming, category, **kwds) 325 return super(FreeNilpotentLieAlgebra, cls).__classcall__( 326 cls, R,r, s, names=tuple(names), naming=naming, --> 327 category=category, **kwds) 328 329 def __init__(self, R, r, s, names, naming, category, **kwds): /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/misc/cachefunc.pyx in sage.misc.cachefunc.CachedFunction.__call__ (build/cythonized/sage/misc/cachefunc.c:6301)() 1003 return self.cache[k] 1004 except KeyError: -> 1005 w = self.f(*args, **kwds) 1006 self.cache[k] = w 1007 return w /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/structure/unique_representation.py in __classcall__(cls, *args, **options) 1025 True 1026 """ -> 1027 instance = typecall(cls, *args, **options) 1028 assert isinstance( instance, cls ) 1029 if instance.__class__.__reduce__ == CachedRepresentation.__reduce__: /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/misc/classcall_metaclass.pyx in sage.misc.classcall_metaclass.typecall (build/cythonized/sage/misc/classcall_metaclass.c:2090)() 494 TypeError: Argument 'cls' has incorrect type (expected type, got classobj) 495 """ --> 496 return (type).tp_call(cls, args, kwds) 497 498 # Class for timing:: /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/algebras/lie_algebras/nilpotent_lie_algebra.pyc in __init__(self, R, r, s, names, naming, category, **kwds) 398 for X_ind, X in basis_by_deg[dx]: 399 for Y_ind, Y in basis_by_deg[dy]: --> 400 Z = L[X, Y] 401 if not Z.is_zero(): 402 s_coeff[(X_ind, Y_ind)] = {W_ind: Z[W.leading_support()] /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/algebras/lie_algebras/lie_algebra.pyc in __getitem__(self, x) 573 return x[1].ideal(x[0]) 574 # Otherwise it is the bracket of two elements --> 575 return self(x[0])._bracket_(self(x[1])) 576 return super(LieAlgebra, self).__getitem__(x) 577 /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/algebras/lie_algebras/lie_algebra_element.pyx in sage.algebras.lie_algebras.lie_algebra_element.FreeLieAlgebraElement._bracket_ (build/cythonized/sage/algebras/lie_algebras/lie_algebra_element.c:19132)() 1452 a, b = mr, ml 1453 cr = -cr -> 1454 for b_elt, coeff in self.parent()._rewrite_bracket(a, b).iteritems(): 1455 d[b_elt] = d.get(b_elt, zero) + cl * cr * coeff 1456 if d[b_elt] == zero: /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/misc/cachefunc.pyx in sage.misc.cachefunc.CachedMethodCaller.__call__ (build/cythonized/sage/misc/cachefunc.c:10824)() 1951 return cache[k] 1952 except KeyError: -> 1953 w = self._instance_call(*args, **kwds) 1954 cache[k] = w 1955 return w /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/misc/cachefunc.pyx in sage.misc.cachefunc.CachedMethodCaller._instance_call (build/cythonized/sage/misc/cachefunc.c:10280)() 1827 True 1828 """ -> 1829 return self.f(self._instance, *args, **kwds) 1830 1831 cdef fix_args_kwds(self, tuple args, dict kwds): /opt/sagemath-8.4/local/lib/python2.7/site-packages/sage/algebras/lie_algebras/free_lie_algebra.pyc in _rewrite_bracket(self, l, r) 690 # For a similar reason, we have b >= c. 691 # Compute the left summand --> 692 for m, inner_coeff in iteritems(self._rewrite_bracket(l._right, r)): 693 if l._left == m: 694 continue AttributeError: 'sage.algebras.lie_algebras.lie_algebra_element.Lie' object has no attribute '_right'  2019-03-14 09:25:34 -0500 received badge ● Scholar (source) 2019-03-14 09:25:29 -0500 edited answer Free nilpotent Lie algebra of step 3 with 11 generators This has been submitted as trac ticket #27018 and subsequently it was fixed. 2019-03-13 12:39:43 -0500 answered a question multiplicity of a point in a scheme Yes, this is a bug. Typing Q.multiplicity?? we see that it calls X.multiplicity(Q) and X.multiplicity?? shows that it changes the monomial ordering to negdegrevlex and then tries to use the interface to Singular. However the interface to Singular is broken for univariate polynomial rings with "local" monomial ordering (it always uses the global ordering instead); I reported this as trac ticket #27479. As a temporary workaround I guess you can add a variable and set it to zero: sage: A2. = AffineSpace(QQ,2) sage: X = A2.subscheme([x^1789+x,y]) sage: Q = X([0,0]) sage: Q.multiplicity() 1  2019-03-12 11:21:52 -0500 received badge ● Nice Answer (source) 2019-03-12 04:11:26 -0500 answered a question How to compute the sums of squares of elements of a quotient ring? You can do it like this (also simplifying the notation a bit): sage: R. = PolynomialRing(GF(3)) sage: S. = R.quotient(x^2 + 1) sage: sum(p^2 for p in S) 0  We can explain the result (for this choice of $f$) as follows: $p(i) = a+bi$ implies $p(i)^2 = a^2 - b^2 + 2abi$ and the sum runs over all $a$ and $b$, so the "real part" of the sum is $$\sum_{a,b} a^2 - b^2 = \sum_a a^2 - \sum_b b^2 = 0,$$ and the "imaginary part" is $$\sum_{a,b} 2ab = \sum_{a,b} ab + \sum_{a,b} ab = \sum_{a,b} ab + \sum_{a,b} (-a)b = \sum_{a,b} ab - ab = 0.$$ 2019-03-11 13:17:22 -0500 commented question Evaluating a derivative of an unknown function This does not quite do what you're asking, but is close: sage: from sage.symbolic.constants import Constant sage: l = Constant('l').expression() sage: g.diff(z)(l) D[0](g)(l)  2019-03-09 05:24:34 -0500 commented question transformation matrix for variable matrix of given jordan type The quotient ring is not so nice because it has zero divisors, but see my answer below. 2019-03-09 05:19:51 -0500 answered a question transformation matrix for variable matrix of given jordan type We can set up the problem as follows: n = 3 jordan_type = (2,1) assert sum(jordan_type) == n import itertools A_coeff_names = {(i,j) : 'a_{}_{}'.format(i,j) for (i,j) in itertools.combinations(range(n),2)} R = PolynomialRing(QQ, names=A_coeff_names.values()) A_coeff = {idx : R(A_coeff_names[idx]) for idx in A_coeff_names} A = Matrix(R, n, A_coeff) J = block_diagonal_matrix([jordan_block(0, s) for s in jordan_type], subdivide=False)  A necessary condition for $gAg^{-1} = J$ is that the rank of $A^k$ equals the rank of $J^k$ for $1 \leq k < n$. (And since $A$ has all zeros as eigenvalues, this will give the right Jordan type.) In particular it is necessary that the $(\operatorname{rank}(J^k)+1)$-minors of $A^k$ vanish; these are polynomial equations upon the coefficients of $A$: A_coeff_eqns = [] J_power = identity_matrix(QQ, n) A_power = identity_matrix(R, n) for k in range(n): r = J_power.rank() # need (r+1) x (r+1) minors to vanish A_coeff_eqns.extend([eqn for eqn in A_power.minors(r+1) if eqn != 0]) J_power *= J A_power *= A  We can solve these equations symbolically, substitute a solution into $A$, change to an exact ring (the fraction field of the polynomial ring in the new variables), calculate the Jordan form (to be sure), and change the names of the variables back to the originals where possible: A_coeff_symb = map(SR, A_coeff_names.values()) A_coeff_eqns_symb = map(SR, A_coeff_eqns) for sol in solve(A_coeff_eqns_symb, A_coeff_symb): # substitute symbolic solution A_new = A.change_ring(SR).apply_map(lambda x: x.subs(sol)) A_new_vars = list(set(sum([list(x.variables()) for x in A_new.list()], []))) # change ring to fraction field of polynomial ring in the new variables A_new_ring = PolynomialRing(QQ, names=A_new_vars).fraction_field() A_new_rat = A_new.change_ring(A_new_ring) if A_new_rat.jordan_form(subdivide=False) != J: continue # change names of variables back to the originals where possible A_final_subs = {A_new_rat[i,j] : A_coeff[(i,j)] for (i,j) in A_coeff if A_new_rat[i,j] in A_new_vars} A_final = A_new_rat.apply_map(lambda x: x.subs(A_final_subs)) print A_final print print A_final.jordan_form(subdivide=False, transformation=True)[1] print print '---' print  Output: [ 0 a_0_1 a_0_2] [ 0 0 0] [ 0 0 0] [ a_0_1 0 0] [ 0 1 1] [ 0 0 a_0_1/(-a_0_2)] --- [ 0 0 a_0_2] [ 0 0 a_1_2] [ 0 0 0] [a_0_2 0 1] [a_1_2 0 0] [ 0 1 0] ---  Beware that this assumes Sage can solve the polynomial equations symbolically, which may not be true for large n. To give Sage an easier time (though the solutions may become less pretty), you can replace the list of equations by a Groebner basis for the same ideal: A_coeff_eqns = R.ideal(A_coeff_eqns).groebner_basis()  2019-03-07 10:17:48 -0500 answered a question Cannot solve differential equation (Lane-Emden equation) numerically I don't know how to fix that code, but as an alternative you can use desolve_system_rk4: var('y0 y1 t') ode_rhs = [y1, -2/t*y1 - y0^(3/2)] points = desolve_system_rk4(ode_rhs,[y0,y1],ics=[0.1,1,0],ivar=t,end_points=10,step=0.1)  To plot y0 (i.e. $x$) against t: ty0_points = [ [i,j] for i,j,k in points] list_plot(ty0_points, plotjoined=True)  To plot the curve in the (y0,y1)-plane (together with the vector field): y0y1_points = [ [j,k] for i,j,k in points] list_plot(y0y1_points, plotjoined=True, color='red') + sum(arrow2d(p[1:], vector(p[1:]) + vector([eqn_rhs.subs({t : p[0], y0 : p[1], y1: p[2]}) for eqn_rhs in ode_rhs]).normalized()*0.03,arrowsize=1,color='blue') for p in points if p[1] >= 0)  The vector field only makes sense when y0 >= 0 (same for the ODE and its "solution" plotted above). I don't know what the solver does past that point; I wouldn't trust it. 2019-03-07 04:46:06 -0500 commented answer Graphs having highest second smallest laplacian eigen value from a collection You're welcome. You can vote and accept using the buttons to the left of the answer. 2019-03-06 01:45:34 -0500 commented answer Graphs having highest second smallest laplacian eigen value from a collection @Kuldeep you're right, I got the definition of algebraic connectivity wrong because Wikipedia was not clear. I updated both Wikipedia and my answer. 2019-03-05 11:27:22 -0500 answered a question Obtaining integers from a linear extension of a poset. This is not my field, but it seems you can do sage: map(lambda x: x.element, p) [1, 4, 2, 3] sage: sigma = Permutation(map(lambda x: x.element, p)) sage: sigma [1, 4, 2, 3] sage: sigma(2) 4  2019-03-05 10:00:57 -0500 answered a question Summation of simbolic variables As part of the symbolic summation you are trying to access e.g. alpha[i] where i is an indeterminate variable; that doesn't work. Instead you should use the ordinary summation. Also, indices start from 0 in Python and SageMath. So instead of e.g. sum(area[i] * sin(2*alpha[i])^2, i, 1, 5)  it should be e.g. sum(area[i] * sin(2*alpha[i])^2 for i in range(5))  Note also that you can do things like n = 5 alpha = vector([var('alpha_%d' % (k+1)) for k in range(n)])  to save yourself some typing. 2019-03-05 06:58:54 -0500 commented answer Graphs having highest second smallest laplacian eigen value from a collection Why do you say it's not correct? When I run G = max((g for g in graphs(8) if g.is_connected() and g.girth() == 5), key=algebraic_connectivity) (or more efficiently G = max((g for g in graphs.nauty_geng('8 -c') if g.girth() == 5), key=algebraic_connectivity)) (or the longer code with analogous modification) then it seems correct. 2019-03-05 03:13:06 -0500 edited question Linear code over a finite ring The command LinearCode(M)  works well on a Finite Field F, but doesn't work for a finite ring. Is there any way to construct a linear code over a finite ring using sage? In other words, could we construct a module over ZZ/nZZ? I try to understand the trac tickets #6452 and #19345 for constructing a linear code over a finite ring. Unfortunately, I failed to understand. 2019-03-04 03:48:50 -0500 answered a question Graphs having highest second smallest laplacian eigen value from a collection For convenience we can define algebraic_connectivity = lambda g: sorted(g.laplacian_matrix().eigenvalues())[1]  Then we can obtain at once G = max((g for g in graphs(7) if g.girth() == 4), key=algebraic_connectivity)  But this assumes there is only one maximum. If we want to be more on the safe side: my_graphs = [g for g in graphs(7) if g.girth() == 4] from collections import defaultdict acs = defaultdict(list) for g in my_graphs: ac = algebraic_connectivity(g) acs[ac].append(g) max_ac = max(acs.keys()) print 'There are {} graphs of (maximal) algebraic connectivity {}'.format(len(acs[max_ac]), max_ac) for g in acs[max_ac]: show(g)  Output: There are 1 graphs of (maximal) algebraic connectivity 3  2019-03-03 16:50:12 -0500 answered a question What software does Sage use to solve linear equations ? From the documentation of solve: algorithm - string (default: ‘maxima’);  Sage uses Maxima by default. Also, more precisely, the source code of solve includes: if len(s) == 0: # if Maxima's solve gave no solutions, try its to_poly_solve try: s = m.to_poly_solve(variables) except Exception: # if that gives an error, stick with no solutions s = [] if len(s) == 0: # if to_poly_solve gave no solutions, try use_grobner try: s = m.to_poly_solve(variables,'use_grobner=true') except Exception: # if that gives an error, stick with no solutions s = []  So try Maxima's solve and to_poly_solve, and try the latter with use_grobner=true. 2019-03-03 06:49:44 -0500 received badge ● Good Answer (source) 2019-03-03 03:37:38 -0500 commented question transformation matrix for variable matrix of given jordan type An alternative approach is to write $A = g^{-1}Jg$ with $g$ symbolic, and solve for $g$ such that $A$ is strictly upper triangular ($n(n+1)/2$ polynomial equations for the $n^2$ entries of $g$). Would that suffice, or not? What quotient ring did you want to work in? 2019-03-03 03:13:26 -0500 received badge ● Nice Answer (source) 2019-03-03 02:40:48 -0500 answered a question Iterator: amazing behavior The generator is defined correctly, but you're not using it in the right way. Expression count() creates a new instance of the generator; so the first time you call next() it will give zero. Instead, you should create it once, and then call next() several times on the resulting object, e.g.: sage: mycount = count() sage: [mycount.next() for i in range(18)] [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17] sage: mycount = count() sage: [next(mycount) for i in range(18)] [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17] sage: next(mycount) 18 sage: next(mycount) 19  2019-03-02 07:01:08 -0500 commented question transformation matrix for variable matrix of given jordan type Can you give a more explicit example? 2019-03-02 01:41:50 -0500 received badge ● Nice Answer (source) 2019-03-01 10:38:40 -0500 answered a question Please can you advise how I can account allocate nominal account codes This forum is for questions about the SageMath mathematics software, not the Sage accounting software. For questions about the accounting software, see e.g. their support page. 2019-03-01 04:40:15 -0500 received badge ● Nice Answer (source) 2019-02-28 17:17:39 -0500 commented answer Eigenvalues over symbolic ring You're welcome! Yes, the characteristic polynomial can be computed in either ring, and when this can be applied it is never harmful. The matrix Matrix(SR,[[0,0,0,0,1], [1,0,0,0,1], [0,1,0,0,0], [0,0,1,0,0], [0,0,0,1,0]]) is 'multiplication by $\alpha$' w.r.t. the basis $1,\alpha,\alpha^2,\alpha^3,\alpha^4$ of $\mathbb{Q}(\alpha)$ where $\alpha$ is a root of $x^5 - x - 1$. So its characteristic polynomial is $x^5 - x - 1$ with Galois group $S_5$ which is not solvable, so the roots cannot be expressed in terms of radicals. This is a special case where there are no variables; then one should factor over $\mathbb{Q}$ instead of a polynomial ring. The polynomial won't factor in this case, and the subsequent attempt at a symbolic solution also won ... 2019-02-28 15:38:58 -0500 answered a question Eigenvalues over symbolic ring You can change the base ring of the matrix to the ring of rational functions (as you did), factor the characteristic polynomial, and set the factors equal to zero symbolically again (and solve). R = PolynomialRing(QQ, names='a,b,c,d') M_rat = M.change_ring(R.fraction_field()) S. = R[] f = M_rat.characteristic_polynomial(x) f_factors = f.factor() linear_factors = [(g,e) for (g,e) in f_factors if g.degree() == 1] rational_roots = [(solve(SR(g)==0,SR(x),solution_dict=True)[0][SR(x)],e) for (g,e) in linear_factors] nonlinear_factors = [(g,e) for (g,e) in f_factors if g.degree() > 1] nonrational_roots = sum([map(lambda eqn: (eqn.rhs(), e), solve(SR(g)==0,SR(x))) for (g,e) in nonlinear_factors], []) eigenvalues = rational_roots + nonrational_roots  This runs in 33 seconds on my machine. The result is [(3*(a^2 - b^2 + 2*b*c - c^2)*d + 6*a + 1, 1), (-3*(a^2 - b^2 - 2*a*c + c^2)*d + 6*b + 1, 1), (-3*(a^2 - 2*a*b + b^2 - c^2)*d + 6*c + 1, 1), ((3*a^2 + 2*a*b + 3*b^2 + 2*(a + b)*c + 3*c^2)*d - 2*a - 2*b - 2*c + 1, 2), ((3*a^2 - 3*b^2 - 2*b*c - 3*c^2)*d - 2*a + 1, 3), (-(3*a^2 - 3*b^2 + 2*a*c + 3*c^2)*d - 2*b + 1, 3), (-(3*a^2 + 2*a*b + 3*b^2 - 3*c^2)*d - 2*c + 1, 3), ((3*a^2 - 2*a*b + 3*b^2 - 2*(a + b)*c + 3*c^2)*d + 2*a + 2*b + 2*c - 4*sqrt((a^2*b^2 + (a^2 - a*b + b^2)*c^2 - (a^2*b + a*b^2)*c)*d^2 + a^2 - a*b + b^2 - (a + b)*c + c^2 + (a^2*b + a*b^2 + (a + b)*c^2 + (a^2 - 6*a*b + b^2)*c)*d) + 1, 1), ((3*a^2 - 2*a*b + 3*b^2 - 2*(a + b)*c + 3*c^2)*d + 2*a + 2*b + 2*c + 4*sqrt((a^2*b^2 + (a^2 - a*b + b^2)*c^2 - (a^2*b + a*b^2)*c)*d^2 + a^2 - a*b + b^2 - (a + b)*c + c^2 + (a^2*b + a*b^2 + (a + b)*c^2 + (a^2 - 6*a*b + b^2)*c)*d) + 1, 1)]  Maybe an algorithm like this could be included into Sage in the future. (The hard part is to detect when this can be applied.) 2019-02-28 11:01:30 -0500 received badge ● Good Answer (source) 2019-02-28 06:19:05 -0500 received badge ● Nice Answer (source) 2019-02-28 04:53:27 -0500 received badge ● Nice Answer (source)