New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
adding jitted scalar maximization routine, first build #416
Merged
Merged
Changes from 4 commits
Commits
Show all changes
12 commits
Select commit
Hold shift + click to select a range
152c067
adding jitted scalar maximization routine, first build
jstac 4384579
add args arguments
natashawatkins c9933d3
fix typo
natashawatkins 22dce1d
update install requires in setup.py
mmcky e890f71
added status_flag to scalar maximizer
jstac f722b25
modified return value to include number of function calls
jstac 70be953
modified return value to accommodate feedback and changed function na…
jstac 44df5df
add optimize subpackage to docs
mmcky 7526895
add optimize to top level index and setup for installation
mmcky 7f24600
add tracked generated sources
mmcky d221eb1
fix index generation
mmcky 09b5531
add tracked generated source files
mmcky File filter
Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,6 @@ | ||
""" | ||
Initialization of the optimize subpackage | ||
""" | ||
|
||
from .scalar_maximization import maximize_scalar | ||
|
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,138 @@ | ||
import numpy as np | ||
from numba import jit, njit | ||
|
||
@njit | ||
def maximize_scalar(func, a, b, args=(), xtol=1e-5, maxiter=500): | ||
""" | ||
Uses a jitted version of the maximization routine from SciPy's fminbound. | ||
The algorithm is identical except that it's been switched to maximization | ||
rather than minimization, and the tests for convergence have been stripped | ||
out to allow for jit compilation. | ||
|
||
Note that the input function `func` must be jitted or the call will fail. | ||
|
||
Parameters | ||
---------- | ||
func : jitted function | ||
a : scalar | ||
Lower bound for search | ||
b : scalar | ||
Upper bound for search | ||
args : tuple, optional | ||
Extra arguments passed to the objective function. | ||
maxiter : int, optional | ||
Maximum number of iterations to perform. | ||
xtol : float, optional | ||
Absolute error in solution `xopt` acceptable for convergence. | ||
|
||
Returns | ||
------- | ||
fval : float | ||
The maximum value attained | ||
xf : float | ||
The maximizer | ||
|
||
Example | ||
------- | ||
|
||
``` | ||
@njit | ||
def f(x): | ||
return -(x + 2.0)**2 + 1.0 | ||
|
||
fval, xf = maximize_scalar(f, -2, 2) | ||
``` | ||
|
||
""" | ||
|
||
maxfun = maxiter | ||
|
||
sqrt_eps = np.sqrt(2.2e-16) | ||
golden_mean = 0.5 * (3.0 - np.sqrt(5.0)) | ||
|
||
fulc = a + golden_mean * (b - a) | ||
nfc, xf = fulc, fulc | ||
rat = e = 0.0 | ||
x = xf | ||
fx = -func(x, *args) | ||
num = 1 | ||
fmin_data = (1, xf, fx) | ||
|
||
ffulc = fnfc = fx | ||
xm = 0.5 * (a + b) | ||
tol1 = sqrt_eps * np.abs(xf) + xtol / 3.0 | ||
tol2 = 2.0 * tol1 | ||
|
||
while (np.abs(xf - xm) > (tol2 - 0.5 * (b - a))): | ||
golden = 1 | ||
# Check for parabolic fit | ||
if np.abs(e) > tol1: | ||
golden = 0 | ||
r = (xf - nfc) * (fx - ffulc) | ||
q = (xf - fulc) * (fx - fnfc) | ||
p = (xf - fulc) * q - (xf - nfc) * r | ||
q = 2.0 * (q - r) | ||
if q > 0.0: | ||
p = -p | ||
q = np.abs(q) | ||
r = e | ||
e = rat | ||
|
||
# Check for acceptability of parabola | ||
if ((np.abs(p) < np.abs(0.5*q*r)) and (p > q*(a - xf)) and | ||
(p < q * (b - xf))): | ||
rat = (p + 0.0) / q | ||
x = xf + rat | ||
|
||
if ((x - a) < tol2) or ((b - x) < tol2): | ||
si = np.sign(xm - xf) + ((xm - xf) == 0) | ||
rat = tol1 * si | ||
else: # do a golden section step | ||
golden = 1 | ||
|
||
if golden: # Do a golden-section step | ||
if xf >= xm: | ||
e = a - xf | ||
else: | ||
e = b - xf | ||
rat = golden_mean*e | ||
|
||
if rat == 0: | ||
si = np.sign(rat) + 1 | ||
else: | ||
si = np.sign(rat) | ||
|
||
x = xf + si * np.maximum(np.abs(rat), tol1) | ||
fu = -func(x, *args) | ||
num += 1 | ||
fmin_data = (num, x, fu) | ||
|
||
if fu <= fx: | ||
if x >= xf: | ||
a = xf | ||
else: | ||
b = xf | ||
fulc, ffulc = nfc, fnfc | ||
nfc, fnfc = xf, fx | ||
xf, fx = x, fu | ||
else: | ||
if x < xf: | ||
a = x | ||
else: | ||
b = x | ||
if (fu <= fnfc) or (nfc == xf): | ||
fulc, ffulc = nfc, fnfc | ||
nfc, fnfc = x, fu | ||
elif (fu <= ffulc) or (fulc == xf) or (fulc == nfc): | ||
fulc, ffulc = x, fu | ||
|
||
xm = 0.5 * (a + b) | ||
tol1 = sqrt_eps * np.abs(xf) + xtol / 3.0 | ||
tol2 = 2.0 * tol1 | ||
|
||
if num >= maxfun: | ||
break | ||
|
||
fval = -fx | ||
|
||
return fval, xf | ||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,58 @@ | ||
""" | ||
Tests for scalar maximization. | ||
|
||
""" | ||
import numpy as np | ||
from numpy.testing import assert_almost_equal | ||
from numba import njit | ||
|
||
from quantecon.optimize import maximize_scalar | ||
|
||
@njit | ||
def f(x): | ||
""" | ||
A function for testing on. | ||
""" | ||
return -(x + 2.0)**2 + 1.0 | ||
|
||
def test_maximize_scalar(): | ||
""" | ||
Uses the function f defined above to test the scalar maximization | ||
routine. | ||
""" | ||
true_fval = 1.0 | ||
true_xf = -2.0 | ||
fval, xf = maximize_scalar(f, -2, 2) | ||
assert_almost_equal(true_fval, fval, decimal=4) | ||
assert_almost_equal(true_xf, xf, decimal=4) | ||
|
||
@njit | ||
def g(x, y): | ||
""" | ||
A multivariate function for testing on. | ||
""" | ||
return -x**2 + y | ||
|
||
def test_maximize_scalar_multivariate(): | ||
""" | ||
Uses the function f defined above to test the scalar maximization | ||
routine. | ||
""" | ||
y = 5 | ||
true_fval = 5.0 | ||
true_xf = -0.0 | ||
fval, xf = maximize_scalar(g, -10, 10, args=(y,)) | ||
assert_almost_equal(true_fval, fval, decimal=4) | ||
assert_almost_equal(true_xf, xf, decimal=4) | ||
|
||
|
||
if __name__ == '__main__': | ||
import sys | ||
import nose | ||
|
||
argv = sys.argv[:] | ||
argv.append('--verbose') | ||
argv.append('--nocapture') | ||
nose.main(argv=argv, defaultTest=__file__) | ||
|
||
|
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -1,4 +1,4 @@ | ||
""" | ||
This is a VERSION file and should NOT be manually altered | ||
""" | ||
version = '0.3.7' | ||
version = '0.3.8' | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. this will be set in |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think a status flag (
ierr
inscipy.fminbound
) should also be returned, which inscipy.fminbound
is0
if converged and1
ifmaxiter
is reached.