Clean parameterized hfun interface

docs/pounders.rst

@@ -88,16 +88,28 @@ Python
 ^^^^^^^^
 .. mat:autofunction:: pounders.m.pounders
 
-:math:`\hfun` Functions
-^^^^^^^^^^^^^^^^^^^^^^^
-The following :math:`\hfun` functions are available for use with both the Python
-and |matlab| implementations of |pounders|.  While they are presented through
-their integration into the Python package, the documentation is valid for the
-|matlab| version of these functions, which are located in
+General :math:`\hfun` Functions
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+The following :math:`\hfun` functions are available for immediate use with both
+the Python and |matlab| implementations of |pounders|.  While they are presented
+through their integration into the Python package, the documentation is valid
+for the |matlab| version of these functions, which are located in
 ``pounders/m/general_h_funs``.
 
 .. autofunction:: ibcdfo.pounders.h_leastsquares
 .. autofunction:: ibcdfo.pounders.h_neg_leastsquares
 .. autofunction:: ibcdfo.pounders.h_identity
 .. autofunction:: ibcdfo.pounders.h_emittance
-.. autofunction:: ibcdfo.pounders.h_squared_diff_from_mean
+
+Parameterized :math:`\hfun` Functions
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+|pounders| specifies parameterized :math:`\hfun` functions, which it can only
+use once users have chosen a single set of parameter values for formulating a
+specific :math:`\hfun` function and, therefore, a single related problem.  The
+following routines can be used to create a single ``hfun`` and ``combinemodels``
+matched pair for a single set of desired parameter values.  While these routines
+are presented through their integration into the Python package, the
+documentation is valid for the |matlab| version of these routines, which are
+located in ``pounders/m/general_h_funs``.
+
+.. autofunction:: ibcdfo.pounders.create_squared_diff_from_mean_functions

pounders/m/general_h_funs/Testcreatefunctions.m

@@ -0,0 +1,75 @@
+% From MATLAB and the directory containing this file, execute
+%     >> runtests
+%
+% If you would like to run this from a different folder that includes these
+% tests as a subfolder, then from that folder execute
+%     >> runtests("IncludeSubfolders", true)
+%
+% To execute the test suite with coverage enabled and to generate an HTML-format
+% coverage report, execute from /path/to/IBCDFO/pounders/m
+%     >> runtests("IncludeSubfolders", true, "ReportCoverageFor", pwd)
+%
+
+% TODO: If this is a good idea, then test the combinemodel function in the
+% same way.
+classdef Testcreatefunctions < matlab.unittest.TestCase
+    methods (Test)
+
+        function badArguments(testCase)
+            testCase.verifyError( ...
+                @()create_squared_diff_from_mean_functions([]), ...
+                'MATLAB:validators:mustBeNonempty' ...
+            );
+
+            testCase.verifyError( ...
+                @()create_squared_diff_from_mean_functions([1.1 2.2]), ...
+                'MATLAB:validators:mustBeScalarOrEmpty' ...
+            );
+
+            testCase.verifyError( ...
+                @()create_squared_diff_from_mean_functions("hello"), ...
+                'MATLAB:validators:mustBeReal' ...
+            );
+
+            for bad = [inf -inf nan]
+                testCase.verifyError( ...
+                    @()create_squared_diff_from_mean_functions(bad), ...
+                    'MATLAB:validators:mustBeFinite' ...
+                );
+                testCase.verifyError( ...
+                    @()create_squared_diff_from_mean_functions([bad]), ...
+                    'MATLAB:validators:mustBeFinite' ...
+                );
+            end
+        end
+
+        function actualArguments(testCase)
+            F = [1.1 2.2 3.3 4.4];
+            expected_0 = 6.05;
+            expected_2 = -9.075;
+
+            alpha = 0.0;
+            [hfun_0, cm_0] = create_squared_diff_from_mean_functions(alpha);
+            h_0 = hfun_0(F);
+            testCase.assertTrue(abs(h_0 - expected_0) <= 5 * eps);
+            testCase.assertEqual(h_0, hfun_0(F));
+
+            % Confirm that we get same result as for actual scalar.
+            [hfun_tmp, cm_tmp] = create_squared_diff_from_mean_functions([alpha]);
+            testCase.assertEqual(h_0, hfun_tmp(F));
+
+            % Intentionally alter variable used to create *_0 functions
+            alpha = 2.0;
+            [hfun_2, cm_2] = create_squared_diff_from_mean_functions(alpha);
+            h_2 = hfun_2(F);
+            testCase.assertEqual(h_2, expected_2);
+            testCase.assertEqual(h_2, hfun_2(F));
+            testCase.assertNotEqual(h_0, h_2);
+            % Confirm that changes to actual alpha argument used to construct
+            % hfun_0 and cm_0 do not alter those functions.  This check is
+            % motivated by technical subtleties seen with Python.
+            testCase.assertEqual(h_0, hfun_0(F));
+        end
+
+    end
+end

pounders/m/general_h_funs/create_squared_diff_from_mean_functions.m

@@ -0,0 +1,38 @@
+function [hfun, combinemodels] = create_squared_diff_from_mean_functions(alpha)
+    % Please refer to the documentation for the Python version of this
+    % h function.
+    arguments
+        alpha {mustBeScalarOrEmpty, mustBeNonempty, mustBeReal, mustBeFinite}
+    end
+
+    function [h] = h_squared_diff_from_mean(F, alpha)
+        h = sum((F - 1 / length(F) * sum(F)).^2) - alpha * (1 / length(F) * sum(F))^2;
+    end
+
+    function [G, H] = combine_squared_diff_from_mean(Cres, Gres, Hres, alpha)
+        [n, ~, m] = size(Hres);
+
+        m_sumF = mean(Cres);
+        m_sumG = 1 / m * sum(Gres, 2);
+        sumH = sum(Hres, 3);
+
+        G = zeros(n, 1);
+        for i = 1:m
+            G = G + (Cres(i) - m_sumF) * (Gres(:, i) - m_sumG);
+        end
+        G = 2 * G - 2 * alpha * m_sumF * m_sumG;
+
+        H = zeros(n, n);
+        for i = 1:m
+            H = H + (Cres(i) - m_sumF) * (Hres(:, :, i) + sumH) + (Gres(:, i) - m_sumG) * (Gres(:, i) - m_sumG)';
+        end
+        H = 2 * H;
+
+        H = H - (2 * alpha / m) * m_sumF * sumH - (2 * alpha) * m_sumG * m_sumG';
+
+        % [grad, Hess] = matlab_symbolic_grad(Cres,Gres,Hres);
+    end
+
+    hfun = @(F) h_squared_diff_from_mean(F, alpha);
+    combinemodels = @(Cres, Gres, Hres) combine_squared_diff_from_mean(Cres, Gres, Hres, alpha);
+end

pounders/m/tests/benchmark_pounders.m

@@ -56,8 +56,7 @@
             combinemodels = @combine_leastsquares;
         elseif hfun_cases == 2
             ALPHA = 0;
-            hfun = @(F) h_squared_diff_from_mean(F, ALPHA);
-            combinemodels = @(Cres, Gres, Hres) combine_squared_diff_from_mean(Cres, Gres, Hres, ALPHA);
+            [hfun, combinemodels] = create_squared_diff_from_mean_functions(ALPHA);
         elseif hfun_cases == 3
             if m ~= 3 % Emittance is only defined for the case when m == 3
                 continue

pounders/m/tests/test_bounds_and_sp1.m

@@ -65,8 +65,7 @@
                 combinemodels = @combine_leastsquares;
             elseif hfun_cases == 2
                 ALPHA = 0;
-                hfun = @(F) h_squared_diff_from_mean(F, ALPHA);
-                combinemodels = @(Cres, Gres, Hres) combine_squared_diff_from_mean(Cres, Gres, Hres, ALPHA);
+                [hfun, combinemodels] = create_squared_diff_from_mean_functions(ALPHA);
             elseif hfun_cases == 3
                 hfun = @h_neg_leastsquares;
                 combinemodels = @combine_neg_leastsquares;

pounders/py/__init__.py

@@ -12,10 +12,10 @@
     combine_identity, h_identity,
     combine_leastsquares, h_leastsquares,
     combine_neg_leastsquares, h_neg_leastsquares,
-    combine_emittance, h_emittance,
-    combine_squared_diff_from_mean, h_squared_diff_from_mean
+    combine_emittance, h_emittance
 )
 # fmt: on
+from .create_squared_diff_from_mean_functions import create_squared_diff_from_mean_functions
 
 # -- Python unittest-based test framework
 # Used for automatic test discovery by main package

pounders/py/create_squared_diff_from_mean_functions.py

@@ -0,0 +1,67 @@
+import numbers
+
+import numpy as np
+
+
+def create_squared_diff_from_mean_functions(alpha):
+    r"""
+    Create an :math:`\hfun` function and its associated combinemodel function,
+    both using the given parameter value :math:`\alpha`, for constructing and
+    using the |pounders| objective function
+
+    .. math::
+
+        f(\psp; \alpha) = \hfun\left(\Ffun(\psp); \alpha\right)
+          = \sum_{i=1}^{\nd} \left(\Ffuncomp{i}(\psp) - \overline{\Ffun}(\psp)\right)^2
+            - \alpha \overline{\Ffun}(\psp)^2
+
+    where
+
+    .. math::
+
+        \overline{\Ffun}(\psp) = \frac{1}{\nd}\sum_{i=1}^{\nd} \Ffuncomp{i}(\psp)
+
+    is the average value of all components in :math:`\Ffun(\psp)`.   This
+    objective, therefore, prefers vectors close to their average.
+
+    :param: :math:`\alpha` is a problem-specific parameter that specifies how
+        much the objective function should penalize small (large) averages for
+        :math:`\alpha` positive (negative).
+    :return: (hfun, combinemodels) constructed with the given :math:`\alpha`
+    """
+    # NOTE: Don't alter alpha anywhere in this function.
+
+    if not isinstance(alpha, numbers.Real):
+        raise TypeError("alpha parameter must be floating point number")
+
+    # Both of these definitions assume that alpha is a variable in the local
+    # scope of this function.  In this case, it's the function's argument.
+    def hfun(F):
+        F_avg = np.mean(F)
+        return np.sum((F - F_avg) ** 2) - alpha * F_avg**2
+
+    def combinemodels(Cres, Gres, Hres):
+        n, _, m = Hres.shape
+
+        m_sumF = np.mean(Cres)
+        m_sumG = 1 / m * np.sum(Gres, axis=1)
+        sumH = np.sum(Hres, axis=2)
+
+        G = np.zeros(n)
+        for i in range(m):
+            G = G + (Cres[i] - m_sumF) * (Gres[:, i] - m_sumG)
+        G = 2 * G - 2 * alpha * m_sumF * m_sumG
+
+        H = np.zeros((n, n))
+        for i in range(m):
+            H = H + (Cres[i] - m_sumF) * (Hres[:, :, i] + sumH) + np.outer(Gres[:, i] - m_sumG, Gres[:, i] - m_sumG)
+
+        H = 2 * H
+
+        H = H - (2 * alpha / m) * m_sumF * sumH - (2 * alpha) * np.outer(m_sumG, m_sumG)
+
+        return G, H
+
+    # The value of the given actual alpha argument is fixed into the returned
+    # functions only at this point.
+    return hfun, combinemodels

pounders/py/general_h_funs.py

@@ -108,64 +108,3 @@ def h_emittance(F):
     """
     assert len(F) == 3, "Emittance must have exactly 3 inputs"
     return F[0] * F[1] - F[2] ** 2
-
-
-def h_squared_diff_from_mean(F, alpha):
-    r"""
-    :math:`\hfun` function for constructing the |pounders| objective function
-
-    .. math::
-
-        f(\psp; \alpha) = \hfun\left(\Ffun(\psp); \alpha\right)
-          = \sum_{i=1}^{\nd} \left(\Ffuncomp{i}(\psp) - \overline{\Ffun}(\psp)\right)^2
-            - \alpha \overline{\Ffun}(\psp)^2
-
-    where
-
-    .. math::
-
-        \overline{\Ffun}(\psp) = \frac{1}{\nd}\sum_{i=1}^{\nd} \Ffuncomp{i}(\psp)
-
-    is the average value of all components in :math:`\Ffun(\psp)`.   This
-    objective, therefore, prefers vectors close to their average.
-
-    :param: :math:`\alpha` is a problem-specific parameter that specifies how
-        much the objective function should penalize small (large) averages for
-        :math:`\alpha` positive (negative).
-
-    Users can create |pounders|-compatible versions of this function and its
-    combine models function for a particular :math:`\alpha` value with code such
-    as
-
-    .. code:: python
-
-        import functools
-        ALPHA = X.Y
-        hfun = functools.partial(h_squared_diff_from_mean, alpha=ALPHA)
-        combinemodels = functools.partial(combine_squared_diff_from_mean, alpha=ALPHA)
-    """
-    F_avg = np.mean(F)
-    return np.sum((F - F_avg) ** 2) - alpha * F_avg**2
-
-
-def combine_squared_diff_from_mean(Cres, Gres, Hres, alpha):
-    n, _, m = Hres.shape
-
-    m_sumF = np.mean(Cres)
-    m_sumG = 1 / m * np.sum(Gres, axis=1)
-    sumH = np.sum(Hres, axis=2)
-
-    G = np.zeros(n)
-    for i in range(m):
-        G = G + (Cres[i] - m_sumF) * (Gres[:, i] - m_sumG)
-    G = 2 * G - 2 * alpha * m_sumF * m_sumG
-
-    H = np.zeros((n, n))
-    for i in range(m):
-        H = H + (Cres[i] - m_sumF) * (Hres[:, :, i] + sumH) + np.outer(Gres[:, i] - m_sumG, Gres[:, i] - m_sumG)
-
-    H = 2 * H
-
-    H = H - (2 * alpha / m) * m_sumF * sumH - (2 * alpha) * np.outer(m_sumG, m_sumG)
-
-    return G, H

pounders/py/tests/TestPoundersExtensive.py

@@ -6,7 +6,6 @@
 import unittest
 
 import ibcdfo
-import functools
 import numpy as np
 import scipy as sp
 from calfun import calfun
@@ -76,8 +75,7 @@ def Ffun_batch(Y):
                     hfun_name = combinemodels.__name__
                 elif hfun_cases == 2:
                     ALPHA = 0.0
-                    hfun = functools.partial(ibcdfo.pounders.h_squared_diff_from_mean, alpha=ALPHA)
-                    combinemodels = functools.partial(ibcdfo.pounders.combine_squared_diff_from_mean, alpha=ALPHA)
+                    hfun, combinemodels = ibcdfo.pounders.create_squared_diff_from_mean_functions(ALPHA)
                     hfun_name = "combine_squared_diff_from_mean"
                 elif hfun_cases == 3:
                     if m != 3:  # Emittance is only defined for the case when m == 3