Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP

Loading…

ENH: optimize/minpack: determine leastsq and fsolve epsfcn by dtype (fixes Trac #920) #427

Merged
merged 1 commit into from

2 participants

@pv
Owner
pv commented

Fix the behavior of minpack optimizers with float32. If the objective function returns float32, then the default value for epsfcn is too small, as it corresponds to float64 machine epsilon.

This PR makes the numerical differentiation epsilon be determined by function return type.

Closes Trac #920

@pv pv ENH: optimize/minpack: determine leastsq and fsolve epsfcn by dtype (…
…fixes Trac #920)

In particular, if the objective function returns float32, use a large
enough epsilon for gradient estimation.
703dc23
@dlax dlax merged commit 38ec1a2 into scipy:master
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Commits on Feb 9, 2013
  1. @pv

    ENH: optimize/minpack: determine leastsq and fsolve epsfcn by dtype (…

    pv authored
    …fixes Trac #920)
    
    In particular, if the objective function returns float32, use a large
    enough epsilon for gradient estimation.
This page is out of date. Refresh to see the latest.
Showing with 40 additions and 9 deletions.
  1. +18 −8 scipy/optimize/minpack.py
  2. +22 −1 scipy/optimize/tests/test_minpack.py
View
26 scipy/optimize/minpack.py
@@ -5,7 +5,8 @@
from numpy import atleast_1d, dot, take, triu, shape, eye, \
transpose, zeros, product, greater, array, \
- all, where, isscalar, asarray, inf, abs
+ all, where, isscalar, asarray, inf, abs, \
+ finfo, inexact, issubdtype, dtype, sqrt
from .optimize import Result, _check_unknown_options
error = _minpack.error
@@ -27,12 +28,16 @@ def _check_func(checker, argname, thefunc, x0, args, numinputs, output_shape=Non
else:
msg += "."
raise TypeError(msg)
- return shape(res)
+ if issubdtype(res.dtype, inexact):
+ dt = res.dtype
+ else:
+ dt = dtype(float)
+ return shape(res), dt
def fsolve(func, x0, args=(), fprime=None, full_output=0,
col_deriv=0, xtol=1.49012e-8, maxfev=0, band=None,
- epsfcn=0.0, factor=100, diag=None):
+ epsfcn=None, factor=100, diag=None):
"""
Find the roots of a function.
@@ -137,7 +142,7 @@ def fsolve(func, x0, args=(), fprime=None, full_output=0,
return res['x']
def _root_hybr(func, x0, args=(), jac=None,
- col_deriv=0, xtol=1.49012e-08, maxfev=0, band=None, eps=0.0,
+ col_deriv=0, xtol=1.49012e-08, maxfev=0, band=None, eps=None,
factor=100, diag=None, full_output=0, **unknown_options):
"""
Find the roots of a multivariate function using MINPACK's hybrd and
@@ -182,7 +187,9 @@ def _root_hybr(func, x0, args=(), jac=None,
n = len(x0)
if type(args) != type(()):
args = (args,)
- _check_func('fsolve', 'func', func, x0, args, n, (n,))
+ shape, dtype = _check_func('fsolve', 'func', func, x0, args, n, (n,))
+ if epsfcn is None:
+ epsfcn = sqrt(finfo(dtype).eps)
Dfun = jac
if Dfun is None:
if band is None:
@@ -241,7 +248,7 @@ def _root_hybr(func, x0, args=(), jac=None,
def leastsq(func, x0, args=(), Dfun=None, full_output=0,
col_deriv=0, ftol=1.49012e-8, xtol=1.49012e-8,
- gtol=0.0, maxfev=0, epsfcn=0.0, factor=100, diag=None):
+ gtol=0.0, maxfev=0, epsfcn=None, factor=100, diag=None):
"""
Minimize the sum of squares of a set of equations.
@@ -347,11 +354,14 @@ def leastsq(func, x0, args=(), Dfun=None, full_output=0,
n = len(x0)
if type(args) != type(()):
args = (args,)
- m = _check_func('leastsq', 'func', func, x0, args, n)[0]
+ shape, dtype = _check_func('leastsq', 'func', func, x0, args, n)
+ m = shape[0]
if n > m:
raise TypeError('Improper input: N=%s must not exceed M=%s' % (n,m))
+ if epsfcn is None:
+ epsfcn = sqrt(finfo(dtype).eps)
if Dfun is None:
- if (maxfev == 0):
+ if maxfev == 0:
maxfev = 200*(n + 1)
retval = _minpack._lmdif(func, x0, args, full_output, ftol, xtol,
gtol, maxfev, epsfcn, factor, diag)
View
23 scipy/optimize/tests/test_minpack.py
@@ -4,7 +4,8 @@
from __future__ import division, print_function, absolute_import
from numpy.testing import assert_, assert_almost_equal, assert_array_equal, \
- assert_array_almost_equal, TestCase, run_module_suite, assert_raises
+ assert_array_almost_equal, TestCase, run_module_suite, assert_raises, \
+ assert_allclose
import numpy as np
from numpy import array, float64
@@ -135,6 +136,11 @@ def test_wrong_shape_fprime_function(self):
deriv_func = lambda x: dummy_func(x, (3,3))
assert_raises(TypeError, optimize.fsolve, func, x0=[0,1], fprime=deriv_func)
+ def test_float32(self):
+ func = lambda x: np.array([x[0] - 1000, x[1] - 10000], dtype=np.float32)**2
+ p = optimize.fsolve(func, np.array([1, 1], np.float32))
+ assert_allclose(func(p), [0, 0], atol=1e-3)
+
class TestRootHybr(TestCase):
def test_pressure_network_no_gradient(self):
"""root/hybr without gradient, equal pipes -> equal flows"""
@@ -242,6 +248,21 @@ def test_wrong_shape_Dfun_function(self):
deriv_func = lambda x: dummy_func(x, (3,3))
assert_raises(TypeError, optimize.leastsq, func, x0=[0,1], Dfun=deriv_func)
+ def test_float32(self):
+ # From Track ticket #920
+ def func(p,x,y):
+ q = p[0]*np.exp(-(x-p[1])**2/(2.0*p[2]**2))+p[3]
+ return q - y
+
+ x = np.array([ 1.475,1.429,1.409,1.419,1.455,1.519,1.472, 1.368,1.286,
+ 1.231], dtype=np.float32)
+ y = np.array([0.0168,0.0193,0.0211,0.0202,0.0171,0.0151,0.0185,0.0258,
+ 0.034,0.0396], dtype=np.float32)
+ p0 = np.array([1.0,1.0,1.0,1.0])
+ p1, success = optimize.leastsq(func, p0, args=(x,y))
+
+ assert_(success in [1,2,3,4])
+ assert_((func(p1,x,y)**2).sum() < 1e-4 * (func(p0,x,y)**2).sum())
class TestCurveFit(TestCase):
def setUp(self):
Something went wrong with that request. Please try again.