Lecture 11: Univariate and Multivariate Polynomials
===================================================

There may not be much difference at the low level of expressions
where we store parameters as ordinary variables, at the mathematical level
polynomials in several variables are different from polynomials in one
variable that have symbols as coefficients.
In this lecture we make the important connection between factoring
and root finding, symbolically as well as numerically.

Polynomials as Expressions
--------------------------

Let us start with a definition in plain words.
A univariate polynomial is a finite sum of terms,
where every term is a coefficient multiplied with a monomial,
and where every monomial is a power of the variable,
by default, call this variable ``x``.
All coefficients in the polynomial are of the same type.

::

   p = x^4 - 4*x^2 - 7*x + 9
   print(p, 'has type', type(p))
   print('the degree :', p.degree(x))

The degree of an expression returns the highest power of the variable
given as an argument of the degree method applied to the expression.
Can we also select the coefficients of the expression?

::

   print('the coefficient of x^2 :', p.coefficient(x,2))
   print('all coefficients :', [p.coefficient(x,k) for k in range(p.degree(x)+1)])
   print('the leading coefficient :', p.leading_coefficient(x))

If we select all coefficients of ``p``, then we get a list of lists.
Each list gives the coefficient and the corresponding power.

::

   print('all coefficients :', p.coefficients())
   print('all coefficients in x :', p.coefficients(x))

If we want to see the terms of the polynomial symbolically,
then we can ask for the operands.

::

   print('the operands of p :', p.operands())
   print('the topmost operator of p :', p.operator())

We see that we may view the expression p as a sum.
If we ask for the roots, then see what happens!

::

   print(p.roots())

We see the roots expressed with the ``sqrt()`` function.
If we are not happy with this output, then we have to change the ring.

::

   print('the real roots :', p.roots(ring=RR))
   print('the complex roots :', p.roots(ring=CC))

Our polynomial has two real roots and two complex conjugated roots.

Univariate Polynomials
----------------------

We can explicitly declare a polynomial ring and coerce our expression
for ``p`` into this ring.

::

   P.<x> = PolynomialRing(QQ)
   q = P(p)
   print(q, 'has type', type(q))

Another, shorter, and perhaps more natural way of writing this conversion is

::

   R = QQ[x]
   r = R(p)
   print(r, 'has type', type(r))

If we were no longer interested in our polynomial ``q``,
then we could just as well pick any random quartic polynomial.

::

   print(QQ[x].random_element(degree=4))

But let us continue with our ``q``.
If we compute the degree and coefficients,
then we may not give the argument ``x`` anymore.

::

   print(q.degree())
   print('coefficients of nonzero terms :', q.coefficients())
   print('all coefficients :', q.coeffs())
   print('the list of all coefficients :', q.list())
   print('a dictionary representation :', q.dict())

The keys of the dictionary are the exponents of those monomials
that appear with nonzero coefficient.  The dictionary representation
could be very convenient for *sparse* polynomials, when only relatively
few monomials appear with nonzero coefficient, relative to the degree
of the polynomial.

We can ask if the polynomial is irreducible.

::

   print(q.is_irreducible())
   print(factor(q))

and build a field extension to add a root of ``q`` to ``QQ``

::

   K.<a> = QQ.extension(q)
   print(K)

Now we can view the original expression as a polynomial over the new field,
extended with ``a`` root of ``q``.

::

   q2 = K[x](q)
   print(q2, 'has type', type(q2))

See if ``q2`` factors or not.

::

   print(q2.is_irreducible())
   print(factor(q2))

We see that adding one algebraic number is not enough to factor 
completely in linear factors.  Observe that we must explicitly
coerce the nonlinear factor into a polynomial with coefficients
in ``K``.

::

   f2 = q2/(x-a); print(f2, type(f2), type(K[x](f2)))
   L.<b> = K.extension(K[x](f2))
   print(L)

We take the original polynomial in the new extended coefficient 
ring and factor.

::

   q3 = L[x](q)
   print(factor(q3))

To select the factor to continue with our field extensions,
we convert the factorization to a list.

::

   lf = list(factor(q3))
   print(lf)

Note that the elements in the list are tuples, 
so we must remove the multiplicities.

::

   f3t = lf[2]; print(f3t, 'has type', type(f3t))
   f3 = L[x](f3t[0]); print(f3, 'has type', type(f3))

Now make our last field extension!

::

   M.<c> = L.extension(f3)
   print(M)

As the polynomial now factors completely, we have solved the polynomial
symbolically, expressing the roots as algebraic numbers a, b, and c.

::

   q4 = M[x](q)
   print(factor(q4))

The polynomial ``q`` is now factored completely as a product of linear
factors, as we see
``(x - c) * (x - b) * (x - a) * (x + c + b + a)``.
Observe the connection between factoring and root finding.

::

   print(q4.roots())

and we see ``[(c, 1), (b, 1), (a, 1), (-c - b - a, 1)]``.
We can still find the roots numerically of our quartic.

::

   crts = q4.roots(ring=CC)
   print(crts)

Because the coefficient in the term with \ :math:`x^3` is zero,
the sum of the roots in both the symbolic and the numeric representation
must be zero.  We can see this easily from the symbolic representation,
but let us verify this on the numerical representation.

::

   print(sum(n for (n,1) in crts))

To recapitulate, we distinguish between two main forms of root finding,
one is symbolic, the other numeric.
The default numeric field is the field of complex numbers,
whereas symbolically we extend the field of rational numbers
with sufficiently many symbols to represent all roots.

Multivariate Polynomials
------------------------

Polynomials in several variables are declared similarly as polynomials 
in one variable.
The quotes around the names are needed if we have not used
or declared them explicitly before as variables.  
We take a random polynomial of degree 4, and with terms we can
give an upper bound on the number of terms in the polynomial.

::

   R.<x,y> = PolynomialRing(QQ)
   p = R.random_element(degree=4,terms=10)
   print(p, 'has type', type(p))

We can select the monomials, get its variables, and its degrees.

::

   print('the monomials :', p.monomials())
   print('corresponding coefficients :', p.coefficients())
   print('the variables :', p.variables())
   print('the degree in x :', p.degree(x))
   print('the degree in y :', p.degree(y))
   print('the degree of p :', p.degree())

The order of the monomials is important.

::

   print(R.term_order())

The default order appears to be ``Degree reverse lexicographic term order``.
To change the ordering of the monomials in the polynomial,
we coerce p into a another ring.
In a lexicographic order, all monomials in which ``x`` occurs come first.

::

   Rlex.<x,y> = PolynomialRing(QQ,order = 'lex')
   print(Rlex); print(Rlex.term_order())
   print(Rlex(p))

We can view a polynomial in several variables as a polynomial in one variable
by collecting terms.  Because of the type of argument of ``polynomial()``,
the selection of the tuple of the outcome of ``p.variables()`` is needed.

::

   print('as polynomial in x :', p.polynomial(p.variables()[0]))
   print('as polynomial in y :', p.polynomial(p.variables()[1]))

Assignments
-----------

1. Execute the following sequence of commands:
   ``var('x'); p = prod([x-k for k in range(20)]); print p; q = p.expand()``
   and then type
   ``print q.roots()`` and ``print q.roots(ring=CDF).``
   Compare the differences between the output of the two ``roots``.
   Did you expect to see those differences?

2. Declare ``x`` as a polynomial variable over the rational number
   with the statement ``x = polygen(QQ)``.

   1. Give the Sage command(s) to
      compute the greatest common divisor of the polynomials
      \ :math:`p = 2 x^5 + 11 x^4 + 14 x^3 + 11 x^2 + 12 x` and
      \ :math:`q = 2 x^5 + 5 x^4 + 7 x^3 + 8 x^2 + 5 x + 3`.

   2. How can Sage compute the cofactors \ :math:`k`
      and \ :math:`\ell` so
      that \ :math:`{\rm gcd}(p,q) = k p + \ell q`?

   3. Finally, give the Sage commands to verify the relation
      \ :math:`{\rm gcd}(p,q) = k p + \ell q` for 
      the \ :math:`k` and the \ :math:`\ell` that were found.

   Hint: do ``help(xgcd)``.

3. Do ``x = polygen(RR); print factor(x^2 - 2.25)``
   and ``x = polygen(QQ); print factor(x^2 - 9/4)``.
   Explain the differences in the outcomes of the two factor commands.

4. Consider the polynomial
   \ :math:`p = 2 x^5 + 9 x^4 + 16 x^3 + 15 x^2 + 12 x + 9`.
   Write \ :math:`p` as a product of linear factors:

   1. symbolically, by adding sufficiently many formal roots; and

   2. numerically, by finding all complex roots of \ :math:`p`. 

5. Consider the polynomial
   \ :math:`p = 2 x z^4 + x z^3 + 2 y z^2`.
   Give the Sage commands to bring \ :math:`p` in the forms

   (a) \ :math:`2 z^4 x + z^3 x + 2 z^2 y`,

   (b) \ :math:`2 x z^4 + x z^3 + 2 z^2 y`, and

   (c) \ :math:`2 y z^2 + 2 x z^4 + x z^3`.

6. Consider the polynomial :math:`p = x^3 + 4 x + 7` 
   over a finite field of 17 elements.

   (a) Give the Sage commands to show that :math:`p`
       is irreducible over this finite field.

   (b) Add sufficiently many formal roots to this finite field so 
       that :math:`p` factors as a product of linear polynomials.  

   (c) Give all relevant Sage commands.
       Write the final factorization of :math:`p`.

7. Consider the statements ``var('x')`` and ``QQ['x']``.

   What is the main difference between the roles of ``'x'``
   after these statements?
   Start your answer with a precise description on the effect 
   of each statement.  Illustrate the difference.

   When should you use ``var('x')``, when ``QQ['x']``?

   What can you do after ``QQ['x']`` but not after ``var('x')``?