eqSolve#

Purpose#

Solves a system of nonlinear equations.

Format#

{ x, retcode } = eqSolve(&F, start)#

Parameters:

&F (scalar) – a pointer to a procedure which computes the value at x of the equations to be solved and returns them in a column vector.
start (Kx1 vector) – starting values.

Returns:

x (Kx1 vector) – solution.

retcode (scalar) –

the return code:

1	Norm of the scaled function value is less than __Tol. x given is an approximate root of \(F(x)\) (unless __Tol is too large).
2	The scaled distance between the last two steps is less than the step-tolerance (_eqs_StepTol). x may be an approximate root of \(F(x)\), but it is also possible that the algorithm is making very slow progress and is not near a root, or the step-tolerance is too large.
3	The last global step failed to decrease `norm2(F(x))` sufficiently; either x is close to a root of \(F(x)\) and no more accuracy is possible, or an incorrectly coded analytic Jacobian is being used, or the secant approximation to the Jacobian is inaccurate, or the step-tolerance is too large.
4	Iteration limit exceeded.
5	Five consecutive steps of maximum step length have been taken; either `norm2(F(x))` asymptotes from above to a finite value in some direction or the maximum step length is too small.
6	x seems to be an approximate local minimizer of `norm2(F(x))` that is not a root of \(F(x)\). To find a root of \(F(x)\), restart `eqSolvefrom()` a different region.

Global Input#

The following are set by eqSolveSet:

_eqs_JacobianProc#: pointer to a procedure which computes the analytical Jacobian. By default, eqSolve() will compute the Jacobian numerically.

_eqs_MaxIters#: scalar, the maximum number of iterations. Default = 100.

_eqs_StepTol#: scalar, the step tolerance. Default = \(\_\_macheps^{\frac{2}{3}}\).

_eqs_TypicalF#: Kx1 vector of the typical \(F(x)\) values at a point not near a root, used for scaling. This becomes important when the magnitudes of the components of \(F(x)\) are expected to be very different. By default, function values are not scaled.

_eqs_TypicalX#: Kx1 vector of the typical magnitude of x, used for scaling. This becomes important when the magnitudes of the components of x are expected to be very different. By default, variable values are not scaled.

_eqs_IterInfo#: scalar, if nonzero, iteration information is printed. Default = 0.

The following are set by gausset():

__Tol#: scalar, the tolerance of the scalar function \(f = 0.5*\left|\left|F(x)\right|\right|^2\) required to terminate the algorithm. Default = 1e-5.

__altnam#: Kx1 character vector of alternate names to be used by the printed output. By default, the names X1, X2, X3... or X01, X02, X03... (depending on how __vpad is set) will be used.

__output#: scalar. If non-zero, final results are printed.

__title#: string, a custom title to be printed at the top of the iterations report. By default, only a generic title will be printed.

__vpad#

scalar. If __altnam is not set, variable names are automatically created. Two types of names can be created:

0	Variable names are not padded to give them equal length. For example, `X1, X2, ..., X10, ...`
1	Variable names are padded with zeros to give them an equal number of characters. For example, `X01, X02, ..., X10, ...` This is useful if you want the variable names to sort properly.

Examples#

// Reset all eqSolve global variables
eqSolveSet();

// Procedure to compute nonlinear equations
proc (1) = f(x);
   local f1, f2, f3;

   // Function one
   f1 = 3*x[1]^3 + 2*x[2]^2 + 5*x[3] - 10;

   // Function two
   f2 = -x[1]^3 - 3*x[2]^2 + x[3] + 5;

   // Function three
   f3 = 3*x[1]^3 + 2*x[2]^2 - 4*x[3];

   retp(f1|f2|f3);
endp;

// Procedure to compute Jacobian
proc (1) = fjc(x);
   local fjc1, fjc2, fjc3;

   fjc1 = 9*x[1]^2 ~ 4*x[2] ~ 5;
   fjc2 = -3*x[1]^2 ~ -6*x[2] ~ 1;
   fjc3 = 9*x[1]^2 ~ 4*x[2] ~ -4;

   retp(fjc1|fjc2|fjc3);
endp;

// Starting values
start = { -1, 12, -1 };

// Set pointer to Jacobian procedure
_eqs_JacobianProc = &fjc;

// Perform estimation and print report
{ x, tcode } = eqSolve(&f, start);

=========================================================
 EqSolve Version 14.0.1
=========================================================

||F(X)|| at final solution:                   0.93699762
---------------------------------------------------------
Termination Code = 1:

Norm of the scaled function value is less than __Tol;
---------------------------------------------------------

---------------------------------------------------------
VARIABLE     START          ROOTS               F(ROOTS)
---------------------------------------------------------
X1          -1.00000       0.54144351      4.4175402e-006
X2          12.00000        1.4085912     -6.6263102e-006
X3          -1.00000        1.1111111      4.4175402e-006
---------------------------------------------------------

Remarks#

The equation procedure should return a column vector containing the result for each equation. For example, consider a two-equation system given by:

\[x_1^2 + x_2^2 - 2 = 0\]

\[e^{x_1-1} + x_2^3 - 2 = 0\]

proc (1) = f(var);
   local x1, x2, eqns;

   // Set x1
   x1 = var[1];

   // Set x2
   x2 = var[2];

   // Pre-allocate output vector
   eqns = zeros(2,1);

   // Equation 1
   eqns[1] = x1^2 + x2^2 - 2;

   // Equation 2
   eqns[2] = exp(x1-1) + x2^3 - 2;

   retp(eqns);
endp;

Note that the first equation in the system is contained in the first row of eqns and the second equation is in the second row of eqns.

Source#

eqsolve.src