Working with polynomials and calculus

Polynomials are among the simplest functions in mathematics and are defined as a sum:

x represents a placeholder to be substituted, and a_i is a number. Since polynomials are simple, they provide an excellent means for a brief introduction to calculus. Calculus concerns the differentiation and integration of functions. Integration is, roughly speaking, anti-differentiation, in the sense that first integrating and then differentiating yields the original function.

In this recipe, we will define a simple class that represents a polynomial and write methods for this class to perform differentiation and integration.

Getting ready

Geometrically, the derivative, obtained by differentiating, of a function is its gradient, and the integral, obtained by integrating, of a function is the area that lies between the curve of the function and the x axis, accounting for whether the curve lies above or below the axis. In practice, differentiating and integrating are done symbolically, using a set of rules and standard results that are particularly simple for polynomials.

There are no additional packages required for this recipe.

How to do it...

The following steps describe how to create a class representing a polynomial and implement differentiation and integration methods for this class:

1. Let's start by defining a simple class to represent a polynomial:

class Polynomial:
    """Basic polynomial class"""

    def __init__(self, coeffs):
        self.coeffs = coeffs

    def __repr__(self):
        return f"Polynomial({repr(self.coeffs)})"

    def __call__(self, x):
        return sum(coeff*x**i for i, coeff 
                in enumerate(self.coeffs))

1. Now that we have defined a basic class for a polynomial, we can move on to implement the differentiation and integration operations for this Polynomial class to illustrate how these operations change polynomials. We start with differentiation. We generate new coefficients by multiplying each element in the current list of coefficients without the first element. We use this new list of coefficients to create a new Polynomial instance that is returned:

    def differentiate(self):
        """Differentiate the polynomial and return the derivative"""
        coeffs = [i*c for i, c in enumerate(self.coeffs[1:], start=1)]
        return Polynomial(coeffs)

1. To implement the integration method, we need to create a new list of coefficients containing the new constant (converted to a float for consistency) given by the argument. We then add to this list of coefficients the old coefficients divided by their new position in the list:

    def integrate(self, constant=0):
        """Integrate the polynomial, returning the integral"""
        coeffs = [float(constant)] 
        coeffs += [c/i for i, c in enumerate(self.coeffs, start=1)]
        return Polynomial(coeffs)

1. Finally, to make sure these methods work as expected, we should test these two methods with a simple case. We can check this using a very simple polynomial, such as x² - 2x + 1:

p = Polynomial([1, -2, 1])
p.differentiate()
# Polynomial([-2, 2])
p.integrate(constant=1)
# Polynomial([1.0, 1.0, -1.0, 0.3333333333])

How it works...

Polynomials offer an easy introduction to the basic operations of calculus, but it isn't so easy to construct Python classes for other general classes of functions. That being said, polynomials are extremely useful because they are well understood and, perhaps more importantly, calculus for polynomials is very easy. For powers of a variable x, the rule for differentiation is to multiply by the power and reduce the power by 1, so that xⁿ becomes nx^n-1.

Integration is more complex, since the integral of a function is not unique. We can add any constant to an integral and obtain a second integral. For powers of a variable x, the rule for integration is to increase the power by 1 and divide by the new power, so that xⁿ becomes xⁿ⁺¹/(n+1), so to integrate a polynomial, we increase each power of x by 1 and divide the corresponding coefficient by the new power.

The Polynomial class that we defined in the recipe is rather simplistic, but represents the core idea. A polynomial is uniquely determined by its coefficients, which we can store as a list of numerical values. Differentiation and integration are operations that we can perform on this list of coefficients. We include a simple __repr__ method to help with the display of Polynomial objects, and a __call__ method to facilitate evaluation at specific numerical values. This is mostly to demonstrate the way that a polynomial is evaluated.

Polynomials are useful for solving certain problems that involve evaluating a computationally expensive function. For such problems, we can sometimes use some kind of polynomial interpolation, where we "fit" a polynomial to another function, and then use the properties of polynomials to help solve the original problem. Evaluating a polynomial is much "cheaper" than the original function, so this can lead to dramatic improvements in speed. This usually comes at the cost of some accuracy. For example, Simpson's rule for approximating the area under a curve approximates the curve by quadratic polynomials over intervals defined by three consecutive mesh points. The area below each quadratic polynomial can be calculated easily by integration.

There's more...

Polynomials have many more important roles in computational programming than simply demonstrating the effect of differentiation and integration. For this reason, a much richer Polynomial class is provided in the NumPy package, numpy.polynomial. The NumPy Polynomial class, and the various derived subclasses, are useful in all kinds of numerical problems, and support arithmetic operations as well as other methods. In particular, there are methods for fitting polynomials to collections of data.

NumPy also provides classes, derived from Polynomial, that represent various special kinds of polynomials. For example, the Legendre class represents a specific system of polynomials called the Legendre polynomials. The Legendre polynomials are defined for x satisfying -1 ≤ x ≤ 1 and form an orthogonal system, which is important for applications such as numerical integration and the finite element method for solving partial differential equations. The Legendre polynomials are defined using a recursive relation. We define

and for each n ≥ 2, we define the nth Legendre polynomial to satisfy the recurrence relation,

There are several other so called orthogonal (systems of) polynomials, including Laguerrepolynomials, Chebyshev polynomials, and Hermite polynomials.

See also

Calculus is certainly well documented in mathematical texts, and there are many textbooks that cover from the basic methods all the way to the deep theory. Orthogonal systems of polynomials are also well documented among numerical analysis texts.