**Homotopy method using arc length**Consider a set of

nonlinear algebraic equations that depend on a parameter

. In vector notation, these equations are represented as

, (5)

where the components of the vectors

and

are

and

.

Suppose we have found a solution to equation (5), say by Newton's method for

. We represent the solution as

, so that

.

We assume that the solution

is an analytic function of

, that is,

is continuous and differentiable. Then, given a solution

, we can find a neighboring solution

on the solution curve by constructing a Taylor series expansion about

. If the solution curve is a multivalued function of

, it is convenient to parameterize the solution in terms of a suitable parameter along the solution curve. The preferred parameter in this study is the arc length

along the solution curve so that we can write equation (5) as

.

The total derivative of

along the solution curve is then:

, (6)

which is a set of

equations in

unknowns, where the unknowns are

.

Recall that

is the familiar Jacobian used in Newton's method. It is convenient to write equation (6) in matrix notation

, (7)

where

and

. Thus the vector

is tangent to the solution curve

.

The homogeneous system of equations given by equation (7) is an underdetermined system of equations and will have a nontrivial solution if the rank

of

has full row rank, namely

. The nontrivial solution is not unique, however. In order to obtain a unique solution, we need to add an additional equation to equation (7). The appropriate equation is given by the definition of arc length:

. (8)

Equations (7) and (8) represent

equations in terms of the

unknowns. Although equation (7) is linear in the unknowns, equation (8) is not. Moreover, the coefficient matrix in equation (7) depends on

and

, usually in a nonlinear way, which means that we need to solve equations (7) and (8) using an appropriate ODE solver. Matters are complicated further by the fact that

is singular at bifurcation and turning points. However, the augmented matrix

has full rank at turning points. This feature allows us to devise a numerical algorithm that can integrate through turning points.

Generally one is not really interested in how

and

depend on

, and thus equation (8) is not really needed. The homogeneous system of equations defined by equation (7) can be solved using Cramer's rule. Thus if we define

, then the vector

obeys the system of equations:

, (9)

where

is the matrix with the

column removed. In the homotopy literature, Garcia and Zangwill [3] call equation (9) the basic differential equation. The advantage of using it instead of the full system of equations (7) and (8) is that equation (9) is linear in the derivatives. Note that the independent variable

is not the arc length along the solution curve, but is related to arc length by

.

Hence once the solution

has been found, a simple integration gives the dependence of

on

. The right-hand side of equation (9) can be determined readily when

is not too large using

*Mathematica *as shown in the recent paper by Higgins and Binous [2].

[1] A. Uppal, W. H. Ray, and A. B. Poore, "On the Dynamic Behavior of Continuous Strirred Tank Reactors,"

* Chemical Engineering Science,* **29**(4), 1974 pp. 967–985.

[2] B. G. Higgins and H. Binous, "A Simple Method for Tracking Turning Points in Parameter Space,"

*Journal of Chemical Engineering of Japan*,

**43**(12), 2010 pp. 1035–1042.

[3] C. B. Garcia and W. I. Zangwill,

*Pathways to Solutions, Fixed Points, and Equilibria*, Englewood Cliffs, NJ: Prentice–Hall, 1981.