\chapter{The Elastica}
\label{elastica-chapter}
The word \emph{spline} originally refers to a thin strip of wood used
to smoothly interpolate control points. A full understanding of
splines requires an examination of the theory of such elastic
strips. This theory, that of the \emph{elastica,} has a long and rich
history, much of which will be explored in the next chapter, but the
literature is missing a concise, accessible treatment of the
topic. This chapter fills that gap.
One difficulty of the literature is that there are many different ways
to approach the elastica: as an equilibrium of moments in a static
physical system; as an equilibrium of forces in the same physical
system (analyzed in terms of a finite element approximation
composed of rigid struts and pivot springs); as the curve that
minimizes the total bending energy of that system; as the solution to
an abstract differential equation; and through the analogue with
another physical system, the pendulum, which happens to share
the same differential equation.
Each approach sheds light on a different facet of the elastica. For
instance, treating it as a statics problem is the easiest way to
understand the effect of different constraints operating on the curve,
and of course the physics interpretation is the only way to think
about the differences between a real elastic lamina and the
mathematical ideal -- especially, the effect of friction at the
constraint points. The variational formulation is most revealing when
examining the minimization of some other metric than total bending
energy. And the pendulum analogy is helpful for developing an
intuitive understanding of the underlying differential equation,
particularly the fact that it is periodic.
My approach is to pull together \emph{all} of these
approaches, letting each tell its own story. For a given problem
involving the elastica, any one of these approaches may
offer the simplest, most intuitive solution. Yet, the underlying
mathematics is the same. In the equilibrium of moments formulation,
the equation takes on a remarkably simple form: $\kappa = \lambda y$,
where $\kappa$ is the curvature, $y$ is the Cartesian ordinate, with
$x$ the abscissa, and $\lambda$ is a constant.
\section{Statement of the elastica problem}
The elastica models a thin strip of flexible material (traditionally
wood, but spring metal will do). In the spirit of idealizing the
problem, ``thin'' means that the thickness is negligible compared
with the amount of bending. The mathematical idealization further
assumes that the centerline of this strip lies entirely on a plane,
and that there is no twisting, stretching, or compression.
Thus, the curve traced by the strip is characterized by its curvature
at each point along its length. A convenient definition of curvature
is the derivative of tangent angle $\theta$ with respect to arc length
$s$. Unless otherwise stated, in this chapter the $\theta'$ notation will
represent $d\theta/ds$, so curvature $\kappa$ is written simply $\theta'$.
Bending such a strip induces potential energy, much like that of a
spring. According to the theory of elasticity, the bending energy is
proportional to the square of the curvature. Thus, the total bending
energy of a strip of length $l$ is
\begin{equation}
\label{elastica-def}
E[\kappa(s)] = \int_0^l \kappa(s)^2\ ds\:.
\end{equation}
When completely unconstrained, the elastica will assume the shape of a
straight line, in which the curvature is everywhere zero, and thus the
total bending energy is also zero. When constrained, the bending
energy will tend to the minimum possible under the constraints.
\section{Variational study}
The elastica equation \ref{elastica-def} is particularly well suited
to study by the calculus of variations \cite{Elsgolts}, a technique for
transforming problems specified in terms of infinite-dimensional
functions into tractable ordinary differential equations. The calculus
of variations was founded at the end of the 17th century, and was
developed by Johann and Jacob Bernoulli, Leonhard Euler, and others.
In particular, Equation \ref{elastica-def} yields readily to the
Euler--Lagrange equation.
\subsection{The Euler--Lagrange equation}
If there exists a minimum value of the functional
\begin{equation}
\label{var-statement}
v[y(x)] = \int_{x_0}^{x_1} F(x, y(x), y'(x))\:dx\: ,
\end{equation}
\noindent subject to the $n$ constraints, each of the form $1 \leq i \leq n$,
\begin{equation}
\label{var-constraints}
\int_{x_0}^{x_1} c_i(y(x))\:dx = k_i\: ,
\end{equation}
\noindent then this minimum is characterized by a tuple of values $\lambda_i$ such that
\begin{equation}
\label{var-solution}
\frac{\partial F}{\partial y} - \frac{d}{dx} \frac{\partial F}{\partial y'} + \sum_i \lambda_i \frac{d c_i}{dy} = 0\:.
\end{equation}
Therefore, the problem of finding a {\em function} minimizing the cost
functional $v$ in an infinite-dimensional space of functions is
reduced to that of finding a tuple of values minimizing that
functional. In general, such an approach yields easily to numerical
techniques for both integrating the differential equation
\ref{var-solution} and finding a set of values $\lambda_i$ minimizing
the cost functional.
\subsection{Application of the Euler--Lagrange equation to the elastica}
There are several ways to apply the Euler--Lagrange equation to find
curves satisfying the elastica problem as defined by Equation
\ref{elastica-def}. The main choice is that of independent and
dependent variables ($x$ and $y$ in the formulation of Equation
\ref{var-statement}). At least three such pairs appear in the
literature of the elastica: arc length and tangent angle, arc length and
signed curvature, and Cartesian coordinates.
Of these, the choice of arc length $s$ as the independent variable and
tangent angle $\theta$ as the dependent yields the simplest
variational solution. The most general endpoint constraints are
specification of curve length, endpoint positions, and endpoint
tangents. Of these, the endpoint position constraint requires Lagrange
multipliers in the form of Equation \ref{var-constraints}. Assume
without loss of generality that the starting point is $(0, 0)$ and the
endpoint is $(x, y)$ (simple translation will allow arbitrary starting
point). Then, the endpoint constraint is expressed as the integral of
the unit speed vector of direction $\theta(s)$:
\begin{equation}
\begin{array}{l}
\label{elastica-pos}
\int_{0}^{l}\cos \theta(s)\: ds = x\:, \\
\int_{0}^{l}\sin \theta(s)\: ds = y\:.
\end{array}
\end{equation}
In terms of arc length and tangent, these constraints are easily
expressible in the form of equation \ref{var-constraints}.
It's now easy to see why representing the curve in terms of tangent
angle as a function of arc length yields most readily to the
variational approach. In terms of curvature $\kappa$, the endpoint
constraint is a double integral, which does not fit neatly into a
formulation of constraints as a single integral of some
function. Alternatively, Cartesian coordinates would make the endpoint
constraints even simpler, but the variational
function becomes more complex, requiring second derivatives as well as
the first derivatives of the Euler--Lagrange equation. The first
solution was indeed done using Cartesian coordinates, but this
approach is definitely simpler.
Applying the variational technique to minimizing the bending energy of
a strip, Equation \ref{elastica-def}, we arrive at the solution,
\begin{equation}
\label{elastica-eq}
\theta'' + \lambda_1 \sin \theta + \lambda_2 \cos \theta = 0\: .
\end{equation}
Equation \ref{elastica-eq} is the simplest form of the solution, but
other forms can be readily derived, and often appear in the
literature. Taking the derivative with respect to $s$,
\begin{equation}
\label{elastica-deriv}
\theta''' + \lambda_1 \theta' \cos \theta
- \lambda_2 \theta' \sin \theta = 0\: .
\end{equation}
Conversely, multiplying both sides of \ref{elastica-eq} by $d \theta / ds$ and then taking
the integral,
\begin{equation}
\label{elastica-int}
\frac{1}{2}(\theta')^2 - \lambda_1 \cos \theta
+ \lambda_2 \sin \theta + C = 0\: .
\end{equation}
Combining \ref{elastica-deriv} and \ref{elastica-int}, we can
eliminate the Lagrange multipliers, retaining the single constant of
integration $C$.
\begin{equation}
\theta''' + {\frac{1}{2}}(\theta')^3 + C\theta' = 0\: .
\end{equation}
Since $\kappa = \theta'$, we can reformulate this differential
equation entirely in terms of curvature:
\begin{equation}
\label{elastica-kappa}
\kappa'' + \frac{1}{2} \kappa^3 + C\kappa = 0\: .
\end{equation}
Another formulation in the literature (for example, Eq. 12 of Lee and
Forsythe \cite{ForsytheLee}) multiplies this last equation by
$1/\kappa$, yielding
\begin{equation}
\label{fl-elastica-eq}
\frac{1}{\kappa}\kappa'' + \frac{1}{2} \kappa^2 + C = 0\: .
\end{equation}
%In this formulation, the units correspond to forces acting on the thin
%beam (with the modulus of elasticity normalized away); $(1 /
%\kappa)(d^2 kappa / ds^2)$ is the longitudinal force on the beam.
%Present equations of equilibrium, and show that variational solution
%is same as equilibrium.
\section{A family of solutions}
\label{elastica-fam-sect}
The equations derived in the previous section for the shape of the
elastica don't have a single solution. Rather, they describe a family
of solutions, characterized by a single scalar parameter. All such
solutions are realizable by a physical elastica; none are purely an
artifact of the mathematical approach.
For mathematical convenience, assign $\kappa = 1$ at the midpoint of
the curve, and use $\lambda_1$ as the parameter in Equation
\ref{elastica-eq}, setting $\lambda_2 = 0$. If we write $\lambda_1 =
\lambda \cos \phi$ and $\lambda_2 = \lambda \sin \phi$, then all other
instances of the elastica can be transformed to meet these conditions
by scaling and rotating the entire system so that $\phi = 0$. Since
the assumption of $\phi = 0$ obviates the need to distinguish multiple
$\lambda$ parameters, it is convenient to drop the subscript ${}_1$
and simply write $\lambda$, resulting in a simplified equation.
\begin{equation}
\label{elastica-eq-simpl}
\theta'' + \lambda \sin \theta = 0\: .
\end{equation}
Figure \ref{elastica-fam} shows a representative selection from
this family, for various values of $\lambda$.
\begin{figure*}
\begin{center}
\includegraphics[scale=.8]{figs/elastica-fam}
\caption{\label{elastica-fam}The family of elastica solutions.}
\end{center}
\end{figure*}
The remainder of this chapter is devoted to exploring these solutions
and their physical meaning, but a few highlights are immediately
apparent. All solutions have a mirror (even) symmetry at each
curvature maximum. All solutions save $\lambda = 0.25$ are
periodic. The solutions for values larger than $0.25$ have
inflections, and have odd symmetry around each inflection point. For
values smaller than $0.25$, there is no inflection; the curvature
always remains the same sign. Note that the curve has a mirror
symmetry for all values of $\lambda$.
As $\lambda$ approaches $0$, the elastica tends to a circle. As
$\lambda$ approaches $\infty$, the elastica tends to a sine wave. The
solution at $\lambda = 0.5$ is special, and is known as the
\emph{rectangular elastica}. The tangents at the inflection points of
this curve are perpendicular to the horizontal axis; the curve is tangent
to the rectangle formed by inflections and curvature extrema.
\begin{figure*}
\begin{center}
\includegraphics[scale=.8]{figs/elastica-fam-k}
\caption{\label{elastica-fam-k}Curvature plots for the family of elastica solutions.}
\end{center}
\end{figure*}
Figure \ref{elastica-fam-k} shows the curvature plots of the
instances of the elastica family shown in Figure \ref{elastica-fam},
plotted with curvature on the vertical axis and arc length on the
horizontal. As is usual, inflection points on the curve correspond to
the curvature crossing (or touching) the horizontal axis.
\section{Equilibrium of moments}
\begin{figure*}
\begin{center}
\includegraphics[scale=.9]{figs/moment}
\caption{\label{moment}Equilibrium of moments.}
\end{center}
\end{figure*}
Analyzing the elastica using the theory of elasticity yields three
physical quantities which must be in equilibrium along its length:
tension (longitudinal force), shearing force, and bending moment. The
tension and shearing forces are useful when working out the effect of
constraints, but it is possible to derive the equation for the shape
of the elastica from moments alone. Indeed, computing the moments is
by far the easiest derivation of the elastica equation.
This section uses elementary principles of statics to carry
out the derivation of the elastica from equilibrium of moments. A
classic textbook introducing these concepts is \emph{Statics,} by
Goodman and Warner \cite{Statics}.
Consider the elastica in Figure \ref{moment}. Two forces of magnitude
$\mathbf{F}$ and opposite directions push the ends of the elastica
together. Because these forces are equal in magnitude, opposite in
direction, and have parallel lines of action, they form a
\emph{couple} \cite[p. 65]{Statics} (the line of action is defined as
a vector in the same direction as the force, emerging from the point
where the force is applied).
At any point along the elastica, then, the \emph{resultant moment} is
simply the force of the couple times the distance from the line of
action: $M = Fy$, where $F$ is the magnitude of the force
vector $\mathbf{F}$.
The simple theory of elasticity says that the relation between bending
moment and curvature is linear. Thus, collecting all the constants
together into a single $\lambda$ yields this remarkably simple formulation
of the shape of a general elastica,
\begin{equation}
\label{eq-simple-elastica}
\kappa = \lambda y\: .
\end{equation}
This equation is readily shown to be equivalent to the variational
solution above. Take the derivative with respect to arc length $s$, and
note that $y' = \sin \theta$. Then,
\begin{equation}
\kappa' = \lambda \sin \theta,
\end{equation}
\noindent and this equation is equivalent to Equation
\ref{elastica-eq-simpl}, assuming without loss of generality that the line
of action is parallel to the horizontal axis.
All elastica solutions are equivalent to the system shown in Figure
\ref{moment}, but in non-inflectional cases the line of action of the
couple is virtual, in that it doesn't touch the curve.
\section{Pendulum analogy}
The differential equation \ref{elastica-eq-simpl} for the shape of the
elastica is mathematically equivalent to that of the dynamics of a
simple swinging pendulum. This kinetic analog is useful for developing
intuition about the solutions of the elastica equation. In particular,
it's easy to see that all solutions are periodic, and that the family
of solutions is characterized by a single parameter (modulo scaling,
translation, and rotation of the system).
\begin{figure*}
\begin{center}
\includegraphics[scale=.9]{figs/pendulum}
\caption{\label{pendulum}An oscillating pendulum.}
\end{center}
\end{figure*}
Consider, as shown in Figure \ref{pendulum}, a weight of mass $m$ at
the end of a pendulum of length $r$, with the angle from vertical
specified as a function of time $\theta(t)$.
Then the velocity of the mass is $r\theta'$ (where, in this section,
the $'$ notation represents derivative with respect to time), and its
acceleration is $r\theta''$. The net force of gravity, acting with the
constraint of the pendulum's angle is $mg \sin \theta$, and thus we
have
\begin{equation}
F_{\mbox{\scriptsize net}} = ma = mr\theta'' = mg \sin \theta\:.
\end{equation}
With the substitution $\lambda_1 = -g/r$, and the replacement of
arc length $s$ for time $t$, this equation becomes equivalent to
Equation \ref{elastica-eq-simpl}, the equation of the elastica. Note that
angle (as a function of arc length) in the elastica corresponds to
angle (as a function of time) in the pendulum, and the elastica's
curvature corresponds to the pendulum's angular momentum.
% According to Love\cite{Love}, the systematic application of the
%kinetic analogue was worked out by W. Hess in 1885.
The swinging of a pendulum is perhaps the best-known example of a
periodic system. Further, it is intuitively easy to grasp the
parameter space of the system. Transformations of scaling (varying the
pendulum length) and translation (assigning different phases of the
pendulum swing to $t = 0$) leave the solutions essentially unchanged.
Only one parameter, how high the pendulum swings, changes the
fundamental nature of the solution.
The solutions of the motion of the pendulum, as do the solutions of
the elastica equation \ref{elastica-eq-simpl}, form a family characterized
by a single scalar parameter, once the trivial transforms of rotation
and scaling are factored out. Without loss of generality, let $t = 0$
at the bottom of the swing (i.e. maximum velocity) and let the
pendulum have unit length. The remaining parameter is, then, the ratio
of the total energy of the system (which remains unchanged over time)
to the potential energy of the pendulum at the top of the swing (the
maximum possible), in both cases counting the potential energy at the
bottom of the swing as zero. Thus, the total system energy is also
equal to the kinetic energy at the bottom of the swing.
For mathematical convenience, define the parameter $\lambda$ as one
fourth the top-of-swing potential energy divided by the total
energy. This convention is equivalent to varying the force of gravity
while keeping the maximum kinetic energy fixed, which may be justified
by noting that the zero-gravity case (which would require infinite
kinetic energy if we were assuming unit gravity) is far more relevant
in the context of splines (it corresponds, after all, to a circular
arc, so any elastica-based spline meeting the criterion of roundness
must certainly have this solution) than the opposite of zero kinetic
energy.
The potential energy at the top of the pendulum is $2mgr$. The kinetic
energy at the bottom of the swing is $\frac{1}{2} mv^2$, where, as
above, $v = r\theta'$. Thus, according to the definition above,
\[
\lambda = \frac{1}{4}\frac{2mgr}{\frac{1}{2m(r\theta')^2}} =
\frac{g}{r(\theta')^2}\: .
\]
In other words, assuming (without loss of generality) unit pendulum
length and unit velocity the bottom of the swing, $\lambda$ simply
represents the force of gravity.
Figure \ref{elastica-fam-k}, the curvature plots for solutions of the
elastica equation, thus also shows the velocity of the pendulum as a
function of time for various solutions of the motion of the
pendulum. For each value of $\lambda$, velocity ($\theta'$) is plotted
against time. The bottom of the swing (maximum velocity) is the center
of the plot, and $25$ units of time are plotted to either side.
Note that the period has a singularity at $\lambda = 0.25$. For
values of $\lambda$ greater than $0.25$, the pendulum swings back
and forth periodically, velocity reaching zero at the highest point in
the swing. For values less than $0.25$, the pendulum winds round
and round, slowing at the top of the cycle, but never actually
reaching zero. At exactly $0.25$, the kinetic energy at the
bottom is equal to the potential energy at the top, so the velocity is
zero at the top.
\begin{figure*}
\begin{center}
\includegraphics[scale=.9]{figs/pendex}
\caption{\label{pendulum-ex}Examples of pendulum analogy.}
\end{center}
\end{figure*}
In the kinetic analogy with the elastica, velocity of the pendulum
corresponds to curvature. Thus, values of $\lambda$ above
$\frac{1}{4}$ have inflection points. Figure \ref{pendulum-ex} shows
the relationship graphically, including the right angle made by the
pendulum at the top of its swing for the $\lambda = 0.5$ rectangular
elastica case.
Note that while $\lambda = 0.3027$ is a special value for the elastica
(it loops over itself and occupies finite space), there is no special
physical meaning in the pendulum analogy.
\section{Equilibrium of forces: Finite element approach}
When an elastic strip is in a minimal energy configuration, it is in
mechanical equilibrium, meaning that at all points along its length,
the sum of forces on that point is zero.
\begin{figure*}[tbh]
\begin{center}
\includegraphics{figs/chain}
\caption{\label{chain}Finite element approximation of elastica as
chain of struts and pivots.}
\end{center}
\end{figure*}
The standard derivation of the equations of mechanical equilibrium
from the theory of elasticity \cite{Love} expresses relations between
shear forces, tension (longitudinal) forces, and flexural couples. The
equations are straightforward, but the relevant physics
requires appeal to results in the theory of elasticity.
This section, therefore, presents a simplified finite-element model
with only bending and tension forces.
\begin{figure*}[tb]
\begin{center}
\includegraphics[scale=.8]{figs/strutfig}
\caption{\label{strutfig}Forces on a single strut.}
\end{center}
\end{figure*}
The finite element model is a chain consisting of alternating rigid
struts and pivots with torsion springs. Such a chain, with spring
torque indicated schematically at each pivot, is shown in Figure
\ref{chain}.
Each strut generates a tension
force at both ends, in the same direction as the strut's length. In
Figure~\ref{strutfig}, this tension force is indicated by the quantity
$T$. By convention, compression force is denoted with a negative sign.
\begin{figure*}[tb]
\begin{center}
\includegraphics[scale=.8]{figs/pivotfig}
\caption{\label{pivotfig}Forces on a single pivot.}
\end{center}
\end{figure*}
At each pivot is a spring generating forces tending to straighten out
the chain, i.e. towards a turning angle of zero. In
Figure~\ref{pivotfig}, the turning angle is denoted by $\theta$ and
the the force generated by the moment (or spring torque) by $M$. At
the endpoints, the force generated by the moment is normal to the
strut joining the center to the endpoint. According to Newton's third
law, the sum of all forces arising from the element must be zero, so
there is a balancing force at the center point equal to the vector sum
of the forces at the endpoints, pointing in the opposite direction.
The spring at the pivot point is assumed to be purely elastic, which
is to say that the moment (torque) is proportional to the turning
angle. Using $c$ for this constant of proportionality (which, in a
real physical model, would be Young's modulus of elasticity times the
second moment of the section of the beam about the axis of bending),
\begin{equation}
\label{pivot-elastic}
M_i = c\Delta \theta_i\: .
\end{equation}
\begin{figure*}[tb]
\begin{center}
\includegraphics[scale=.8]{figs/chainfig}
\caption{\label{chainfig}Forces on a single pivot point in a chain.}
\end{center}
\end{figure*}
As shown in Figure~\ref{chainfig}, when these primitive elements are
linked together in a chain, there are five force vectors acting on a
single pivot point: the tension force from each adjoining strut ($T_0$
and $-T_1$), the moment forces from the previous and next pivots in
the chain ($M_0$ and $M_2$), and the balancing force $-2M_1\cos
(\Delta \theta_1/2$).
In flatland, the equilibrium condition that these vector forces
sum to zero induces two scalar equations. When the axes are chosen to
be parallel and perpendicular to one of the struts, these equations have
a particularly simple form. Arbitrarily choosing the strut labeled
$1$ in Figure~\ref{chainfig}, the parallel component is
\begin{equation}
-T_1 + \cos\Delta \theta_1\;T_0 - \sin\Delta \theta_1(M_0 - M_1) = 0\:,
\end{equation}
\noindent and the perpendicular component is
\begin{equation}
\sin\Delta \theta_1\;T_0 + \cos\Delta \theta_1(M_0 - M_1) + M_2 - M_1 = 0\:.
\end{equation}
In the continuous limit of this finite element model, the individual
turning angles at the pivots $\Delta \theta_i$ tend to zero, justifying the
usual linear approximations $\sin\Delta \theta \approx \Delta \theta$ and
$\cos\Delta \theta \approx 1$. The result of using these approximations, and substituting
the elastic pivot assumption of Equation~\ref{pivot-elastic} is
\begin{equation}
-T_1 + T_0 - c\Delta \theta_1(\Delta \theta_0 - \Delta \theta_1) = 0\:,
\end{equation}
\begin{equation}
\Delta \theta_1\;T_0 + c(\Delta \theta_0 + \Delta \theta_2 - 2\Delta \theta_1) = 0\:.
\end{equation}
In the limit as $\Delta \theta$ approaches zero (as the number of links in
the chain increases), these difference equations become the
differential equations
\begin{equation}
\label{elastica-diffeq-fst}
-\frac{dT}{ds} + c\kappa\frac{d\kappa}{ds} = 0\:, \\
\end{equation}
\begin{equation}
\label{elastica-diffeq-snd}
\kappa T + c\frac{d^2\kappa}{ds^2} = 0\:.
\end{equation}
Equation \ref{elastica-diffeq-fst} is readily integrated, yielding
\begin{equation}
T = \frac{c}{2}\kappa^2 + C\:.
\end{equation}
Substituting this into Equation \ref{elastica-diffeq-snd}, and
dividing by the common constant $c$, yields a restatement of
Equation~\ref{elastica-kappa}.
\begin{equation}
\kappa'' + \frac{1}{2} \kappa^3 + C\kappa = 0\:.
\end{equation}
The fact that the constant $c$, related to the coefficient of
elasticity, factors out has the physical interpretation that the
elastica takes the same shape regardless of the value of this
coefficient, as long as it is consistent across the length of the
elastica.
It is satisfying but not surprising that the variational and the
mechanical equilibrium treatments of the elastica problem result in
the same equations. The theory of conservative systems states that at
a stable equilibrium, the energy is at a local minimum. In addition to
confirming the variational derivation, the mechanical equilibrium
equations prove a useful physical interpretation for the constant $C$
in terms of curvature and tension (longitudinal force). In particular,
if $C$ is zero, then the tension vanishes at inflection points.
\section{Solutions with length and endpoint constraints}
Now that we have the mathematical solution to the elastica equation,
including plots for a variety of values of the parameter $\lambda$, we
turn our attention to particular solutions corresponding to real elastic
strips under given constraints. Can real elastic strips take on all
values of $\lambda$, or is some of the solution space a mathematical
artifact?
Without loss of generality, assume that the length of the strip is
$1$, and that the endpoints lie on the horizontal axis (all other
cases can be recovered through scaling, rotation, and
translation). Three parameters remain: the tangent angles at the
endpoints with respect to the horizontal axis, and the distance
between the endpoints (chord length).
Figure \ref{lengraph1-fig} shows one linear section of this three
dimensional parameter space: the endpoint tangents are fixed at $\pi /
2$ from horizontal, and the chord length varies from $0$ to $1$.
\begin{figure*}
\begin{center}
\includegraphics[scale=.75]{figs/lengraph1}
\caption{\label{lengraph1-fig}Values of $\lambda$ for varying chord lengths.}
\end{center}
\end{figure*}
Several properties are apparent: $\lambda$ varies smoothly from about
$0.34$ to a maximum of $0.5$, then falls quickly to $0$, where the
solution is a precise half-circle, then increases asymptotically to
$0.25$ as the strip is pulled tighter and tighter towards its full
length. For these endpoint angles, the $\lambda$ maximum at $0.5$ also
corresponds to zero curvature (inflection points) at the endpoints.
\begin{figure*}
\begin{center}
\includegraphics[scale=.65]{figs/multi_mec}
\caption{\label{multimec-fig}Multiplicity of solutions for same tangent constraints.}
\end{center}
\end{figure*}
Note that the above discussion only reflects the \emph{principal
solution} of the elastica with respect to the constraints, i.e. the
one with the minimum total arc length from the original curve. The
principal solution always has at most one inflection point. The
periodic nature of the elastica implies that additional solutions may
also exist, with more than one inflection point. In general, these
additional solutions are not interesting for the purpose of generating
spline curves, but they are valid solutions. Figure \ref{multimec-fig}
shows a number of these solutions for the end constraints $\theta_0 =
\pi/4,\ \theta_1= \pi/4,\ \lambda = 1/2$. It is interesting that some
of these solutions lack symmetry even though the constraints
themselves are symmetrical.
\section{Elastica with general constraints}
Working towards the application of an elastica to splines, this
section considers the effect of general constraints on an
elastica. A goal is to develop physical intuition, so the constraints
are expressed in terms of forces.
There are four basic types of constraints. Each constraint, can either
rotate freely or fix the tangent angle of the elastica. Further, it
can either allow the elastica to slide freely, or can fix the relative
arc length of the elastica. A segment between two constraints that fix
relative arc length has fixed arc length, even if there are sliding
constraints in between. In all these types of constraints, the
position of the constraint is fixed, and the elastica is constrained
to go through that point.
A position constraint, freely rotating, and with the elastica free to
slide through the point, exerts only a normal force on the
elastica. If the tangent angle is fixed, the constraint exerts a normal
force on the elastica. Similarly, if the relative arc length is fixed,
the constraint exerts a longitudinal force.
\newtheorem{result}{Result}
\begin{result}
A normal force acting on an elastica preserves $G^2$-continuity across
the point where the force is applied, but creates a discontinuity in
the first derivative of curvature.
\end{result}
\begin{result}
A shear force (torque) acting on an elastica preserves
$G^1$-continuity of the elastica but creates a discontinuity in
curvature.
\end{result}
\begin{result}
A segment of an elastica with zero longitudinal force at the endpoints
takes the form of a rectangular elastica. Nonzero longitudinal force
induces values of $\lambda$ other than $0.5$.
\end{result}
Proofs of these results are not included here, as they are known in
the literature of the elastica. An accessible source for them, in the
context of interpolating splines, is Lee and Forsythe
\cite{ForsytheLee}.
\begin{figure*}
\begin{center}
\includegraphics[scale=.8]{figs/rectel}
\caption{\label{rectel-fig}Forces acting on rectangular elastica.}
\end{center}
\end{figure*}
Considering a full cycle of the rectangular elastica, as shown in
Figure \ref{rectel-fig}, it is clear that the longitudinal forces are
zero at the endpoints, because the curve is perpendicular to the force
$\mathbf{F}$ acting on it.
\section{Angle constraints}
The general constraint approach of the previous section yields a
specific result relevant to using the elastica as a spline.
\begin{result}
An elastica segment with an angle constraint at each endpoint, but
allowed to slide freely through those constraints, is a section cut from
the rectangular elastica.
\end{result}
Equivalently, the minimum energy curve (of arbitrary length) with
given endpoints and given angles at the endpoints is a section cut from
the rectangular elastica.
\begin{figure*}
\begin{center}
\includegraphics[scale=.8]{figs/mecrange}
\caption{\label{mecrange-fig}Angle-constrained minimum energy curves.}
\end{center}
\end{figure*}
Figure \ref{mecrange-fig} shows the full parameter space of
angle-constrained minimum energy curves. For sufficiently large
angles, no such solution exists; the envelope for which the MEC
exists is shown as well. This envelope was determined by implementing
a very robust numerical solver and measuring the boundary at which
the solver did not converge. As far as is known, there is no
analytical formula for the envelope's shape. Only half of the complete
parameter space is shown, to avoid clutter. The envelope is
symmetrical around the $\theta_{\mbox{\scriptsize \textit{left}}} =
\theta_{\mbox{\scriptsize \textit{right}}}$ diagonal.
\begin{figure*}
\begin{center}
\includegraphics[scale=.8]{figs/mecrangek}
\caption{\label{mecrangek-fig}Curvature plots, angle-constrained
minimum energy curves.}
\end{center}
\end{figure*}
Similarly, Figure \ref{mecrangek-fig} shows the curvature plots of all
angle-constrained MECs within their parameter space. It is worth
noting that, for small angles, the curvature is very nearly
linear with arc length.
Note that the curvatures at the endpoints are fully determined by
the angles at the endpoints.
\section{Elastica as an interpolating spline}
\label{mec-interp-sec}
The classical formulation of an interpolating spline is the minimum
energy curve passing through a sequence of control points.
Equivalently, it is an elastica passing through a sequence of
constraints, each of which rotates freely and allows the elastica to
slide freely through it.
Because each of these constraints exerts only a normal force, the
curve is $G^2$ continuous. Since the elastica is unbent past
either endpoint, the curvature at the endpoints is zero (this is also
known as the ``natural'' endpoint condition).
Each segment of the spline, bounded by two control points, is an
instance of an angle-constrained MEC as described in the previous
section.
To find an interpolating spline, it is necessary and
sufficient to find a global solution assigning a tangent angle to each
point so that the resulting curvatures are continuous across each
point (and zero at the endpoints).
\section{SI-MEC}
\label{simec-section}
To address some of the shortcomings of the MEC as a spline, Moreton
and S\'equin
proposed the Scale Invariant Minimum Energy Curve (SI-MEC) \cite{Moreton93}. It is
defined as the curve minimizing the functional
\begin{equation}
\label{simec-eq}
\mbox{SIMEC}[\kappa(s)] = l\int_0^l \kappa(s)^2\ ds\:.
\end{equation}
The SI-MEC derives its name from the fact that the functional is
invariant with respect to scaling of the curve. The MEC, by contrast,
is inversely proportional to the scale factor. However, this scaling
property has little, if any, physical meaning, as the resulting shape
is invariant to scaling of the entire physical system in either case.
The SI-MEC functional of Equation \ref{simec-eq} is equivalent to the
bending energy of a lamina with length normalized to unity. Thus, a
SI-MEC with angle constraints at the endpoints is equivalent to this
physical system: the endpoint constraints slide on a rail, but have
their angles fixed. The lamina is of fixed length and is not allowed
to slide through the endpoint constraints.
The SI-MEC will, for any given constraints, be shorter than the
corresponding MEC. Otherwise, the SI-MEC would have even lower MEC
energy than the MEC curve, which is impossible by definition.
Any elastica meeting the given angle constraints can be evaluated in
its MEC energy and its SI-MEC energy. In short, the MEC energy is the
bending energy when the chord length is held fixed and the arc length
varies, and the SI-MEC energy is the bending energy when the arc length
is held fixed.
\begin{figure*}
\begin{center}
\includegraphics[scale=.75]{figs/lengraph2}
\caption{\label{lengraph2-fig}MEC and SI-MEC energies for varying
lengths ($\pi / 2$).}
\end{center}
\end{figure*}
\begin{figure*}
\begin{center}
\includegraphics[scale=.75]{figs/lengraph3}
\caption{\label{lengraph3-fig}MEC and SI-MEC energies for varying
lengths ($\pi / 3$).}
\end{center}
\end{figure*}
Consider first the case of symmetrical angle constraints. Figure
\ref{lengraph2-fig} shows the MEC and SI-MEC energy for all elastica
meeting the angle constraints $(\pi/2, \pi/2)$, of chord length 0 to 1
(assuming unit arc length). The MEC energy is minimized, as expected,
for $\lambda = 0.5$, the rectangular elastica. The SI-MEC energy is
minimized for the circle, $\lambda = 0$. Figure \ref{lengraph3-fig}
shows a similar plot for angle constraints $(\pi/3, \pi/3)$, showing
that as the angles become less extreme, the difference in chord length
between MEC and SI-MEC becomes less pronounced.
\begin{figure*}
\begin{center}
\includegraphics[scale=0.8]{figs/simecplot}
\caption{\label{simecplot-fig}SI-MEC and MEC curve energy for $(\pi/2,
-\pi/2)$ angle constraints.}
\end{center}
\end{figure*}
Figure \ref{simecplot-fig} plots the MEC and SI-MEC energy as a function of
the chord length / arc length ratio, given angle constraints of
$(\pi/2, -\pi/2)$.
Along with each representative elastica curve is shown the continuation
of the elastica curve (with dashed lines), and the line of action (of
which the length is proportional to the force of action applied by the
constraints).
Note that the line of action for the curve minimizing the SI-MEC energy
is vertical. In fact, the line of action is always perpendicular to
the chord for a SI-MEC with angle constraints at the endpoints, as can
be deduced from the equilibrium of forces in the physical model. If
the action had nonzero component along the chord, the constraints
would slide closer together or farther apart along the rail.
Several interesting properties follow immediately from the vertical
line of action. Based on equilibrium of moments, the curvature varies
linearly along the chord length. For symmetrical angle constraints, the
curvature is also symmetrical, and thus must be equal at the two
endpoints. Thus, the curvature is constant, and the curve is a circle.
For the angle constraints $(\pi/2, -\pi/2)$, the SI-MEC is a segment of
the lemnoid elastica. This is also true for $\pi/2$ and the lemnoid
elastica angle ($\approx .71$ radians).
Further, all such angle-constrained SI-MEC curves have $\lambda \leq
0.5$, for the $\lambda > 0.5$ cases are single-valued in $x$ for any
given $y$ (where the line of action is drawn as vertical), and of
course both endpoints must lie on the $x$ axis.
Because curvature is linear along the chord length, for angle constraints
with absolute value $\leq \pi/2$, curvature is monotonic with arc
length. Contrast with the MEC, which has a curvature peak for all
symmetric and nearly-symmetric angle constraints, even as the angles
approach zero. As described in Section \ref{monotone-curvature-sec},
monotonic curvature is a desirable property for curve segments used in
interpolating splines.
The domain where the angle-constrained SI-MEC has a solution is about
double the size of the MEC's. In the symmetrical case, it reaches a
singularity for the full circle $(\pi, \pi)$ (exactly double that of
the MEC), and in the anti-symmetrical case the singularity is a full
figure-eight around the lemnoid elastica, approximately $(5.57,
-5.57)$, as opposed to the MEC's $(2.22, -2.22)$.
The SI-MEC has some problems when generalizing from angle
constraints at the endpoints to multiple position constraints. In
particular, it does not form an extensional spline.
\section{Real splines}
\begin{figure*}
\begin{center}
\includegraphics[width=5in]{figs/draftweight02}
\caption{\label{real-spline}A real spline.}
\end{center}
\end{figure*}
As mentioned, the MEC is a mathematical idealization
of a real wooden spline, the kind shown in Figure \ref{real-spline}
that has been used in shipbuilding for hundreds of years.
For small angles, the two agree closely. But for large angles, how
accurately does the MEC model a real spline? The fact that the MEC
does not exist for certain angle constraints is
particularly suspicious.
The key to this question is, I believe, tension forces and
friction. In a real spline, small tension forces are countered by
friction at the control points.
Indeed, shipbuilders use traditional spline weights in two different
ways. In one such way, the spline weight pushes against the side of
the spline, and the spline is free to slide with fairly low friction.
\begin{figure*}
\begin{center}
\includegraphics[width=4in]{figs/draftweight07}
\caption{\label{real-high-friction}High-friction spline constraints on a real
spline.}
\end{center}
\end{figure*}
Another common way to use spline weights is to place the point of the
pin on the spline curve, as shown in Figure \ref{real-high-friction}.
In this way, the weights effectively constrain the arc length as well
as the position, leaving only rotation entirely free.
It is understood that applying such constraints indiscriminately
yields curves of poor fairness (applying tension at the constraints
yields curvature discontinuity). Thus, designers use an iterative
technique to improve fairing. Once the curve is laid out roughly, each
weight is lifted, the spline makes small, fairly localized adjustments
when it is free from the weight, and then the weight is replaced, in
turn. Of course, if the resulting curve diverges too much from the
desired shape, the spline is pushed when the weight is replaced, or a
new weight is added.
Iterative local adjustment is reminiscent of successive
over-relaxation, a common numerical technique for solving global
systems of equations. Note that, as a weight is lifted, the total arc
length bounded by the two neighboring weights remains constant. Thus,
this process does not inevitably converge to the MEC solution, which
has longer arc length than a spline under some tension. Solutions
similar to the SI-MEC, or to circular arcs, remain stable.
\begin{figure*}
\begin{center}
\includegraphics[width=4in]{figs/matthew_real_spline}
\caption{\label{matthew-spline}Bricks as high friction constraints on a real
spline.}
\end{center}
\end{figure*}
Friction, then, seems an appealing technique for improving the
robustness of the pure MEC, and for making it conform to circles more
closely. Shipbuilders have built up a tradition for using splines
effectively. For the spline, they often prefer poplar or white pine
for the material, about a quarter inch thick (but thinner for more
detailed work involving smaller radii, or thicker for large boat
lines). For the weights, they usually prefer a pound or two so it
holds the spline firmly, a rounded shape resembling a duck or whale
for easier handling, and a configuration of the pin so that the weight
itself is separated from the spline. Even so, the traditional weight
has no different effect on the spline than would a brick (and simple
bricks are pressed into service by do-it-yourself boatbuilders, as
well, as in Figure \ref{matthew-spline}). Aside from shape of the
weight, though, all these factors affect the degree to which the real
spline models a MEC, especially the friction at the constraints.
A number of researchers have been tempted by the possibility of using
tension to make elastica-based splines better than the pure MEC.
However, there is no clear way to accomplish this. Real splines use
spline weight and wood thickness to set the free parameters, but on a
computer these parameters would either need to be set explicitly by
the user, or in some other way -- at worst, a computer representation
of physical quantity such as ``one pound weight'' or ``quarter inch
thick white pine'' are ``voodoo constants,'' which are well known for
limiting the scope of usefulness of a system.
All these problems are solved, though, when the elastica is left
behind, and the Euler spiral is used in its stead. As described in
Chapters \ref{hist-elast-chapter} and \ref{hist-euler-chapter},
dedicated to the history of these curves, their fates are closely
intertwined. In fact, the complete mathematical characterization of
both curves can be traced to the same 1744 treatise by Leonhard Euler.
But the details will wait for Chapter \ref{hist-elast-chapter}. The
next chapter explores the properties of $G^2$-continuous interpolating
splines in general, and quantifies exactly why the Euler spiral is
optimal as a curve primitive for interpolation.