In the previous post, I said that Lissajous curves are the result of plotting a curve whose x and y coordinates come from (undamped) harmonic oscillators. If we add a little bit of dampening, multiplying our cosine terms by negative exponentials, the resulting curve is called a harmonograph.

Here’s a bit of Mathematica code to draw harmonographs.

harmonograph[f1_, f2_, ph1_, ph2_, d1_, d2_, T_] := 
    ParametricPlot[{
        Cos[f1 t + ph1] Exp[-d1 t], 
        Cos[f2 t + ph2] Exp[-d2 t]
    }, {t, 0, T}]

For example,

harmonograph[3, 4, 0, Sqrt[2], 0.01, 0.02, 100]

is a slight variation on one of the curves from the previous post, adding dampening terms exp(-0.01t) and exp(-0.01t). The result is the following.

You can get a wide variety of graphs by playing with the parameters. A couple years ago I wrote about an attempt to reproduce the background of the image in the poster for Alfred Hitchcock’s movie Vertigo. The curve is another example of a harmonograph.

Related posts

• Lissajous curves and Lissajous knots
• Recreating the Vertigo poster
• Four part series on mechanical vibrations

Suppose that over time the x and y coordinates of a point are both given by a harmonic oscillator, i.e.

x(t) = cos(n_x t + φ_x)
y(t) = cos(n_y t + φ_y)

Then the resulting path is called a Lissajous curve.

If you add a z coordinate also given by harmonic oscillator

z(t) = cos(n_z t + φ_z)

then the result is a Lissajous knot. This means that if we project a Lissajous knot onto a coordinate plane, we get a Lissajous curve.

For example, here is a Lissajous knot. See the code at the bottom for parameter details.

Here are its projections onto the xy plane,

the yz plane,

and the xz plane.

If you’ve ever seen an oscilloscope, you may have seen Lissajous curves on its screen.

In order to get a knot, we require that the frequencies, the n‘s, are pairwise relatively prime. We also require that the curve not intersect itself, which imposes a constraint on the phases, the φ’s, as well.

Last week I used torus knots as an example of a family of knots. The Lissajous knots are a distinct family of knots: no torus knot is also a Lissajous knot.

Here’s the Mathematica code to produce the images above.

x[t_] := Cos[3 t]
y[t_] := Cos[4 t + Sqrt[2]]
z[t_] := Cos[5 t + Sqrt[3]]

ParametricPlot3D[{x[t], y[t], z[t]}, {t, 0, 2 Pi}]

ParametricPlot[{x[t], y[t]}, {t, 0, 2 Pi}]
ParametricPlot[{y[t], z[t]}, {t, 0, 2 Pi}]
ParametricPlot[{x[t], z[t]}, {t, 0, 2 Pi}]

Take a calculator and enter any number. Then press the cosine key over and over. Eventually the numbers will stop changing. You will either see 0.99984774 or 0.73908513, depending on whether your calculator was in degree mode or radian mode. This is an example of a fixed point, a point that doesn’t change when you apply a function.

The example above is actually two examples, one for cosine of x degrees and one for cosine of x radians. These are two different functions, and they have different fixed points. Note that the two fixed points are not simply related to each other by converting between degrees and radians.

Contraction mapping theorem

The functions

f(x) = cos(x)

and

g(x) = cos(πx/180)

are both contraction mappings, the former corresponding to radians and the latter to degrees. A function h(x) is a contraction mapping if for any two points x and y,

|h(x) – h(y)| ≤ c|x – y|

for some constant c < 1. You can use the mean value theorem to show that c = sin(1) for the function f, and c = π sin(π/180) for the function g.

The contraction mapping theorem says that if a function h is a contraction mapping on a closed interval, then h has a unique fixed point. You can generalize this from working on closed interval to working in any complete metric space.

Weak contraction mappings

The constant c above is important. Without it, you don’t necessarily have a fixed point. For example, let

f(x) = x + 1/x

on the interval [1, ∞). Then you can show that

|f(x) – f(y)| < |x – y|

for x ≠ y but the function f has no fixed point. To see this, suppose f(x) = x, then 1/x = 0, but the reciprocal of a positive number is positive.

If you start with a point x ≥ 1 and repeatedly apply f, you’ll get a sequence of points that move closer to each other, but not closer to any particular number. The sequence of iterates diverges.

A function satisfying the inequality above is called weak contraction mapping. In general, a weak contraction mapping need not have a fixed point, as in the example above. But a weak contraction on a compact space will have a fixed point.

Commentary

As discussed in the link about Kepler’s equation below, something like the contraction mapping theorem was used in practice long before it was formalized in theory.

When I was in grad school, I studied nonlinear PDEs. After seeing a long list of existence theorems and getting lost in the proofs, I wondered what was at the heart of all these theorems. Under all the layers of operator and function space definitions, somewhere there had to be core theorems that said something exists, and I went back over my lecture notes to dig for what these theorems were. In every case it boiled down to a fixed point theorem, not the contraction mapping theorem per se but analogous theorems.

Related posts

• Kepler and the contraction mapping theorem
• Fixed points of the Cliff random number generator

Suppose you have a polynomial p(x) and in interval [a, b] and you want to know how many distinct real roots the polynomial has in the interval. You can answer this question using Sturm’s algorithm.

Let p₀(x) = p(x) and letp₁(x) be its derivative p‘(x).

Then define a series of polynomials for i ≥ 1

p_i+1(x) = – p_i-1(x) mod p_i(x)

until you reach a constant. Here f mod g means the remainder when f is divided by g.

This sequence of polynomials is called the Sturm series. Count the number of sign changes in the Sturm series at a and at b. Then the difference between these two counts is the number of distinct roots of p in the interval [a, b].

Example with Mathematica

As an example, let’s look at the number of real zeros for

p(x) = x⁵ – x – c

for some constant c. We’ll let Mathematica calculate our series.

    p0[x_] := x^5 - x - c
    p1[x_] := D[p0[x], x]
    p2[x_] := -PolynomialRemainder[p0[x], p1[x], x]
    p3[x_] := -PolynomialRemainder[p1[x], p2[x], x]

This works out to

p₀(x) = x⁵ – x – c
p₁(x) = 5x⁴ – 1
p₂(x) = (4/5)x + c
p₃(x) = 1 – 3125c⁴/256

Now suppose c = 3 and we want to find out how many zeros p has in the interval [0, 2].

Evaluating our series at 0 we get -3, -1, 3, -3109/16. So the pattern is – – + -, i.e. two sign changes.

Evaluating our series at 2 we get 27, 79, 23/5, -3109/16. So the pattern is + + + -, i.e. one sign change.

This says x⁵ – x – 3 has one real root between 0 and 2.

By the way, you can multiply the polynomials in the sequence by any positive constant you like if that makes calculations easier. This multiplies subsequent polynomials by the same amount and doesn’t change the signs.

Fine print

Note that the algorithm counts distinct roots; multiple roots are counted once.

You can let the end points of the interval be infinite, and so count all the real roots.

I first tried using Mathematica’s PolynomialMod function, but

   PolynomialMod[5 x^4 - 1, 4 x/5 + c]

gave the unexpected result 5x⁴ – 1. That’s because PolynomialMod does not let you specify what the variable is. It assumed that 4x/5 + c was a polynomial in c. PolynomialRemainder is explicit about the variable, and that’s why the calls to this function above have x as the last argument.

Related post

• Few coefficients, few roots
• Interlaced zeros of ODEs