# How to plot a bicycle with square wheels

I want to draw a bicycle with square wheels similar to this picture, but I can't plot the trajectory along the curve. enter link description here

(*https://pastebin.com/3UbbfG6W*)

corner[x_] :=
Module[{θ = N@SawtoothWave[{0, 2 Pi}, x/(2 Pi)]},
Piecewise[{{{-1 - Cos[θ], Sin[θ]},
0 <= θ < Pi/2}, {{-Sqrt[2] Cos[θ - Pi/4],
Sqrt[2] Sin[θ - Pi/4]},
Pi/2 <= θ < Pi}, {{1 + Sin[θ - Pi],
Cos[θ - Pi]},
Pi <= θ < 3 Pi/2}, {{2, 0}, θ >=
3 Pi/2}}] + {4 Floor[x/(2 Pi)], 0}];

frame[θ_] :=
Module[{t}, t = If[θ <= 3.5 Pi, θ, 7 Pi - θ];
Show[Graphics[{RGBColor[0.8, 0.3, 0.2],
Polygon[{corner[t], corner[t + Pi/2] - {1, 0},
corner[t + Pi] - {2, 0}, corner[t + 3 Pi/2] - {3, 0}}], Black,
PointSize[Medium], Point[corner[t]]},
PlotRange -> {{-2, 5}, {-0.1, 1.5}}],
If[θ >= 0.5 Pi,
ParametricPlot[corner[x], {x, t, 0.5 Pi+ 0},
PlotStyle -> Black], Graphics[]]]];

Manipulate[frame[θ], {θ, 0, 7 Pi}]
• See Stan Wagon, Mathematica in Action, 3e, electronic version, RollingSquare[]. – Michael E2 Jan 16 at 13:31
• I haven't bought this paper book with CD at present,How can I get the RollingSquare[] code content in the DVD of this e-book? – A little mouse on the pampas Jan 17 at 1:07
• Thanks a lot,I've found a link to download this book's CD-ROM from mathematica.stackexchange.com/questions/87284/… – A little mouse on the pampas Jan 17 at 1:55
• You're welcome. BTW, I think both the current answers are wrong. For one thing the trajectory of a corner of a square should have a 90º angle at the point of contact with the catenary surface it rolls on. (The trajectory of a point is normal to the surface at the point of contact, and the surface has a 90º corner where the vertex of the square contacts the surface.) – Michael E2 Jan 17 at 2:13
• Before I clicked the book link I thought I was about to discover that "Stan Wagon" was the name for a wagon with square wheels. – IPoiler Jan 21 at 22:31

## 6 Answers

This question is too interesting to resist, so I'll talk about how to analyze the problem.

Take a look at sketch above. It describes an arbitrary moment during the rolling. From the kinematics view, $$P$$ is the "instant center of rotation". From the energy view, the square's center of mass $$O$$ keeps its height, thus the potential of the square doesn't change, means it must be in balance. Either way we arrive to the same conclusion that $$\overline{OP}$$ is perpendicular to the trajectory of $$O$$ (the horizontal red dash-line).

Suppose the side length of the square is $$2$$, and the equation of the questioned curve is $$\boldsymbol{r}(s):=\left(x(s),y(s)\right)$$, where parameter $$s$$ is the length of $$\overline{CP}$$, which should be equal to the arc-length of $$\overset{\mmlToken{mo}{⏜}}{C'P\,}$$ due to the slipping-less rolling. It's straightforward to see tangent vector $$\frac{\mathrm{d}\boldsymbol{r}}{\mathrm{d}s}$$ at $$P$$ is parallel to $$\overline{CP}$$, so ($$\dot{F}$$ is a short-form for $$\frac{\mathrm{d}F}{\mathrm{d}s}$$ for any $$F$$)

$$\frac{-\dot{y}}{\dot{x}}=\frac{\mathrm{length}_\overset{\rightharpoonup}{CP}}{\mathrm{length}_\overset{\rightharpoonup}{OC}}=s\implies s\dot{x}+\dot{y}=0\;\text{.}$$

Additionally, because $$s$$ is an arc-length parameter, we have

$$\dot{x}^2+\dot{y}^2=1\;\text{.}$$

We setup the coordinate frame so the trajectory of $$O$$ lies on x-axis and $$C'$$ lies on y-axis. Solving the system is a one-liner

DSolve[{
s x'[s] + y'[s] == 0,
x'[s]^2 + y'[s]^2 == 1,
x[0] == 0,
y[0] == -1
}, {x, y}, s]

$$\left\{\left\{x\to-\sinh^{-1}(s), y\to\sqrt{s^2+1}-2\right\}, \left\{x\to\sinh^{-1}(s), y\to-\sqrt{s^2+1}\right\}\right\}$$

Selecting the solution with positive $$\dot{x}$$, we have

\left\{ \begin{align} x&=\sinh^{-1}(s) \\ y&=-\sqrt{s^2+1} \\ \end{align} \right.\;\text{,}

or

Block[{\$Assumptions = x \[Element] Reals},
y == -Sqrt[1 + s^2] /. Solve[x == ArcSinh[s], s] // FullSimplify
]

i.e.

$$y=-\cosh(x)$$

At last the animation:

ClearAll[catenaryGround, origcube, point, perp]
catenaryGround =
Plot[-Cosh[x], {x, -ArcSinh[1], ArcSinh[1]}, PlotRange -> All,
AspectRatio -> Automatic] // Cases[#, _Line, Infinity] & //
First;
origcube = {
{EdgeForm[GrayLevel[0.3]], FaceForm[GrayLevel[0.9]], Cuboid[{-1, -1}, {1, 1}]},
{GrayLevel[0.3], Line[{{0, 0}, {0, -1}}]}
};
point = {EdgeForm[{Hue[0., 1., 0.66], Thick}], FaceForm[GrayLevel[0.9]], Disk[{0, 0}, .04]};
perp = Line[{{1, 0}, {1, 1}, {0, 1}}];
ClearAll[cubeTF]
cubeTF[x_] := RotationTransform[ArcTan[1, -Sinh[x]]] /* TranslationTransform[{x, 0}]
ClearAll[periodLen, totalPeriod]
periodLen = 2 ArcSinh[1];
totalPeriod = 5;
DynamicModule[{period = 1, xshift, xC = -(periodLen/2), x, tf, center, contact, bottom},
DynamicWrapper[
Deploy@Graphics[{
{EdgeForm[GrayLevel[0.3]], FaceForm[GrayLevel[0.9]], Translate[FilledCurve@catenaryGround, {(# - 1) periodLen, 0} & /@ Range[totalPeriod]]}
, Dynamic@GeometricTransformation[origcube, tf]
, {Hue[0., 1., 0.66], Dashed, InfiniteLine[{0, 0}, {1, 0}]}
, {Hue[0.54, 1., 0.66], Dashed, Line[Dynamic@{center, contact}]}
, {Hue[0.54, 1., 0.66],
Dynamic@GeometricTransformation[perp, RightComposition[
ScalingTransform[1/8 {1, 1}],
RotationTransform[Pi/2 (<|-1 -> 2, 0 -> 2, 1 -> 3|>@Sign[x])],
TranslationTransform[center]
]]}
, {GrayLevel[0], Dynamic@GeometricTransformation[perp, RightComposition[
ScalingTransform[1/10 {1, 1}],
RotationTransform[Pi/2 (<|-1 -> 1, 0 -> 1, 1 -> 0|>@Sign[x])],
TranslationTransform[{0, -1}], tf
]]}
, {Black, AbsoluteThickness[4], CapForm[None], Line[Dynamic@{bottom, contact}]}
, {Black, AbsoluteThickness[4], CapForm[None],
Line@Dynamic[Function[s, {ArcSinh[s] + xshift, -Sqrt[1 + s^2]}] /@ N[Rescale[Rescale[Range[100]], {0, 1}, Sort@{0, Sinh[x]}]]]
}
, Dynamic@Translate[point, {center, contact}]
, Text[Style["O", Italic, 12], Dynamic[center], {0, -1}]
, Text[Style["P", Italic, 12], Dynamic[contact], Dynamic@{-Sign[x] 2, 0}]
, Text[Style["C", Italic, 12], Dynamic[bottom], Dynamic@{Sign[x] 2, -1}]
}
, ImageSize -> 800, PlotRange -> {{-1, 2 totalPeriod - 1} periodLen/2 + {-1, 1} Sqrt[2], {-1, 1} Sqrt[2]}, PlotRangePadding -> None
]
,
xC = -Cos[2 Clock[Pi, 10]] // Rescale[#, {-1, 1}, {-1, 2 totalPeriod - 1} periodLen/2] &
; center = {xC, 0}
; period = Round[xC/periodLen] + 1
; xshift = (period - 1) periodLen
; x = xC - xshift
; contact = {x, -Cosh[x]} + {xshift, 0}
; tf = cubeTF[x] /* TranslationTransform[{xshift, 0}]
; bottom = tf@{0, -1}
]
]

• Another correct answer! And such a beautiful animation! – Michael E2 Jan 21 at 12:34
• @MichaelE2 Thanks! Glad to hear my animation is liked! – Silvia Jan 22 at 15:26

Your curves show a bit strange case when central point keeps the same altitude. However, such situation can be simulated by the following way:

R = 1;
ϕ = π/4;

center = {#, 1} &;

crn = {center@# + {-R Cos@(# + ϕ), R Sin@(# + ϕ)},
center@# + {-R Cos@(# + ϕ + π/2),
R Sin@(# + ϕ + π/2)},
center@# + {-R Cos@(# + ϕ + π),
R Sin@(# + ϕ + π)},
center@# + {-R Cos@(# + ϕ + (3 π)/2),
R Sin@(# + ϕ + (3 π)/2)}} &;

trjc = Table[center@f, {f, 0, 5 π, π/10}];
trjcorn = Table[crn@f, {f, 0, 5 π, π/10}];
Manipulate[
Show[
ListPlot[{trjc, trjcorn[[All, 1]], trjcorn[[All, 2]],
trjcorn[[All, 3]], trjcorn[[All, 4]]}, AspectRatio -> 2/(5 π),
PlotStyle -> {Gray, Red, Orange, Green, Blue}, ImageSize -> 800,
Joined -> True],
Graphics[{
Black, PointSize[0.01], Point@center@f,
PointSize[0.01], Red, Point[crn[f][[1]]], Orange,
Point[crn[f][[2]]], Green, Point[crn[f][[3]]], Blue,
Point[crn[f][[4]]],
Black, Line@crn@f, Line[crn[f][[{1, 4}]]]},
PlotRange -> {{-π/2, 5 π}, All}]
]
, {f, 0, 5 π, π/10}]

Thus, the trjc and trjcorn contain points of trajectories of center and corners correspondingly

• I need to scroll forward along the curve drawn by the function Plot[-(Cosh[Mod[x, 2 Log[Sqrt[2] + 1], -Log[Sqrt[2] + 1.]]]) + Sqrt[ 2], {x, 0, 10}],that's the point. – A little mouse on the pampas Jan 17 at 5:26
• @Gowiththewind, and what's the problem with this? Just replace my ListPlot block with your Plot and enjoy fun :) – Rom38 Jan 17 at 5:58
• Thanks a lot,I have basically solved this problem in the following answers. – A little mouse on the pampas Jan 17 at 8:57
• The problem is that the vertices of a square-wheeled bicycle do not trace cycloids as shown in this answer. What this answer shows is the trajectory of the square on a circular bicycle wheel, which is entirely different. Replacing the ListPlot with @Go's Plot produces even worse results. – Michael E2 Jan 18 at 11:32

To make it a bit shorter

n = 4 (*number of corners*)
crn = {Cos[#/n 2 Pi], Sin[#/n 2 Pi]} & /@ Range[n];
cen = {0, 1};
v[t_] = {t, 0};(*linear velocity*)
w[t_] = -0.2 t 2 Pi;(*angular velocity*)
tmax = 10; dt = 0.1;
trj = Table[(cen + v[t] + RotationMatrix[w[t]].#) & /@ crn, {t, 0, tmax, dt}];

ListAnimate[Table[Graphics[{Gray, Polygon[t],
Black, Point[Mean[t]],(*Centre*)
Table[{Hue[i/n], Line[trj[[All, i]]]}, {i, n}],
Black, Dashed, Line[Table[Mean[t], {t, trj}]]}]
,{t, trj}]]

You can play with different polygon and velocity. For an arbitrary polygon, you have to define its corners yourself in crn.

## Cycling on an arbitrary path

n = 4 (*number of corners*)
crn = {Cos[#/n 2 Pi], Sin[#/n 2 Pi]} & /@ Range[n];
cen = {0, 1};
v[t_] = {t, -(Cosh[Mod[t, 2 Log[Sqrt[2] + 1], -Log[Sqrt[2] + 1.]]]) +   Sqrt[2]};
(*parametric linear velocity*)
w[t_] = -0.2 t 2 Pi;(*angular velocity*)
tmax = 10; dt = 0.1;
trj = Table[(cen + v[t] + RotationMatrix[w[t]].#) & /@ crn, {t, 0, tmax, dt}];

ListAnimate[Table[Graphics[{Gray,Polygon[trj[[1]]], Polygon[trj[[-1]]],
Polygon[t], Black, Point[Mean[t]],(*Centre*)
Table[{Hue[i/n], Line[trj[[All, i]]]}, {i, n}], Black, Dashed,
Line[Table[Mean[t], {t, trj}]]}, PlotRange -> {{-1, 11}, {0, 3}}], {t, trj}]]

Note that I use a PlotRange here to fix the frame. Otherwise your animation might be shaky. A related post is Why is my GIF shaky?

• I need to scroll forward along the curve drawn by the function Plot[-(Cosh[Mod[x, 2 Log[Sqrt[2] + 1], -Log[Sqrt[2] + 1.]]]) + Sqrt[ 2], {x, 0, 10}] – A little mouse on the pampas Jan 17 at 5:24
• With your friendly help, I have basically solved this problem in the following answers.Thank you very much. – A little mouse on the pampas Jan 17 at 9:07

By the following proof, we can know that the horizontal ordinate of the center of the rectangle is the same as that of the tangent point.

(x + (1 - ArcLength[-Cosh[t] + Sqrt[2], {t, 0, x}]) Cos[
ArcTan[D[-(Cosh[x]) + Sqrt[2],
x]]] + (x - (1 + ArcLength[-Cosh[t] + Sqrt[2], {t, 0, x}]) Cos[
ArcTan[D[-(Cosh[x]) + Sqrt[2], x]]] -
2*Sin[ArcTan[D[-(Cosh[x]) + Sqrt[2], x]]]) // Simplify)

The result 2x that we get shows this.

Animate[Show[
Graphics[{Rotate[
Translate[{EdgeForm[Black], LightGray,
Rectangle[{-1, Sqrt[2] - 1}, {1, Sqrt[2] + 1}], Black,
Line[{{-1, Sqrt[2] - 1}, {1, Sqrt[2] + 1}}],
Line[{{-1, Sqrt[2] + 1}, {1, Sqrt[2] - 1}}]}, {x, 0}],
Evaluate[
Simplify[
ArcTan[D[-(Cosh[
Mod[t, 2 Log[Sqrt[2] + 1], -Log[Sqrt[2] + 1.]]]) + Sqrt[
2], t]], t >= 0] /. t -> x], {x, Sqrt[2]}]}, Axes -> True,
AxesOrigin -> {0, 0}],
Plot[-(Cosh[Mod[t, 2 Log[Sqrt[2] + 1], -Log[Sqrt[2] + 1.]]]) + Sqrt[
2], {t, 0.001, x}, AxesOrigin -> {0, 0}], AxesOrigin -> {0, 0},
PlotRange -> {{0, 10}, {0, 3}}], {x, 0, 10},
AnimationRunning -> False, AnimationDirection -> ForwardBackward,
DefaultDuration -> 2]

• (+1) Finally, a correct answer! :) – Michael E2 Jan 18 at 2:53

Silvia already gave a pretty rigorous derivation (compare this with the treatment by Hall and Wagon), so let me show how to plot the desired trajectory of the rolling polygon's corner. One could certainly modify the code I gave here for this, but I will instead adapt this solution I previously wrote in OpenGL to Mathematica:

With[{n = 4}, (* number of sides *)
DynamicModule[{ang = π/n, apo, ha, pr, sa, traj, th, road},
apo = Cos[ang]; ha = ArcSinh[Tan[ang]];
pr = N[apo ha, 35]; sa = π/2 - ang; (* some extra precision needed *)

th[x_] := -sa - Gudermannian[ha Mod[x/pr, 2, -1]] + 2 ang Floor[(pr - x)/(2 pr)];
road[x_] := -apo Cosh[ha Mod[x/pr, 2, -1]];

traj = Table[{u, {u, 0} + AngleVector[th[u]]}, {u, -6 pr, 6 pr, pr/5}];

Manipulate[Show[Plot[road[x], {x, -6 pr, 6 pr}, AspectRatio -> Automatic,
Axes -> None, Background -> RGBColor[0., 0.169, 0.212],
Frame -> True, PlotRange -> {-2, 2},
PlotStyle -> RGBColor[0.149, 0.545, 0.824]],
Graphics[{{Gray, InfiniteLine[{0, 0}, {1, 0}]},
{Directive[FaceForm[],
EdgeForm[RGBColor[0.522, 0.6, 0.]]],
RegularPolygon[{t, 0}, {1, th[t]}, n]},
{Directive[AbsolutePointSize[5],
RGBColor[0.863, 0.196, 0.184]],
Point[Select[traj, First[#] <= t &][[All, -1]]]},
{Directive[AbsolutePointSize[3],
RGBColor[0.827, 0.212, 0.51]],
Point[{t, 0}]}}]],
{t, -6 pr, 6 pr}, SaveDefinitions -> True]]]

Here's how it looks for n = 6:

• Thanks for the Hall and Wagon reference! The Fourier series section is especially interesting! – Silvia Jan 22 at 15:28
• @Silvia, you're welcome; the paper is definitely a fascinating read. – J. M.'s discontentment Jan 24 at 2:53

This code is written by Shenlan

rotationPoint[φ_] := Line[{Cos[#] - φ, Sin[#] + 1} & /@ Range[φ, φ + 2 π, π/2]]
f[φ_] := 1 - Sqrt[2]/2 Sec[Mod[Pi/4 - φ, Pi/2, -Pi/4]]
Manipulate[
Show[{
Graphics[
{{Red, Circle[{φ, 1}]}, rotationPoint[-φ]},
PlotRange -> {{0, 2 Pi}, {0, 2}}],
Plot[f[φ1], {φ1, 0, φ}]}], {φ, 0.001, 2 Pi, 0.1},
TrackedSymbols :> {φ}
]

• This answer is wrong. (1) The vertices do not trace cycloids, which is what they do here. (2) The secant graph is approximately right, but zoom in shows it is not exactly right; the actual curve the square rolls on is a catenary. -- Your other answer is right. – Michael E2 Jan 18 at 11:40
• This is to use the first order approximation of catenary to solve this problem. The result is acceptable. Thank you for your guidance. – A little mouse on the pampas Jan 18 at 11:43
• The square's motion is not even approximately correct. Which answer to accept is your prerogative, of course, but accepting such a wrong solution does not seem helpful to others. – Michael E2 Jan 18 at 11:51
• I accept your suggestion and change my choice.Thanks a lot. – A little mouse on the pampas Jan 18 at 12:11