Interpreting multivariable root-finding results from Matlab's fsolve algorithm

Question

Edit: So I was able to get the same value of r that's given, when coding up the sum of squares of function values directly in the script file, rather than on the Command Window. So, maybe there's a round-off error from doing it manually. Thanks.

I'm learning multivariable root finding, starting with some simple sets of nonlinear equations in several variables -- solving for F = 0. I'm making up easy equations in which I know some obvious multivariable roots, and I am choosing starting points just a bit away from the roots, in order to get solutions from Matlab's fsolve algorithm. At the moment, I'm trying to learn how to interpret the algorithm's feedback, since my results differ from Matlab's, which is here:

Equation solved.

The sum of squared function values, r = 6.701089e-19, is less than sqrt(options.FunctionTolerance) = 1.000000e-03. The relative norm of the gradient of r, 7.127651e-09, is less than options.OptimalityTolerance = 1.000000e-06.

Optimization Metric

relative norm(grad r) = 7.13e-09
r = 6.70e-19

Options

OptimalityTolerance = 1e-06 (default)

sqrt(FunctionTolerance) = 1.0e-03 (default)

When I manually compute the sum of squared function values in the Command Window, I am getting something that's roughly 1e-7, which is nothing close to their r-value of 6.7e-19. What is this value r? Matlab calls it the "sum of squared function values", too. They also give a "relative norm" of grad r. I know of the usual Euclidean norm, p-norms, the infinity norm, and the Frobenius norm, but I've never heard of a relative norm. So, I suspect that Matlab's definitions are proprietary and not standard -- or, manual computations in the Command Window might have round-off errors.

Answer 1

The sum of squares should add up to what MATLAB says it adds up to, that is for $F:\mathbb{R}^m\to\mathbb{R}^n$ $$ r(x) = \sum\limits_{i=1}^{n}(F_i(x))^2$$

It is then easy to compute the gradient of $r$ as $$\nabla r(x) = J(x)^TF(x) $$, where $J$ is the Jacobian of $F$.

Let us now try to reverse engineer if this is indeed what MATLAB means. Since it is not specified in the what is $F$ exactly in the question, I tried fsolve on

$$F = \begin{pmatrix}e^{e^{-(x_1+x_2)}}-x_2*(1+x_1^2) \\x_1*\cos(x_2)+x_2*\sin(x_1)-\frac{1}{2}\end{pmatrix}$$,

which is an example found in the MATLAB documentation page for fsolve:

https://se.mathworks.com/help/optim/ug/fsolve.html

The following code snippet solves $F(x) = 0$ for this example and prints out the sum of squares and the norm of the gradient of r.

f = @(x) [exp(exp(-(x(1)+x(2))))-x(2)*(1+x(1)^2);x(1)*cos(x(2))+x(2)*sin(x(1))-0.5];

x_ = fsolve(f,[0 0]);

r = norm(f(x_))^2 %Sum of squares

%Estimate the gradient of r using finite differences
h = 1e-12;
J_t = (1/h)*[f(x_+[h,0])-f(x_),f(x_+[0,h])-f(x_)];%Transpose of Jacobian
norm_gradr = norm(J_t'*f(x_)) %Norm of gradient of r

Ideally, you would find the Jacobian by hand differentiation. This is a lazy estimation using finite differences, but it does the job. The output is:

r =

   6.8935e-23


norm_gradr =

   9.8783e-12

Equation solved. The sum of squared function values, r = 6.893461e-23, is less than
sqrt(options.FunctionTolerance) = 1.000000e-03. The relative norm of the gradient of r,
9.832419e-12, is less than options.OptimalityTolerance = 1.000000e-06.

Optimization Metric                                         Options
relative norm(grad r) =   9.83e-12              OptimalityTolerance =   1e-06 (default)
r =   6.89e-23                              sqrt(FunctionTolerance) =  1.0e-03 (default)

This short bit of reverse engineering confirms that if $x_*$ is the solution computed by MATLAB, then it means $r(x_*)$ when it says "sum of squares", and $\|\nabla r(x_*)\|_2 = \| J^T(x_*)F(x_*)\|$ " when it says "relative norm(grad r)".

You should check your computation of the sum of squares again.