7.3.2 Nonisolated stationary points

Because the equation is equivalent to a system of equations in unknowns, we expect that in most cases the solution set of this system will consist of a few discrete, isolated points in . However, it is possible for the solution set to contain curves, surfaces, or hypersurfaces as well as points. For example, suppose ; then the solution set of consists of the two lines and .

In this section, we will examine a marching method helpful in tracing curves of critical points. This type of degeneracy can occur if, for example, we are trying to find the minimum of the squared distance between two surfaces which happen to intersect. The method we use to trace out such curves involves setting up a system of differential equations which can be solved numerically using a standard ordinary differential equation system solver (see Sect. 5.8.1).

We begin our discussion by defining the concept of a degenerate critical point [323]:

Definition 7.3.3. A critical point of a function (that is, a point at which ) is called degenerate if the Hessian matrix of evaluated at is singular.

The following well-known theorem of differential geometry relates the concepts of degeneracy and isolation.

Theorem 7.3.5. Suppose a critical point of is nondegenerate, i.e. its Hessian is nonsingular. Then there exists a neighborhood of such that for all , . (A point having this property is called an isolated critical point.)

Proofs of this theorem are given in [323,129]. Unfortunately, the converse is not true; all nondegenerate critical points are isolated, but some isolated critical points are degenerate. For example, the function has an isolated critical point at the origin, but its Hessian matrix there is the zero matrix and therefore singular. Nevertheless, this theorem is practically helpful in eliminating most cases.

We will use this theorem to help us detect curve branches of nonisolated stationary points as follows:

1. Run the IPP algorithm described in Chap. 4 on our initial box of search, with a fairly coarse level of accuracy. For example, if we start with the search box , we may run our root-finding algorithm with a tolerance of or .

2. Check the remaining bounding boxes after this step. If we observe a number of boxes adjacent to one another, there may be a curve in the solution set.

3. Use Newton-Raphson iteration to find roots within these boxes to a high accuracy. Check the Hessian at these points; if it is singular, it is very likely that a curve exists. Accordingly, choose one of these roots as a starting point for our tracing technique.

This approach is by no means infallible; it is theoretically possible for there to be a number of degenerate isolated roots that will fool us into thinking that there is a curve involved. It is also possible that the solution set is a surface or even a higher dimensional entity. However, in practice these exceptions rarely occur.

Let us now assume that we have found a point on this curve, which we will call . We will assume that is a function of the parameter , which takes on a value 0 at our starting point and that has unit speed everywhere (i.e. everywhere). This unit speed condition is imposed to ensure that has an arc length parameterization.

Rewrite the equation as a system of equations in unknowns:

Now, in order to stay on the curve in moving a small distance from to , we need to know what increments have to be added to . Accordingly, we take a Taylor expansion of each equation around the critical point to first order:

Or, rewriting this as a matrix equation, we have

where and is the Hessian of evaluated at .

Because is degenerate, the rank of should be less than . In fact, we anticipate that it will be , since any satisfying (7.28) must be a vector tangent to at , and hence the nullspace of should have dimension . (It is possible on rare occasions that there is a curve passing through but the rank of the Hessian is less than . This problem may occur if the first-order Taylor expansion in (7.27) contains insufficient information due to zero derivatives. Fortunately, these cases are extremely rare in practice; see [323] for more information.)

Because we need to find a tangent vector to at with unit length, we need to find some vector in the nullspace of and then normalize it to length . The Singular Value Decomposition method [410] may be applied to to generate such a vector in a stable manner. Then, after making the vector unit length, we pass it to our ordinary differential equation solver as the tangent to at . We can trace backwards by changing the sign of the tangent vector.

Table 7.1: Curve and surface description (adapted from [461])

Case	Property	Degree(s)	Property	Degree(s)
P/C			rational	9
P/S			rational	7 and 7
C/C	integral	2	rational	2
C/S	rational	2	rational	4 and 3
S/S	rational	1 and 2	rational	1 and 2