Tags

, ,

I haven’t put nearly the energy into this blog as I have my main blog, Logos Con Carne. My intentions are good, but somehow I never seem to get around to posting here. (It’s certainly not due to lack of interest.)

In an attempt to get more in the habit, I thought I’d write about some simple fun I had recently with a class for calculating polynomials. It was inspired by a lesson from a set of really fun Python tutorial YouTube videos.

Firstly, a polynomial is an expression that looks like this:

\displaystyle{f}({x})={c_4}\,{x}^{4}+{c_3}\,{x}^{3}+{c_2}\,{x}^{2}+{c_1}\,{x}+{c_0}

Where each coefficient, c, is any real number. The example above is a polynomial of degree 4, the highest power of x. In fact, in a polynomial all the terms are powers of x, so the same thing fully expressed is:

\displaystyle{f}({x})={c_4}\,{x}^{4}+{c_3}\,{x}^{3}+{c_2}\,{x}^{2}+{c_1}\,{x}^{1}+{c_0}\,{x}^{0}

The coefficient numbering probably makes more sense now; it just matches the exponent. The first version abbreviates because x¹=x and xº=1. Note that any of the coefficients can be zero, which effectively removes that term from the sum. A simple expression such as x² is really:

\displaystyle{f}({x})={1.0}\,{x}^{2}+{0}\,{x}^{1}+{0}\,{x}^{0} ={x}^{2}

Coefficients can also be negative. Here’s a real-life example:

\displaystyle{f}({x})={1.2}\,{x}^{4}+{2.5}\,{x}^{3}-{3.14}\,{x}^{2}-{5.25}\,{x}-{1.66}

Note that all that’s required is a list of coefficients: [1.2, 2.5, -3.14, -5.25, -1.66]. The algorithm can build the correct polynomial based on positional placement. That means the x² function above would have the inputs: [1, 0, 0].

§

Now here’s the video that inspired me (it’s only six minutes):

Lambda expressions are one of my favorite parts of Python. They’re inspired by (but not the same as) the lambda calculus that’s a mathematical equivalent to Turing Machines. Both are used interchangeably to define digital computation. In Python, lambda expressions are just mini-functions, but they’re cool because they’re closures.

Lambda expressions that return lambda expressions especially tickle my fancy. See my post about currying functions for an exploration of functions that return functions.

Computer science aside, the regularity of polynomials makes them a good topic for programming students. Often the assignment is to just parse the polynomial text into tokens. The intent is to validate the polynomial as “well-written” (having correct syntax).

The funny robot in the tutorial uses a lambda expression to create an degree-2 (quadratic) polynomial:

# Her version (mostly)...
def build_quadratic_function (a, b, c):
   """Returns f(x) = ax^2 + bx + c"""
   return lambda x: a*x**2 + b*x + c
#
# My version...
quadfunc = lambda a,b,c: lambda x:((a*x*x)+(b*x)+(c))
#

For line length I shortened the docstring and I made a style change I can’t live without – a space comes after the function name when its defined. I equally insist on no space when the function is called.

My one-line lambda-lambda version is shown below. It includes another style change I feel is important: Always Use Parentheses. Her version works perfectly fine, but I never trust operator precedence. (More I don’t trust my memory or that of other programmers. Parentheses keep us all safe.)

The two versions work identically, but they’re restricted to quadratics.

§

It got me to thinking how one could create a function that returned a polynomial of any degree, not just quadratics. It turned out to be way less tricky than I thought it might be.

Here’s the simple version:

def polyn (coeffs):
    '''Generate a polynomial function.'''
    # Function to create one term...
    poly1 = lambda c,e:lambda x:(c*pow(x,e))
    # Create a list of terms...
    terms = [poly1(c,ix) for ix,c in enumerate(reversed(coeffs))]
    # Return a function that applies x to those terms...
    return lambda x:[f(x) for f in reversed(terms)]
#

It takes a list of coefficients such as I described above. It returns a function that takes a value, x, and returns a list of numeric terms. The sum of those terms is the f(x)=y for that x.

For example:

# Create the polynomial x^2...
f = polyn([1,0,0])
# Calculate x=[0,1,2,3,4]...
terms = [sum(f(x)) for x in range(5)]
print(terms)
#

Which calculates the values of f(0) through f(4) in steps of 1 for the function f(x)=x². When run it prints: [0, 1, 4, 9, 16]. Obviously those are the squares of [0, 1, 2, 3, 4].

A use case might involve a function curve calculator. That requires generating a range of real x values from some minimum, xmin, to some maximum, xmax, in steps of some interval, xstep. An alternate interface specifies the number of steps to divide the range into, rather than the amount to step. The former gives you precise steps, the latter gives you a precise number of steps.

Below is an example of the latter:

def divisions (start, stop, steps):
    '''Return list of numbers with N steps.'''
    s0 = float(start)
    s1 = float(stop)
    sx = int(steps)
    assert (s0 < s1)
    assert (0 < sx)
    interval = (s1 - s0) / float(sx)
    return [s0+(interval*n) for n in range(sx+1)]
#

And here’s an example of the former (it’s somewhat similar to Python’s range function but with floats and required step value):

def intervals (start, stop, interval):
    '''Return list of numbers by X interval.'''
    s0 = float(start)
    s1 = float(stop)
    sx = float(interval)
    assert (s0 < s1)
    assert (0 < sx)
    acc = [s0] # start...
    while True:
        x = acc[-1] + sx
        if s1 < x: break # stop!
        acc.append(x)
    return acc
#

Note that in both functions the stop value is included in the result!

These are simple examples. Many libraries provide more robust and capable versions. The well-known NumPy library, for example, has many such functions.

§ §

If you know me (and you almost surely don’t), I’ll wrap any nice handful of code in a class that encapsulates its data and extends its capabilities (because a nice array of methods can share the same private data; I love class design).

Of course I made a polynomial-generating class:

class PolyN (object):
    def __init__ (self, coeffs):
        ''''Create a new polynomial object.'''
        # Generating function for one term...
        poly = lambda c,e:lambda x:(c*pow(x,e))
        # Generate a list of terms...
        _coeffs = enumerate(reversed(coeffs))
        terms = [poly(c,ix) for ix,c in _coeffs]
        # Keep a copy of the coefficients...
        self.coeffs = coeffs
        # Bind the list of terms in a function...
        self.expr = lambda x:[f(x) for f in reversed(terms)]
        # Evaluate f(0) for default values...
        self(0)

    def __call__ (self, x):
        '''Evaluate f(x); return result.'''
        # Generate a list of values...
        self.vals = self.expr(x)
        # Sum the values to get (and return) a result...
        self.result = sum(self.vals)
        return self.result

    def __str__ (self):
        '''Pretty print the polynomial.'''
        cs = enumerate(reversed(self.coeffs))
        ts = ['(%+.2f)x^%d' % (c,ix) for ix,c in cs]
        return ' '.join(reversed(ts))

    def __repr__ (self):
        '''Return the coefficients as a string.'''
        return str(self.coeffs)
#

Of course, always implement the toString method — in this case Python’s __str__ method (and the __repr__ method for good measure). (Oops. I was going to link to my post about that, but I never got around to publishing it. Need to fix that post toasties!)

(One improvement I can imagine is a method that returns a new PolyN object that is the derivative of the polynomial. That shouldn’t be too hard, actually. Hmmm…)

With the class, it’s simple to implement a function charting application. Just generate a set of xs and ys like this:

# Create the polynomial to graph...
ftn = PolyN(coefficients)
# Generate the xs...
xs = interval(xmin, xmax, incr)
# Generate the ys...
ys = [ftn(x) for x in xs]
#

And then use your favorite plotting library to create a graph of those xs and ys:

That’s a graph of the real-life example from above. It’s an degree-4 polynomial with coefficients: [1.2, 2.5, -3.14, -5.25, -1.66]. The “W” shape is characteristic of degree-4 polynomials. (As with all even-degree polynomials, it’s possible to have no real roots, but there are always complex roots.)

Put it all into a handy bit of code, and you have a nice graphing calculator for polynomials of any degree.

§

A brief note about my naming convention for Python lists. If they’re characterized by having a type, even an abstract one, I’m inclined to name the list as a plural of the type’s name. For example, if writing a poker routine, I might have a list of cards or hands, presuming a card or hand data type.

Here xs is a list of x values, ys is a list of y values. That might seem brief, but it’s a style I inherited from Haskell, and I rather like it. I feel brief variable names are fine if their scope is restricted to only a few lines and the abbreviation is obvious.

Ø