snlls: Fix violation of boundaries due to float-point round-off err…

…ors (#188)

* Jacobian: add rounding to avoid float-point round-off errors

* Jacobian: use exception handling to avoid rounding in all cases

* snlls: shield linear LSQ results against float-point errors

deerlab/snlls.py

@@ -293,6 +293,9 @@ def snlls(y, Amodel, par0, lb=None, ub=None, lbl=None, ubl=None, nnlsSolver='cvx
         elif nnlsSolver == 'cvx':
             linSolver = lambda AtA, Aty: cvxnnls(AtA, Aty, tol=lin_tol, maxiter=lin_maxiter)
         parseResult = lambda result: result
+
+    # Ensure correct formatting and shield against float-point errors
+    validateResult = lambda result: np.maximum(lbl,np.minimum(ubl,np.atleast_1d(result)))
     # ----------------------------------------------------------
 
     # Containers for alpha-update checks
@@ -319,9 +322,9 @@ def linear_problem(A,optimize_alpha,alpha):
 
         # Solve the linear least-squares problem
         result = linSolver(AtA, Aty)
-        linfit = parseResult(result)
-        linfit = np.atleast_1d(linfit)
-        
+        result = parseResult(result)
+        linfit = validateResult(result)
+
         return linfit, alpha
     #===========================================================================
 
@@ -566,4 +569,3 @@ def _plot(subsets,y,yfit,show):
         plt.close()
     return fig
 # ===========================================================================================
-

deerlab/utils/utils.py

@@ -409,7 +409,7 @@ def Jacobian(fcn, x0, lb, ub):
 
     """
     J = opt._numdiff.approx_derivative(fcn,x0,method='2-point',bounds=(lb,ub))
-    J = np.atleast_2d(J)
+    J = np.atleast_2d(J) 
     return J
 #===============================================================================