Merge 7f24600 into 649cc55

docs/qe_apidoc.py

@@ -80,6 +80,15 @@
     :show-inheritance:
 """
 
+optimize_module_template = """{mod_name}
+{equals}
+
+.. automodule:: quantecon.optimize.{mod_name}
+    :members:
+    :undoc-members:
+    :show-inheritance:
+"""
+
 random_module_template = """{mod_name}
 {equals}
 
@@ -211,6 +220,12 @@ def model_tool():
     # Alphabetize
     markov.sort()
 
+    # list file names with optimize
+    optimize_files = glob("../quantecon/optimize/[a-z0-9]*.py")
+    optimize = list(map(lambda x: x.split('/')[-1][:-3], optimize_files))
+    # Alphabetize
+    optimize.sort()
+
     # list file names with random
     random_files = glob("../quantecon/random/[a-z0-9]*.py")
     random = list(map(lambda x: x.split('/')[-1][:-3], random_files))
@@ -230,7 +245,7 @@ def model_tool():
     # Alphabetize
     util.sort()
 
-    for folder in ["game_theory", "markov", "random", "tools", "util"]:
+    for folder in ["game_theory", "markov", "optimize", "random", "tools", "util"]:
         if not os.path.exists(source_join(folder)):
             os.makedirs(source_join(folder))
 
@@ -257,6 +272,13 @@ def model_tool():
             equals = "=" * len(mod)
             f.write(markov_module_template.format(mod_name=mod, equals=equals))
 
+    # Write file for each optimize file
+    for mod in optimize:
+        new_path = os.path.join("source", "optimize", mod + ".rst")
+        with open(new_path, "w") as f:
+            equals = "=" * len(mod)
+            f.write(optimize_module_template.format(mod_name=mod, equals=equals))
+
     # Write file for each random file
     for mod in random:
         new_path = os.path.join("source", "random", mod + ".rst")
@@ -284,24 +306,28 @@ def model_tool():
 
     gt = "game_theory/" + "\n   game_theory/".join(game_theory)
     mark = "markov/" + "\n   markov/".join(markov)
+    opti = "optimize/" + "\n   optimize/".join(optimize)
     rand = "random/" + "\n   random/".join(random)
     tlz = "tools/" + "\n   tools/".join(tools)
     utls = "util/" + "\n   util/".join(util)
     #-TocTree-#
     toc_tree_list = {"game_theory": gt,
                      "markov": mark,
+                     "optimize" : opti,
                      "tools": tlz,
                      "random": rand,
                      "util": utls,
                      }
 
-    for f_name in ("game_theory", "markov", "random", "tools", "util"):
+    for f_name in ("game_theory", "markov", "optimize", "random", "tools", "util"):
         with open(source_join(f_name + ".rst"), "w") as f:
             m_name = f_name
             if f_name == "game_theory":
                 f_name = "Game Theory"                                             #Produce Nicer Title for Game Theory Module
             if f_name == "util":
                 f_name = "Utilities"            #Produce Nicer Title for Utilities Module
+            if f_name == "optimize":
+                f_name = "Optimize"
             temp = split_file_template.format(name=f_name.capitalize(),
                                               equals="="*len(f_name),
                                               files=toc_tree_list[m_name])

docs/source/optimize.rst

@@ -0,0 +1,7 @@
+Optimize
+========
+
+.. toctree::
+   :maxdepth: 2
+
+   optimize/scalar_maximization

docs/source/optimize/scalar_maximization.rst

@@ -0,0 +1,7 @@
+scalar_maximization
+===================
+
+.. automodule:: quantecon.optimize.scalar_maximization
+    :members:
+    :undoc-members:
+    :show-inheritance:

quantecon/__init__.py

@@ -12,6 +12,7 @@
 from . import game_theory
 from . import quad
 from . import random
+from . import optimize
 
 #-Objects-#
 from .compute_fp import compute_fixed_point

quantecon/optimize/__init__.py

@@ -0,0 +1,6 @@
+"""
+Initialization of the optimize subpackage
+"""
+
+from .scalar_maximization import brent_max
+

quantecon/optimize/scalar_maximization.py

@@ -0,0 +1,144 @@
+import numpy as np
+from numba import jit, njit
+
+@njit
+def brent_max(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
+    -------
+    xf : float
+        The maximizer
+    fval : float
+        The maximum value attained
+    info : tuple
+        A tuple of the form (status_flag, num_iter).  Here status_flag
+        indicates whether or not the maximum number of function calls was
+        attained.  A value of 0 implies that the maximum was not hit.  
+        The value `num_iter` is the number of function calls.
+
+    Example
+    -------
+
+    ```
+        @njit
+        def f(x):
+            return -(x + 2.0)**2 + 1.0
+
+        xf, fval, info = maximize_scalar(f, -2, 2)
+    ```
+
+    """
+
+    maxfun = maxiter
+    status_flag = 0
+
+    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
+
+    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
+
+        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:
+            status_flag = 1
+            break
+
+    fval = -fx
+    info = status_flag, num
+
+    return xf, fval, info

quantecon/optimize/tests/test_scalar_max.py

@@ -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 brent_max
+
+@njit
+def f(x):
+    """
+    A function for testing on.
+    """
+    return -(x + 2.0)**2 + 1.0
+
+def test_brent_max():
+    """
+    Uses the function f defined above to test the scalar maximization 
+    routine.
+    """
+    true_fval = 1.0
+    true_xf = -2.0
+    xf, fval, info = brent_max(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_brent_max():
+    """
+    Uses the function f defined above to test the scalar maximization 
+    routine.
+    """
+    y = 5
+    true_fval = 5.0
+    true_xf = -0.0
+    xf, fval, info = brent_max(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__)
+
+

quantecon/version.py

@@ -1,4 +1,4 @@
 """
 This is a VERSION file and should NOT be manually altered
 """
-version = '0.3.7'
+version = '0.3.8'

setup.py

@@ -98,6 +98,7 @@ def write_version_py(filename=None):
                 'quantecon.game_theory',
                 'quantecon.game_theory.game_generators',
                 'quantecon.markov',
+                'quantecon.optimize',
                 'quantecon.random',
                 'quantecon.tests',
                 'quantecon.util',
@@ -113,7 +114,7 @@ def write_version_py(filename=None):
       download_url='https://github.com/QuantEcon/QuantEcon.py/tarball/' + VERSION,
       keywords=['quantitative', 'economics'],
       install_requires=[
-          'numba>=0.36.2',
+          'numba>=0.38',
           'numpy',
           'scipy>=1.0.0',
           'sympy',