Polynomials in python
Posted January 22, 2013 at 09:00 AM | categories: polynomials, math | tags:
Updated February 27, 2013 at 02:53 PM
Polynomials can be represented as a list of coefficients. For example, the polynomial \(4*x^3 + 3*x^2 -2*x + 10 = 0\) can be represented as [4, 3, -2, 10]. Here are some ways to create a polynomial object, and evaluate it.
import numpy as np ppar = [4, 3, -2, 10] p = np.poly1d(ppar) print p(3) print np.polyval(ppar, 3) x = 3 print 4*x**3 + 3*x**2 -2*x + 10
139 139 139
numpy makes it easy to get the derivative and integral of a polynomial.
Consider: \(y = 2x^2 - 1\). We know the derivative is \(4x\). Here we compute the derivative and evaluate it at x=4.
import numpy as np p = np.poly1d([2, 0, -1]) p2 = np.polyder(p) print p2 print p2(4)
4 x 16
The integral of the previous polynomial is \(\frac{2}{3} x^3 - x + c\). We assume \(C=0\). Let us compute the integral \(\int_2^4 2x^2 - 1 dx\).
import numpy as np p = np.poly1d([2, 0, -1]) p2 = np.polyint(p) print p2 print p2(4) - p2(2)
3 0.6667 x - 1 x 35.3333333333
One reason to use polynomials is the ease of finding all of the roots using numpy.roots.
import numpy as np print np.roots([2, 0, -1]) # roots are +- sqrt(2) # note that imaginary roots exist, e.g. x^2 + 1 = 0 has two roots, +-i p = np.poly1d([1, 0, 1]) print np.roots(p)
[ 0.70710678 -0.70710678] [ 0.+1.j 0.-1.j]
There are applications of polynomials in thermodynamics. The van der waal equation is a cubic polynomial \(f(V) = V^3 - \frac{p n b + n R T}{p} V^2 + \frac{n^2 a}{p}V - \frac{n^3 a b}{p} = 0\), where \(a\) and \(b\) are constants, \(p\) is the pressure, \(R\) is the gas constant, \(T\) is an absolute temperature and \(n\) is the number of moles. The roots of this equation tell you the volume of the gas at those conditions.
import numpy as np # numerical values of the constants a = 3.49e4 b = 1.45 p = 679.7 # pressure in psi T = 683 # T in Rankine n = 1.136 # lb-moles R = 10.73 # ft^3 * psi /R / lb-mol ppar = [1.0, -(p*n*b+n*R*T)/p, n**2*a/p, -n**3*a*b/p]; print np.roots(ppar)
[ 5.09432376+0.j 4.40066810+1.43502848j 4.40066810-1.43502848j]
Note that only one root is real (and even then, we have to interpet 0.j as not being imaginary. Also, in a cubic polynomial, there can only be two imaginary roots). In this case that means there is only one phase present.
1 Summary
Polynomials in numpy are even better than in Matlab, because you get a polynomial object that acts just like a function. Otherwise, they are functionally equivalent.
Copyright (C) 2013 by John Kitchin. See the License for information about copying.