search

MITgcm / gcmfaces

gcmfaces is a Matlab / Octave toolbox that handles gridded earth variables in generic fashion. Read more at:
http://gcmfaces.readthedocs.io/en/latest/
MIT License
24 stars 21 forks source link

Replace mldivide with lsqminnorm in diffsmooth2D_div_inv.m #15

Closed owang01 closed 6 months ago

owang01 commented 2 years ago

When solving A xx = yy, diffsmooth2D_div_inv.m uses xx=A\yy. This becomes less robust when matrix A is close to singular. The solution xx may become all zeros.

Different versions of MATLAB may have different behaviors. For instance, when using calc_barostream.m (which calls diffsmooth2D_div_inv.m) to calculate barotropic stream function, gcmfaces using matlab/2017b was able to find a solution that appears fine, while matlab/2021a gives a solution of all zeros.

This Pull Request uses the more robust MATLAB function lsqminnorm to find a solution in a least-squares sense. This method can also handle sparse matrix and is more efficient than pinv (which cannot handle sparse matrix). Warning is turned on to monitor if the matrix is close to singular.

gaelforget commented 2 years ago

I'll double check with latest Matlab -- I had been using 2018b. Do you know when lsqminnorm was introduced? We may need a check on matlab version being used for backward compatibility.

owang01 commented 2 years ago

That is a good point. lsqminnorm was introduced in R2017b.