Ask Your Question

rburing's profile - activity

2019-07-21 06:05:17 -0500 answered a question Is there a way to temporarily turn off sage's type checking?

Personally I would try to find a way to avoid the symbolic ring (with var, etc.) altogether.

If you want to use the symbolic ring, then you have to do the conversion to a LinearFunction somehow.

Assuming this setup:

sage: var('x,y')
sage: expr = 2*x + y
sage: p = MixedIntegerLinearProgram(maximization=False)
sage: pvar = p.new_variable(real=True, nonnegative=True)

You can do it e.g. like this:

def linear_function_from_symbolic(expr, pvar):
    constant_coeff = pvar.mip().base_ring()(expr.subs(dict(zip(expr.variables(), [0]*len(expr.variables())))))
    return sum(expr.coefficient(v)*pvar[v] for v in expr.variables()) + constant_coeff

and use it as follows:

sage: linear_function_from_symbolic(expr, pvar)
2*x_0 + x_1

You can convert constraints too:

def linear_constraint_from_symbolic(expr, pvar):
    return expr.operator()(linear_function_from_symbolic(expr.lhs(), pvar),
                           linear_function_from_symbolic(expr.rhs(), pvar))

for example:

sage: linear_constraint_from_symbolic(x <= 1, pvar)
x_0 <= 1

It is reasonable that there is no coercion from the Symbolic Ring to Linear functions over RDF, because it would not be well-defined with respect to the ordering of variables: e.g. coercion of two functions of different variables would depend on the order in which you do it. The function above "suffers" from the same (but probably you don't care).

2019-07-17 14:43:45 -0500 received badge  Nice Answer (source)
2019-07-16 14:13:00 -0500 answered a question Plotting derivative of bump function

This is definitely a bug; I reported it as trac ticket #28207.

As a workaround, you can plot the derivative manually:

sage: f(x) = piecewise([((-oo, -1), 0), ((-1, 1), exp(-1/(1 - x^2))), ((1, oo), 0)])
sage: df(x) = piecewise([((-oo,-1), 0), ((-1,1), diff(exp(-1/(1 - x^2)), x)), ((1,oo), 0)])
sage: plot(df, (x,-2,2))

derivative of bump function

For linear combinations, you can use linearity of the derivative (as a workaround).

2019-07-16 13:29:49 -0500 answered a question How to define polynomial p(x_i, x_j) while x_i, x_j runs over available variables?

I suppose our multivariate polynomial ring R is defined as follows:

sage: R = PolynomialRing(QQ, 'x_', 3)
sage: x = R.gens()

Now it depends a little bit on how your bivariate polynomial is defined. Is it a Python function? A Sage callable symbolic expression? An element of a polynomial ring? An element of the symbolic ring?

If it is a Python function:

sage: p = lambda a,b: a^2 + b^2
sage: [p(x[i], x[j]) for i in range(len(x)) for j in range(len(x)) if i != j]

If it is a Sage callable expression:

sage: p(a,b) = a^2 + b^2
sage: [R(p(x[i], x[j])) for i in range(len(x)) for j in range(len(x)) if i != j]

If it is an element of the symbolic ring:

sage: var('a,b')
sage: p = a^2 + b^2
sage: [R(p(a=x[i], b=x[j])) for i in range(len(x)) for j in range(len(x)) if i != j]

If it is a polynomial in a polynomial ring:

sage: S.<a,b> = PolynomialRing(QQ)
sage: p = a^2 + b^2
sage: P = p.change_ring(R)
sage: [P(x[i],x[j]) for i in range(len(x)) for j in range(len(x)) if i != j]

All give the same result, a list of elements of R:

[x_0^2 + x_1^2,
 x_0^2 + x_2^2,
 x_0^2 + x_1^2,
 x_1^2 + x_2^2,
 x_0^2 + x_2^2,
 x_1^2 + x_2^2]

You can add other conditions in the list comprehension using and etc.

Also you can use itertools.product or itertools.combinations:

sage: p = lambda a,b: a^2 + b^2
sage: import itertools
sage: [p(x[i], x[j]) for (i,j) in itertools.product(range(len(x)),repeat=2) if i != j]

And you don't have to use indices if your condition doesn't depend on it:

sage: p = lambda a,b: a^2 + b^2
sage: [p(z,w) for (z,w) in itertools.product(x,repeat=2) if z != w]
2019-07-16 12:39:58 -0500 answered a question Create Matrix in RREF with indeterminates

Here's a way. First we define our variables:

k = 2
n = 3
R = PolynomialRing(QQ, 'x',(n-k)*k)

Then we introduce the following function:

def indet_rref(k, n, pivots, ring):
    assert len(pivots) == k and sorted(pivots) == list(pivots)
    M = Matrix(ring, k, n)
    taken = []
    for (row, pivot) in enumerate(pivots):
        M[row, pivot] = 1
        taken.extend((i, pivot) for i in range(k))
        taken.extend((row,j) for j in range(pivot))
    indet_indices = [(i,j) for i in range(k) for j in range(n) if not (i,j) in taken]
    for (idx, y) in zip(indet_indices, ring.gens()):
        M[idx] = y
    return M

The pivots argument is a list of k column indices (one per row), and the generators of ring are taken as the indeterminates.

For example:

sage: indet_rref(k, n, range(k), R)
[ 1  0 x0]
[ 0  1 x1]

We can test all the combinations:

import itertools
for pivots in itertools.combinations(range(n), k):
    M = indet_rref(k, n, pivots, R)
    assert M.rref() == M
    assert M.rank() == k
    print M


[ 1  0 x0]
[ 0  1 x1]

[ 1 x0  0]
[ 0  0  1]

[0 1 0]
[0 0 1]
2019-07-11 06:33:52 -0500 received badge  Nice Answer (source)
2019-07-11 04:16:33 -0500 answered a question How to generate a nonsquare Vandermonde matrix?

Enter matrix.vandermonde?? to see the source code (3 lines long). We can extend the functionality:

def my_vandermonde(v, ncols=None, ring=None):
    def entries(i, j):
        return v[i]**j
    if ncols is None:
        ncols = len(v)
    return matrix(entries, nrows=len(v), ncols=ncols, ring=ring)

The relative inefficiency here (taking powers instead of repeatedly multiplying) is copied from the original.

Then you can do for instance:

sage: my_vandermonde([2,3,4], 4)
[ 1  2  4  8]
[ 1  3  9 27]
[ 1  4 16 64]
2019-07-10 19:49:07 -0500 received badge  Nice Answer (source)
2019-07-10 13:00:34 -0500 commented answer Defining manifolds in a systematic way.

That should work (and it does for me). What version of SageMath are you running?

2019-07-10 11:51:26 -0500 answered a question Defining manifolds in a systematic way.

See the documentation on charts, particularly the arguments coordinates and names.

For example, you can do the following:

sage: x = U.chart(names=tuple('x_%d' % i for i in range(2*p*q)))
sage: x[:]
(x_0, x_1, x_2, x_3)
sage: x[0]
2019-07-09 15:26:03 -0500 answered a question Does the set_embedding command have a bug?

Setting the combinatorial embedding does not immediately affect how the graph is plotted. However, you can call one of the layout algorithms that use this combinatorial embedding, such as planar, with save_pos=True.

G = Graph([[0,1,2,3,5,4],[[0,1],[0,2],[0,3],[0,4],[0,5]]])
G.set_embedding({0: [1,2,3,4,5], 1 : [0], 2 : [0], 3:[0], 4:[0], 5:[0]})
G.layout('planar', save_pos=True)

graph with embedding

Indeed we get clockwise 1, 2, 3, 4, 5.

Edit: Obviously not clockwise but counterclockwise. I opened trac ticket #28152 for this issue.

2019-07-09 14:51:50 -0500 commented question Does the set_embedding command have a bug?

Please include a full code sample, since G, verts and graphare currently not defined. Also include a full code sample of how you want to call this function, what you expect as output, and what you get instead. (It should be possible to run it e.g. on It is very difficult to help you otherwise.

2019-07-09 06:35:31 -0500 commented question Does the set_embedding command have a bug?

Please edit your question: paste the relevant code, select it and press the 101010 button to indent it by four spaces (to get the proper formatting).

2019-07-08 13:41:08 -0500 answered a question How to turn the function into expression?

Suppose we have an ordinary function

sage: def f(x,y,z):
....:     return x^2 - y + z^3

or (which is the same, but shorter)

sage: f = lambda x,y,z: x^2 - y + z^3

Evaluating the function at arguments of the desired type, we can work with the result, e.g. symbolically:

sage: var('x,y,z')
sage: f(x,y,z).coefficient(y)

or with polynomials:

sage: R.<x,y,z> = QQ[]
sage: f(x,y,z).coefficient({y : 1})

However, the latter variant does not work when f is instead a SageMath callable symbolic expression:

sage: f(x,y,z) = x^2 - y + z^3

because it will always output a symbolic expression. So in order to work with polynomials, you have to do:

sage: R(f(x,y,z)).coefficient({y : 1})

instead (or better: define an ordinary function, as above).

In my opinion, this is one of the reasons why having this notation for callable symbolic expressions is bad.

Side note: in SageMath it is often very reasonable to use expressions instead of functions. Many facilities are designed around this, accepting an expression and a list of variables, so the expression will be interpreted as a function of the specified variables. For example, to take a derivative with diff you can specify the variable with respect to which you want to differentiate.

2019-06-25 14:44:03 -0500 answered a question (Update) How to change/set variables?

Here's how you can do substitution into the coefficients of a differential form, in a simple example:

sage: M = Manifold(2, 'M')
sage: U = M.open_subset('U')
sage: c_xy.<x,y> = U.chart()
sage: eU = c_xy.frame()
sage: f = M.diff_form(2)
sage: var('a,b')
sage: f[eU, 0, 1] = a+b
sage: f.display(eU)
(a + b) dx/\dy
sage: f[eU, 0, 1] = f[eU, 0, 1].expr().subs({a : x^2, b : y^2})
sage: f.display(eU)
(x^2 + y^2) dx/\dy
2019-06-25 14:29:44 -0500 commented question (Update) How to change/set variables?

Also note the variables g[(i,j,k)] are not functions of x, so the differential forms have constant coefficients in this frame. That means "compute with symbolic coefficients first, substitute later" doesn't work (in this way) if you're doing any computations involving derivatives. There seem to be no derivatives here though, so I guess it's fine.

2019-06-25 14:06:04 -0500 commented question (Update) How to change/set variables?

Updating the original forms is pointless. What you want to do is a substitution in the components of CH.

2019-06-16 04:16:32 -0500 answered a question solve equation with two variables over RR

This is a bug; I reported it as trac ticket #27998.

2019-06-10 04:34:29 -0500 received badge  Good Answer (source)
2019-06-09 01:23:05 -0500 answered a question roots of polynomial

It is the multiplicity of the respective root. Try factoring the polynomial with f.factor() and you will see.

2019-06-09 01:15:13 -0500 answered a question roots of polynomial

It is the multiplicity of the respective root. Try factoring the polynomial with f.factor() and you will see.

2019-06-07 03:24:14 -0500 answered a question Problems and errors in solve an equation

This is a quartic polynomial in omega^2, so make the substitution and solve:

solve(omega_nf_eq.subs({omega : sqrt(omega_squared)}), omega_squared)

Take plusminus the square roots of these solutions to get the solutions for omega.

2019-06-06 03:13:09 -0500 received badge  Nice Answer (source)
2019-06-04 14:13:51 -0500 answered a question Adding elements to a set

The documentation on sets explains it well:

If X is a list, tuple, Python set, or X.is_finite() is True, this returns a wrapper around Python’s enumerated immutable frozenset type with extra functionality. Otherwise it returns a more formal wrapper.

If you need the functionality of mutable sets, use Python’s builtin set type.

So Set([a,b]) acts more like Python's frozenset (hence it has no method to add an element).

If you want a mutable set, just use Python's set type.

2019-06-01 06:13:32 -0500 answered a question How to handle elements of two different Galois fields simultaneously?

The modulus should be a polynomial over the prime field.

R.<x> = PolynomialRing(GF(2))
G = GF(2^8, name='x', modulus=x^8 + x^5 + x^3 + x + 1)
F = GF(2^3, name='x', modulus=x^3 + x^2 + 1)

for i in range(2^3):
    print G.fetch_int(i).integer_representation(), '=', G.fetch_int(i)
    print F.fetch_int(i).integer_representation(), '=', F.fetch_int(i)

To get the generator of, say, F, use F.gen(). Or give them different names like so (still using x in R):

G.<a> = GF(2^8, modulus=x^8 + x^5 + x^3 + x + 1)
F.<b> = GF(2^3, modulus=x^3 + x^2 + 1)
2019-05-30 10:19:34 -0500 commented question Sagemath not evaluating complicated expression

I reported this (sub)issue as trac ticket #27897.

2019-05-30 09:56:08 -0500 commented question Sagemath not evaluating complicated expression

There seems to be some weird "holding" going on, e.g. sqrt(abs(1+I)^2 + 14) evaluates to sqrt(16) instead of 4.

2019-05-29 15:56:10 -0500 commented question How to perform modulus in polynomial rings..

f and g are not polynomials but rational functions (with denominators), so what do you mean? What quotient and remainder do you expect?

2019-05-23 10:01:51 -0500 received badge  Nice Answer (source)
2019-05-22 17:05:02 -0500 commented answer Formal determinant of symbolic matrix

Easier would be to implement the determinant separately, for matrices with ordinary entries, but as a formal sum.

2019-05-22 15:58:31 -0500 answered a question Formal determinant of symbolic matrix

You can implement your own "ring" as described in How to implement new algebraic structures in Sage.

Here's a start:

class MyFormalSum(sage.structure.element.Element):
    def __init__(self, parent, x):
            x_iter = iter(x)
            self.x = list(x_iter)
            self.x = [(1,x)]
        # simplify:
        self.x = [(a,b) for (a,b) in self.x if not (a == 0 or b == 0)]
        sage.structure.element.Element.__init__(self, parent)
    def _repr_(self):
        return " + ".join('%s*%s' % (a,b) for (a,b) in self.x)
    def _add_(self, other):
        C = self.__class__
        return C(self.parent(), self.x + other.x)
    def _sub_(self, other):
        C = self.__class__
        return C(self.parent(), self.x + [(-a,b) for (a,b) in other.x])
    def _mul_(self, other):
        import itertools
        C = self.__class__
        return C(self.parent(), [(a*c, b*d) for ((a,b),(c,d)) in itertools.product(self.x, other.x)])

class MyFormalSums(sage.structure.unique_representation.UniqueRepresentation, Ring):
    Element = MyFormalSum
    def __init__(self, base):
        Ring.__init__(self, base)
    def _repr_(self):
        return "MyFormalSums(%s)"%repr(self.base())
    def base_ring(self):
        return self.base().base_ring()
    def characteristic(self):
        return self.base().characteristic()

Now you can do:

sage: R.<x,y> = PolynomialRing(QQ)
sage: S = MyFormalSums(QQ)
sage: Matrix(S, [[x,y], [x,y]]).determinant()
1*x*y + -1*x*y
sage: Matrix(S, [[x,0],[x,y]]).determinant()

Note that this is a hack because the thing is not actually a ring, e.g. TestSuite(S).run() fails.

2019-05-22 03:09:40 -0500 commented question What does cohomology_generators actually do?

I said "suggest" because the code obviously doesn't work, as you demonstrated.

2019-05-21 12:46:39 -0500 commented question What does cohomology_generators actually do?

Strictly speaking the documentation does not claim that the list is anyhow minimal. However the stated algorithm does suggest that this is attempted. Do you see anything wrong with the algorithm or its implementation (approx. 40 lines of code)?

2019-05-21 04:44:32 -0500 received badge  Nice Answer (source)
2019-05-20 06:29:43 -0500 commented question Maximum algebraic connectivity from a given collection of graphs
2019-05-20 04:44:52 -0500 commented answer Maximum algebraic connectivity from a given collection of graphs

It's more efficient to let nauty restrict the number of edges: list(graphs.nauty_geng("8 9:9 -c")).

2019-05-17 05:33:29 -0500 commented question Is there any way to plot3d latex package on sagemath?

Please add the code you have so far to your question.

2019-05-15 15:48:18 -0500 received badge  Nice Answer (source)