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

The sum of squares should add up to what MATLAB says it adds up to, that is for $F:mathbb{R}^mtomathbb{R}^n$ $$ r(x) = sumlimits_{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.