Proving that $x := x - yt$; $y := y + xt$ traces an ellipse

Question

I was investigating the iterated function formed by the following JavaScript code, trying to find an invariant that can prove that the system doesn't slowly diverge:

// Every frame:
x -= y * t
y += x * t

Or, in mathematical notation, $$ \left\{\begin{aligned} x' &= x - y t \\ y' &= y + x' t = y + x t - y t^2 \end{aligned}\right. $$

From what I can see empirically, for $-2<t<2$, the system seems to trace an ellipse in "state space", with a different shape depending on $t$—approaching a circle when $t$ is close to $0$, and approaching a pair of lines when $|t|$ close to $2$.

However, I've been stuck on how to prove that the system traces an ellipse—skewing/shearing and scaling never quite got me a true circle. Any ideas?

Answer 1

Substitute into the second equation:

$$ x^{\prime}=x-t y, \quad y^{\prime}=y+t(x-t y)=y+t x-t^2 y . $$

Note it does not preserve cirdes: If $R=x^2+y^2$, then

$$ \begin{aligned} x^{\prime 2}+y^{\prime 2} & =(x-t y)^2+(y+t(x-t y))^2 \\ & =x^2-2 t x y+t^2 y^2+y^2+2 t x y+t^2 x^2+O\left(t^2\right) \\ & =x^2+y^2+O\left(t^2\right) \end{aligned} $$

So although $x^2+y^2$ is unchanged to first order, the linear map

$$ A=\left(\begin{array}{cc} 1 & -t \\ t & 1 \end{array}\right) $$

distorts shapes through

$$ A^T A=\left(\begin{array}{cc} 1 & t \\ -t & 1 \end{array}\right)\left(\begin{array}{cc} 1 & -t \\ t & 1 \end{array}\right)=\left(\begin{array}{cc} 1 & -2 t \\ -2 t & 1 \end{array}\right)+O\left(t^2\right), $$

and

$$ (x, y)\left(\begin{array}{cc} 1 & -2 t \\ -2 t & 1 \end{array}\right)\binom{x}{y}=x^2+y^2-4 t x y, $$

whose mixed $x y$-term is the shear encoded in the ellipse. To account for that shear, consider instead

$$ Q=x^2-t x y+y^2 $$

Substituting $x^{\prime}, y^{\prime}$ back into $Q$ gives

$$ Q\left(x^{\prime}, y^{\prime}\right)=Q(x, y), $$

so it preserves ellipses

$$ x^2-t x y+y^2=Q\left(x_0, y_0\right) . $$

---

Conversely, Start with the ellipse

$$ Q(x, y)=x^2-t x y+y^2=C . $$

Its gradient is

$$ \nabla Q=(2 x-t y, 2 y-t x) . $$

A motion stays on the ellipse exactly when

$$ \frac{d}{d t} Q(x(t), y(t))=0 $$

so

$$ (2 x-t y) \dot{x}+(2 y-t x) \dot{y}=0 . $$

Solve this linear condition for the family of allowed directions. Rearrange:

$$ (2 x-t y) \dot{x}=-(2 y-t x) \dot{y} . $$

A simple choice that satisfies this relation is

$$ \dot{x}=-(2 y-t x), \quad \dot{y}=(2 x-t y) . $$

Now drop the factor 2 into parameter $t$ (scaling both components leaves direction unchanged):

$$ \dot{x}=-\left(y-\frac{t}{2} x\right), \quad \dot{y}=x-\frac{t}{2} y . $$

For small $t$, the dominant terms are

$$ \dot{x} \approx-y, \quad \dot{y} \approx x, $$

Taking an explicit Euler step gives

$$ x^{\prime}=x-t y, \quad y^{\prime}=y+t x^{\prime}, $$

which is exactly the code.

---

And you may ask yourself,
How did I get here?
This is not my beautiful ellipse in standard form.
This is not my principal-axis frame.

And you may ask yourself,
Well… how do I diagonalize this?

Starting with

$$ Q(x, y)=x^2-t x y+y^2=C, $$

write

$$ \begin{gathered} x=\frac{u+v}{\sqrt{2}}, \quad y=\frac{u-v}{\sqrt{2}} \\ x^2-t x y+y^2=\frac{1}{2}\left[(u+v)^2-t(u+v)(u-v)+(u-v)^2\right] \\ =\frac{1}{2}\left[u^2+2 u v+v^2-t\left(u^2-v^2\right)+u^2-2 u v+v^2\right] \\ =\frac{1}{2}\left[(2-t) u^2+(2+t) v^2\right] \\ =\left(1-\frac{t}{2}\right) u^2+\left(1+\frac{t}{2}\right) v^2 . \end{gathered} $$

Hence

$$ \begin{gathered} x^2-t x y+y^2=C \Longleftrightarrow\left(1-\frac{t}{2}\right) u^2+\left(1+\frac{t}{2}\right) v^2=C \\ \Longleftrightarrow \frac{u^2}{C /(1-t / 2)}+\frac{v^2}{C /(1+t / 2)}=1 . \end{gathered} $$

Answer 2

Put it in matrix form, and as derived in the OP, we have

$\begin{bmatrix} x(k+1) \\ y(k+1) \end{bmatrix} = \begin{bmatrix} 1 && - t \\ t && 1 - t^2 \end{bmatrix} \begin{bmatrix} x(k) \\ y(k) \end{bmatrix} $

where $t$ is some real number.

Define $V(k) = \begin{bmatrix} x(k) \\ y(k) \end{bmatrix} $

and $A = \begin{bmatrix} 1 && - t \\ t && 1 - t^2 \end{bmatrix} $

Then,

$ V(k+1) = A V(k) $

And we have $V(0) = \begin{bmatrix} x(0) \\ y(0) \end{bmatrix}$ given.

The solution of this system is

$ V(k) = A^k V(0) $

To calculate $A^k$ we need to diagonalize it. Diagonalization means writing $A$ in the following form:

$ A = P D P^{-1} $

where $D$ is a diagonal matrix, or a Jordan form matrix.

Once this is done, then it is easy to see that

$ A^k = P D^k P^{-1} $

To diagonalize $A$ we need to find its eigenvalues, and its eigenvectors.

The characteristic polynomial is

$ c ( \lambda) = \det( \lambda I - A ) = (\lambda - 1) (\lambda -1 + t^2 ) + t^2 $

From which,

$ c(\lambda) = \lambda^2 + (t^2 - 2) \lambda + 1 $

The next step is to set $c(\lambda) = 0$ and find the roots of this quadratic. Let the two roots be $\lambda_1$ and $\lambda_2$

$\lambda_1$ and $\lambda_2$ are the diagonal elements of matrix $D$.

Matrix $P$ has $2$ columns, which are eigenvectors corresponding to $\lambda_1$ and $\lambda_2$ respectively.

These are complicated computations, and you're encouraged to use library routines for them. These routines are available for example in MATLAB, Mathematica, etc.

Assuming we now have matrices $P$ and $D$, then

$ V(k) = P D^k P^{-1} V(0) $

So that

$ P^{-1} V(k) = D^k P^{-1} V(0) $

If we define

$ W(k) = P^{-1} V(k) $

Then

$ W(k+1) = D^k W(0) $

And this means that $V(k)$ (when expressed in the two-dimensional basis whose vectors are the columns of matrix $P$), is being scaled non-uniformly and the scale factors are determined by the diagonal elements of matrix $D$.

Edit: From the solution of the OP, it was shown that for $ t \in [-2, 2]$, the eigenvalues of $A$ are a conjugate complex pair of unit length. Consider now, the square of the length of vector $W(k)$, it is given by

$ W^{*T}(k) W(k) = W^{*T}(0) (D^{*k}) ( D^k ) W(0) = W^{*T}(0) W(0) = V(0)^T (P^{-1})^{*T} (P^{-1}) V(0) = c$

for some constant $c$. Thus $W(k)$ has an invariant length if $t \in [-2,2]$. Substituting the expression for $W(k)$ in the above equation gives

$ W^{*T}(k) W(k) = V(k)^T G V(k) = c \hspace{20pt}(1)$

where $G = (P^{-1})^{*T} (P^{-1} ) $ and is some real matrix.

Equation (1) is that of an ellipse because it can shown that matrix $G$ is a real symmetric positive definite matrix. The eigen-structure of matrix $G$ (i.e. its eigenvalues and eigenvectors), in addition to the constant $c$ determine the lengths of the axes of the ellipse, and their orientation relative to the coordinate axes in the $xy$ plane.

Example:

Let take $t = 1.2$, which is in $[-2,2]$. Then

$A = \begin{bmatrix} 1 && -1.2 \\ 1.2 && -0.44 \end{bmatrix} $

$A$ can be diagonalized and written as

$ A = P D P^{-1} $

with

$ D = \begin{bmatrix} 0.28 - j 0.96 && 0 \\ 0 && 0.28 + j 0.96 \end{bmatrix}$

and

$ P = \begin{bmatrix} 0.6 - j 0.8 && 0.6 + j 0.8 \\ 1 && 1 \end{bmatrix} $

And now we have $G$ given by

$ G = (P^{-1})^{*T} P^{-1} = \begin{bmatrix} 0.78125 && - 0.46875 \\ - 0.46875 && 0.78125 \end{bmatrix} $

The constant $c$ is given by

$ c = W(0)^T W(0) = V(0)^T G V(0) $

Suppose $V(0) = [5, 3]^T $ , then

$c = 12.5 $

Further, matrix $G$ can be diagonalized and written as

$ G = R \Lambda R^T $

where

$ \Lambda = \begin{bmatrix} 0.3125 && 0 \\ 0 && 1.25 \end{bmatrix} $

And $R$ is a rotation matrix by $45^\circ$, i.e.

$ R = \dfrac{1}{\sqrt{2}} \begin{bmatrix} 1 && -1 \\ 1 && 1 \end{bmatrix} $

And this means, the resulting ellipse will have a major semi-axis of

$ a = \dfrac{ \sqrt{c} }{\sqrt{\Lambda_{11}}} = \dfrac{\sqrt{12.5}}{\sqrt{0.3125}} = 6.324555 $

and a minor semi-axis of

$ b = \dfrac{ \sqrt{c} }{ \sqrt{\Lambda_{22}}} = \dfrac{\sqrt{12.5}}{\sqrt{ 1.25}} = 3.162277 $

The ratio of major to minor axes lengths is $2:1$.

I simulated this system, and got the following figure

Answer 3

Given the system defined by

$$ \left\{\begin{aligned} x' &:= x - y t \\ y' &:= y + x' t = y + x t - y t^2 \end{aligned}\right. \tag1 $$

the question is

However, I've been stuck on how to prove that the system traces an ellipse—skewing/shearing and scaling never quite got me a true circle. Any ideas?

This old problem would be interesting to trace the history of.

Define the quadratic form

$$ Q(x,y) := x^2 - txy + y^2 \tag2 $$

and verify by elementary algebra that

$$ Q(x,y) = Q_(x',-y) = Q(x',y'). \tag3 $$

Generically, the locus of $Q(x,y)=C$ is a centrosymmetric conic section centered at the origin while the locus of $Q(x,-y)=C$ is the same curve reflected about the $y$ or $x$-axis. Thus, given any point $\,(x,y)\,$ on the first curve, the shear map $\,(x,y)\mapsto(x',y)\,$ maps it to a point on the second curvee and then the map $\,(x',y)\mapsto(x',y')\,$ maps that to a point on the first curve. Now the combined map $\,(x,y)\mapsto(x',y')\,$ when iterated produces a sequence of points on the conic section.

Note that if $\,t=0\,$ then the maps are identity maps and in this special case the sequence of points is a constant sequence. For some other special values of $t$ the sequence is periodic, but otherwise the sequence of points is dense in the curve.

Answer 4

Building on Hosam Hajeer's answer:

Define $A = \begin{bmatrix} 1 & - t \\ t & 1 - t^2 \end{bmatrix} $ to be the matrix controlling the system's evolution through time; that is, for some starting value $v_0 = \begin{bmatrix} x_o \\ y_o \end{bmatrix}$, the system's state after $k$ iterations is $\begin{bmatrix} x_k \\ y_k \end{bmatrix} = A^k v_0 $

To show that the system is periodic (though, not necessarily an ellipse), we need to find some value of $k$ for which $A^k = I$.

Diagonalizing the matrix $A$ (thanks Hosam for the tip!), we get $A = P D P^{-1}$ and $A^k = P D^k P^{-1}$, where:

$$ D = \begin{bmatrix} \lambda_1 & 0 \\ 0 & \lambda_2 \end{bmatrix} $$

Where $\lambda_{1,2}$ are the roots of the characteristic polynomial:

$$c(\lambda) = \det(\lambda I - A) = \lambda^2 + (t^2 - 2) \lambda + 1 = 0$$

Applying the quadratic formula, we get:

$$ \lambda_{1,2} = \frac{2 - t^2 \pm \sqrt{t^4 - 4t^2}}{2} = \frac{2 - t^2 \pm t\sqrt{t^2 - 4}}{2} $$

Since $t^2 < 4$, we know that the discriminant is negative and thus $\lambda_{1,2}$ is a pair of complex conjugates, $$ \lambda_{1,2} = \frac{2 - t^2}{2} \pm \frac{t\sqrt{4 - t^2}}{2} i $$

Expressing $\lambda$ in polar form, we have $\lambda_{1,2} = r e^{\pm i\theta}$ for $r = |\lambda_{1,2}|$ and some angle $\theta$. Computing $r$, we get: $$ r = \left(\frac{2 - t^2}{2}\right)^2 + \left(\frac{t\sqrt{4 - t^2}}{2}\right)^2 = \frac{(4 - 4t^2 + t^4) + (4t^2 - t^4)}{4} = 1 $$

Therefore, $\lambda_{1,2} = e^{\pm i\theta}$.

Setting $k = \frac{2\pi}{\theta}$ and evaluating the original equation, we have:

$A^k = P \begin{bmatrix} \lambda_1 & 0 \\ 0 & \lambda_2 \end{bmatrix} ^ k P^{-1} = P \begin{bmatrix} e^{i\theta k} & 0 \\ 0 & e^{-i\theta k} \end{bmatrix} P^{-1} = P \begin{bmatrix} e^{i 2\pi} & 0 \\ 0 & e^{-i 2\pi} \end{bmatrix} P^{-1} = P I P^{-1} = I$

Therefore, the system described in the OP is periodic. $\blacksquare$

Answer 5

Let ($x_n$, $y_n$) denote the coordinates after $n$ iterations. Then the recurrence relation is:

$$(x_{n+1}, y_{n+1}) = (x_n - ty_n, y_n + tx_{n+1})$$ $$(x_{n+1}, y_{n+1}) = (x_n - ty_n, (1 - t^2)y_n + tx_n)$$

Now, assume that all points lie on a curve described by the equation:

$$Ax^2 + Bxy + Cy^2 + Dx + Ey + F = 0$$

From $(x_n, y_n)$ and $(x_{n+1}, y_{n+1})$ both being on this curve, we get:

$$Ax_n^2 + Bx_ny_n + Cy_n^2 + Dx_n + Ey_n + F = 0 \tag{1}$$ $$A(t^2 y_n^2 - 2 t x_n y_n + x_n^2) + B(t^3 y_n^2 - 2 t^2 x_n y_n + t x_n^2 - t y_n^2 + x_n y_n) + C(t^4 y_n^2 - 2 t^3 x_n y_n + t^2 x_n^2 - 2 t^2 y_n^2 + 2 t x_n y_n + y_n^2) + D(x_n - ty_n) + E(-t^2y_n + t x_n + y_n) + F = 0 \tag{2}$$

Subtracting (1) from (2) gives:

$$A(t^2 y_n^2 - 2 t x_n y_n) + B(t^3 y_n^2 - 2 t^2 x_n y_n + t x_n^2 - t y_n^2) + C(t^4 y_n^2 - 2 t^3 x_n y_n + t^2 x_n^2 - 2 t^2 y_n^2 + 2 t x_n y_n) + D(-ty_n) + E(-t^2y_n + tx_n) = 0$$

Note that we can now factor out a $t$ (assumed to be nonzero). Also, for brevity, let's drop the subscripts.

$$A(ty^2 - 2xy) + B(t^2y^2 - 2 t x y + x^2 - y^2) + C(t^3y^2 - 2 t^2 x y + t x^2 - 2 t y^2 + 2x y) + D(-y) + E(-ty + x) = 0$$

We can rearrange this to get:

$$(B + Ct)x^2 + (-2A - 2Bt - 2Ct^2 + 2C)xy + (At + Bt^2 - B + Ct^3 - 2Ct)y^2 + Ex + (-D - Et)y = 0$$

We can make the LHS everywhere zero by setting $D = E = 0$, $B = -At$, and $C = A$. So the equation of the curve becomes:

$$Ax^2 - Atxy + Ay^2 + F = 0$$

The discriminant of this conic section is $A^2(t^2 - 4)$. Since $A^2$ is nonnegative, the sign depends only on the $t^2 - 4$ term. Therefore,

• If $|t| < 2$, the curve is an ellipse. (And if $t = 0$, it's a circle.)
• If $|t| = 2$, the curve is a parabola.
• If $|t| > 2$, the curve is a hyperbola.