Thank you, that was very instructive to see the "right way" to start by using an appropriate ring. I guess you can go on to divide out by the other linear factor and get the quadratic equation, and solve it, but I got that far myself, and then could not work with the solutions of the quadratic equation. But you may be able to do that because then they will belong to some algebraic number field which will come equipped with useful algorithms.
On Thursday, November 15, 2018 at 9:47:24 AM UTC-8, slelievre wrote: > > For the problem at hand here, I would work with polynomials in `x` > with coefficients in the ring of polynomials in `N` and `M` over > the field of algebraic numbers, and do something like the following. > > Define the ring of polynomials over the field of algebraic numbers. > > sage: R.<N, M> = QQbar[] > > Define the ring of polynomials over the above polynomial ring. > > sage: S.<x> = R[] > > Define `a`, `i`, `b`, `c`, `X`, `f` from the question. > > sage: a = QQbar(3).sqrt() / 2 > sage: i = QQbar(I) > sage: b = (x - ~x) / (2 * i) > sage: c = (a * (x + ~x) + b) / 2 > sage: X = (M / 3) * (a + b + c) > sage: f = 24 * (X^2 - N * b * c) * x^2 > > Note that the definition of `b` involves the inverse of `x` (which can be > denoted by `~x` or `x^-1`) and therefore lives in the fraction field of > `S`, > rather than in `S` itself. Therefore, so do `c`, `X`, and `f`. > > Check if `f` however represents an element in `S` as follows: > > sage: f in S > True > > and then construct the corresponding element of `S` (we could call it `f` > but here we call it `ff`). > > sage: ff = S(f) > > Check that `ff` vanishes at `-1`, or equivalently is divisible by `x + 1`. > > sage: ff(-1) > 0 > sage: d = x + 1 > sage: d.divides(ff) > True > > Perform exact division (to stay in the polynomial ring `S` rather than move > to its fraction field). > > sage: g = ff // d > > Define `t` as an algebraic number. > > sage: t = QQbar(exp(-pi*I/3)) > > Check that `g` vanishes at `t` or equivalently is divisible by `x - t`. > > sage: g(t) > 0 > sage: dd = x - t > sage: dd.divides(g) > True > > Perform exact division. > > sage: h = g // dd > > Inspect the result: > > sage: h > ((-1.000000000000000? - 1.732050807568878?*I)*M^2 > + (3.000000000000000? + 5.196152422706632?*I)*N)*x^2 > + ((1.000000000000000? - 1.732050807568878?*I)*M^2 > + (3.000000000000000? - 5.196152422706632?*I)*N)*x > + 2*M^2 + (-6)*N > > Or as a list (the list of coefficients of `x^j`, for `j` from `0` to the > degree): > > sage: h.list() > [2*M^2 + (-6)*N, > (1.000000000000000? - 1.732050807568878?*I)*M^2 > + (3.000000000000000? - 5.196152422706632?*I)*N, > (-1.000000000000000? - 1.732050807568878?*I)*M^2 > + (3.000000000000000? + 5.196152422706632?*I)*N] > > Define somme pretty_print functions to inspect these algebraic numbers: > > def pretty_print_qqbar(z): > a, b = z.real().radical_expression(), z.imag().radical_expression() > return "{} + {}*i".format(a, b) > > def pretty_print_x_coeff(xc): > return " + ".join("({})*{}".format(pretty_print_qqbar(xc[m]), m) > for m in xc.monomials()) > > and print `h` in a nice form: > > print("h = " + > "\n + ".join("({})*x^{}".format(pretty_print_x_coeff(xc), k) > for k, xc in enumerate(g.list()))) > > The result is: > > h = ((-1 + sqrt(3)*i)*M^2 + (3 + -3*sqrt(3)*i)*N)*x^0 > + ((3 + sqrt(3)*i)*M^2 + (-3 + 3*sqrt(3)*i)*N)*x^1 > + ((3 + -sqrt(3)*i)*M^2 + (-3 + -3*sqrt(3)*i)*N)*x^2 > + ((-1 + -sqrt(3)*i)*M^2 + (3 + 3*sqrt(3)*i)*N)*x^3 > > -- You received this message because you are subscribed to the Google Groups "sage-support" group. To unsubscribe from this group and stop receiving emails from it, send an email to sage-support+unsubscr...@googlegroups.com. To post to this group, send email to sage-support@googlegroups.com. Visit this group at https://groups.google.com/group/sage-support. For more options, visit https://groups.google.com/d/optout.