Simplifying Quadratic Forms via Matrices (Handling the Parabolic Case with Infinitesimal Perturbation)

(To follow this post comfortably you should have some background in linear algebra — at least a working understanding of its basic concepts.)

The semester is drawing to a close, and in our Analytic Geometry course we've reached the topic of simplifying quadratic curves. For students without a linear algebra background, however, directly applying the computational formulas for axis rotation and axis translation to carry out the transformation is enough to drive one crazy with the sheer volume of computation involved. And although the "invariant" method is somewhat simpler computationally, it always feels a bit mysterious — one never quite knows where it comes from, and on top of that you have to memorize a pile of formulas. In fact, if you have a linear algebra background, all of this becomes quite easy to understand. I like to pursue a unified method for solving the same class of problem — that is, treating all quadratic forms in a single, uniform way — and of course I'd also like the computation to be as light as possible.

The general model

A general quadratic form can be written as

$$x^T A x + 2 b^T x + c=0$$

where $x,b$ are both n-dimensional column vectors (with entries $x_i$ and $b_i$), A is an n-th order square matrix (with entries $a_{ij}$), and c is a constant. Here we only discuss the cases n = 2 and n = 3. The process of simplifying a quadratic form can be reduced to simplifying the matrix A. more

Most treatments of this problem extend the matrix A into an (n+1)-dimensional matrix $\widetilde{A}$, expressed in block form as:

$$\widetilde{A}=\left[\begin{array}{c c}A & b\\ b^T & c \end{array}\right]$$

while the vector $x$ is likewise extended into an (n+1)-dimensional vector $\widetilde{x}$:

$$\widetilde{x}=\left[\begin{array}{c}x\\1\end{array}\right]=\left[\begin{array}{c}x_1\\x_2\\...\\x_n\\1\end{array}\right]$$

so that the quadratic form becomes

$$\widetilde{x}^T \widetilde{A} \widetilde{x}=0$$

In this case the problem turns into simplifying the matrix $\widetilde{A}$. To my mind, though, while this approach is elegant and concise, it involves too much computation — for instance, a planar quadratic curve gives rise to a third-order square matrix, so that finding the eigenvalues turns what used to be a quadratic equation into a cubic one, and computing determinants gets considerably more troublesome with each increase in order. In this case, the benefit of conciseness feels like it's not worth the cost. So I tried simplifying an n-dimensional quadratic form directly using an n-th order matrix, and arrived at the following result.

My own idealized model

Consider an isometry matrix R. "Isometry" here means that under this transformation, the "distance" of a vector is preserved, i.e. $v^T v =(Rv)^T (Rv)=v^T (R^T R) v$, so we have $R^T R=I$ — the transpose of the matrix equals its inverse. Under such a transformation the shape of the curve is preserved. The most familiar example of an isometry is rotation about the origin.

If the simplified quadratic form takes the shape:

$$x^T M x =k$$

where M is a diagonal matrix and k is a constant, this is the simplest possible form of a quadratic form, and the one most convenient for our purposes.

Consider the composition of an isometry with a translation, i.e.

$$x \to R(x+x_0)$$

$x_0$ is also an n-dimensional constant column vector; substituting into the above expression gives:

$$\begin{eqnarray*} (x+x_0)^T R^T A R(x+x_0)=k \\ x^T(R^T M R)x+x_0^T[(R^T M R)+(R^T M R)^T]x+x_0^T (R^T M R) x_0-k=0 \end{eqnarray*}$$

Matching this against the original general quadratic form, we get

$$\begin{eqnarray*} R^T M R=A \\ \left[(R^T M R)+(R^T M R)^T \right]x_0 =2b \\ x_0^T (R^T M R) x_0 -k=c \end{eqnarray*}$$

where $R^T A R=R^{-1} A R=M$ is a similarity transformation of the matrix; in general its simplest form is the diagonal matrix consisting of all the eigenvalues of A. So it is relatively easy to determine its form. Next we have:

$$\begin{eqnarray*} (A+A^T)x_0=2b \\ x_0^T A x_0-k=c \end{eqnarray*}$$

The key is to find the explicit form of $x_0$, which amounts to solving a system of n linear equations in n unknowns, namely

$$x_0=2(A+A^T)^{-1}b$$

Regarding quadratic forms, there are some quite exciting results (see The Mathematical Bridge for proofs):

1. Every real quadratic form always corresponds to a real symmetric matrix (symmetric about the diagonal) A;
2. Every real symmetric matrix A always has n real eigenvalues.

With these two guarantees in hand, we no longer need to worry about cases where the "eigenvalues" turn out to be complex numbers, and moreover the computation above can be simplified further, because for a symmetric matrix we have $A^T=A$, so it becomes

$$x_0=A^{-1} b$$

A general treatment of the special case

At this point the problem seems entirely solved: first find the eigenvalues to get the matrix M; then compute $x_0=A^{-1} b$, and substitute into $x_0^T A x_0-k=c$ to find the constant term. Especially for quadratic curves, the amount of computation involved is really quite small. But what if the matrix A is not invertible? This can indeed happen — namely, when one of the eigenvalues is zero, which in the plane corresponds to the case of a parabola. Then computing an inverse matrix becomes impossible — does the method simply fail here?

For a parabola, we can actually still think of it as a limiting case of a hyperbola or an ellipse. For example, the parabola $2py+x^2=0$ can be viewed as $2py+x^2+ \varepsilon y^2$, where $\varepsilon \to 0$, and completing the square gives:

$$x^2+ \varepsilon (y+\frac{p}{\varepsilon})^2-\frac{p^2}{\varepsilon}=0$$

which is still a sum of squared terms.

This suggests an approach: we can take the originally non-invertible matrix A and add a small perturbation to it so that it becomes invertible, work out the result, and then take the limit as the perturbation tends to zero. This is a rather interesting idea, and we'll pursue it below.

The key question is where to add the perturbation, and how. In general, one could perturb any single element of the original matrix, but for the sake of convenience and uniformity in the treatment, we choose to perturb the diagonal elements. If there is exactly one zero eigenvalue, we add $\varepsilon_1$ to the first diagonal entry of A; if there are two zero eigenvalues, we further add $\varepsilon_2$ to the second diagonal entry, and so on. The purpose of doing this is to handle the higher-dimensional case uniformly; call the new matrix $\bar{A}$. If there are k zero eigenvalues, then the determinant $\det \bar{A} =(a_{11}+a_{22}+...+a_{nn})\varepsilon_1 \varepsilon_2 ... \varepsilon_k$ (that is, the sum of the diagonal entries times all the $\varepsilon_i$ — readers can easily verify this using facts about the "trace" of a matrix).

Now everything reverts to the original procedure:

The original n eigenvalues $\lambda_i$ become $\lambda_i -\varepsilon_i$, so that no zero eigenvalue appears anymore (with the zero eigenvalues placed first). Then, from $x_0=(\bar{A})^{-1} b$, we solve for the constant vector $x_0$ — of course this is an infinite vector, but that doesn't matter, since we only need to keep terms up to first order in the infinitesimal. We can then find the constant term: $k=x_0^T A x_0-c$. It should take the form:

$$\frac{C_1}{\varepsilon_1}+\frac{C_2}{\varepsilon_2}+...+Const$$

Comparing with

$$x^2+ \varepsilon (y+\frac{p}{\varepsilon})^2-\frac{p^2}{\varepsilon}=0$$

it's easy to see that $2\sqrt{C_i}x_i$ is exactly the term corresponding to the eigenvalue equal to zero — the parabolic term.

Worked examples:

  1. Elliptic (hyperbolic) type.

$$x^2+6xy+y^2+6x+2y-1=0$$

Solution:

The above expression is

$$\left[\begin{array}{c c}1 & 3\\ 3 & 1 \end{array}\right] \left[\begin{array}{c}x\\ y \end{array}\right]+2[3,1]\left[\begin{array}{c c}x\\ y \end{array}\right]-1=0$$

that is

$$A=\left[\begin{array}{c c}1 & 3\\ 3 & 1 \end{array}\right],b=\left[\begin{array}{c}3\\1 \end{array}\right],c=-1$$

Find the two eigenvalues of A: $(1-\lambda)^2-3^2=0$, namely -2 and 4.

Next compute

$$\begin{aligned}x_0=\left[\begin{array}{c c}1 & 3\\ 3 & 1 \end{array}\right]^{-1} \left[\begin{array}{c}3\\1 \end{array}\right]=\left[\begin{array}{c}0\\1 \end{array}\right] \\ k=x_0^T A x_0-c=x_0^T b=[0,1]\left[\begin{array}{c}3\\1 \end{array}\right]+1=2\end{aligned}$$

The final result is

$$\left[\begin{array}{c c}-2 & 0\\0 & 4 \end{array}\right]\left[\begin{array}{c}x\\ y \end{array}\right]=2$$

that is: $-x^2+2y^2=1$

  1. Parabolic type

$$x^2-4xy+4y^2+2x-2y-1=0$$

Solution:

The above expression is

$$\left[\begin{array}{c c}1 & -2\\ -2 & 4 \end{array}\right] \left[\begin{array}{c}x\\ y \end{array}\right]+2[1,-1]\left[\begin{array}{c c}x\\ y \end{array}\right]-1=0$$

that is

$$A=\left[\begin{array}{c c}1 & -2\\ -2 & 4 \end{array}\right],b=\left[\begin{array}{c}1\\-1 \end{array}\right],c=-1$$

The two eigenvalues are 0 and 5.

Perturb the matrix A into $\bar{A}$, i.e.

$$\bar{A}=\left[\begin{array}{c c}1+\varepsilon & -2\\ -2 & 4 \end{array}\right]$$

For a second-order square matrix, computing the inverse is fairly straightforward:

$$\left[\begin{array}{c c}a & c\\ b & d\end{array}\right]^{-1}=\frac{\left[\begin{array}{c c}d & -b\\ -c & a\end{array}\right]}{\det \left[\begin{array}{c c} a & c\\ b & d\end{array}\right]}$$

The approximate expansion of a second-order determinant is likewise fairly simple:

$$\det \left[\begin{array}{c c} a+\varepsilon & c\\ b & d \end{array} \right] =ad-bc+(a+d)\varepsilon$$

Hence

$$\bar{A}^{-1}=\frac{1}{5\varepsilon}\left[\begin{array}{c c}4 & 2\\ 2 & 1 \end{array}\right]$$

(accurate to first order in the infinitesimal)

$$x_0=\bar{A}^{-1}b=\frac{1}{5\varepsilon}\left[\begin{array}{c c}4 & 2\\ 2 & 1 \end{array}\right]\left[\begin{array}{c c}1\\-1 \end{array}\right]=\frac{1}{5\varepsilon}\left[\begin{array}{c}2\\1 \end{array}\right]$$

which gives:

$$k=x_0^T A x_0-c=x_0^T b=\frac{1}{5\varepsilon}[2,1]\left[\begin{array}{c}1\\-1 \end{array}\right]+1=\frac{1}{5\varepsilon}+1$$

This correspondingly produces the term $2\sqrt{\frac{1}{5}}y$, so the quadratic curve simplifies to

$$5x^2 \pm 2\sqrt{\frac{1}{5}}y=0$$

(the sign here doesn't matter — one just opens upward, the other downward.)

For three-dimensional surfaces, or even higher-dimensional cases, this method proceeds through exactly the same steps; naturally the amount of computation increases, but the essence of the method remains unchanged. As for the "invariant" method of simplification — I believe that the invariants themselves are, in fact, produced through the very sequence of steps described above. When transplanted into analytic geometry, it's really just the same procedure stripped of its linear-algebra clothing, with (not necessarily intentionally) a somewhat mystifying gloss added on top.

At the same time, I sincerely hope that in analytic geometry exams, problems should not be restricted to a fixed method — there shouldn't be rigid questions like "simplify using the invariant method..." That kind of prescriptive framing does nothing to improve our abilities. After all, memorizing the methods given in a textbook is not what matters; what matters fundamentally is mastering the mathematical skill needed to handle a given problem.

English translation of a post from 科学空间 | Scientific Spaces by 苏剑林. Original: https://kexue.fm/archives/1841
Translated automatically with claude-sonnet-5; all equations are reproduced verbatim from the source. Copyright remains with the original author.