More general cleanup (#95)

* Docstrings
* Whitespace (fewer blank lines in functions)
* Removed some cruft comments

gpfit/classes.py

@@ -7,8 +7,8 @@ def max_affine(x, ba):
     Evaluates max affine function at values of x, given a set of
     max affine fit parameters.
 
-    INPUTS
-    ------
+    Arguments
+    ---------
         x: 2D array [nPoints x nDim]
             Independent variable data
 
@@ -18,7 +18,7 @@ def max_affine(x, ba):
              [ ....,          ]
              [bk, ak1, ... akk]]
 
-    OUTPUTS
+    Returns
     -------
         y: 1D array [nPoints]
             Max affine output
@@ -28,12 +28,9 @@ def max_affine(x, ba):
     npt, dimx = x.shape
     K = ba.size // (dimx + 1)
     ba = np.reshape(ba, (dimx + 1, K), order="F")  # 'F' gives Fortran indexing
-
-    # augment data with column of ones
-    X = np.hstack((np.ones((npt, 1)), x))
+    X = np.hstack((np.ones((npt, 1)), x))  # augment data with column of ones
     y, partition = np.dot(X, ba).max(1), np.dot(X, ba).argmax(1)
 
-    # The not-sparse sparse version
     dydba = np.zeros((npt, (dimx + 1) * K))
     for k in range(K):
         inds = np.equal(partition, k)
@@ -50,7 +47,8 @@ def softmax_affine(x, params):
     Evaluates softmax affine function at values of x, given a set of
     SMA fit parameters.
 
-    INPUTS:
+    Arguments:
+    ----------
             x:      Independent variable data
                         2D numpy array [nPoints x nDimensions]
 
@@ -59,7 +57,8 @@ def softmax_affine(x, params):
                         [b1, a11, .. a1d, b2, a21, .. a2d, ...
                          bK, aK1, aK2, .. aKd, alpha]
 
-    OUTPUTS:
+    Returns:
+    --------
             y:      SMA approximation to log transformed data
                         1D numpy array [nPoints]
 
@@ -74,14 +73,11 @@ def softmax_affine(x, params):
         return np.inf * np.ones((npt, 1)), np.nan
     K = np.size(ba) // (dimx + 1)
     ba = ba.reshape(dimx + 1, K, order="F")
-
     X = np.hstack((np.ones((npt, 1)), x))  # augment data with column of ones
     z = np.dot(X, ba)  # compute affine functions
-
     y, dydz, dydsoftness = lse_scaled(z, alpha)
 
     dydsoftness = -dydsoftness * (alpha ** 2)
-
     nrow, ncol = dydz.shape
     repmat = np.tile(dydz, (dimx + 1, 1)).reshape(nrow, ncol * (dimx + 1), order="F")
     dydba = repmat * np.tile(X, (1, K))
@@ -97,7 +93,8 @@ def implicit_softmax_affine(x, params):
     Evaluates implicit softmax affine function at values of x, given a set of
     ISMA fit parameters.
 
-    INPUTS:
+    Arguments:
+    ----------
             x:      Independent variable data
                         2D numpy array [nPoints x nDimensions]
 
@@ -106,7 +103,8 @@ def implicit_softmax_affine(x, params):
                         [b1, a11, .. a1d, b2, a21, .. a2d, ...
                          bK, aK1, aK2, .. aKd, alpha1, alpha2, ... alphaK]
 
-    OUTPUTS:
+    Returns:
+    --------
             y:      ISMA approximation to log transformed data
                         1D numpy array [nPoints]
 
@@ -121,10 +119,8 @@ def implicit_softmax_affine(x, params):
     if any(alpha <= 0):
         return np.inf * np.ones((npt, 1)), np.nan
     ba = ba.reshape(dimx + 1, K, order="F")  # reshape ba to matrix
-
     X = np.hstack((np.ones((npt, 1)), x))  # augment data with column of ones
     z = np.dot(X, ba)  # compute affine functions
-
     y, dydz, dydalpha = lse_implicit(z, alpha)
 
     nrow, ncol = dydz.shape

gpfit/fit.py

@@ -20,7 +20,8 @@ def fit(xdata, ydata, K, ftype="ISMA", alpha0=10):
     """
     Fit a log-convex function to multivariate data, returning a FitConstraintSet
 
-    INPUTS
+    Arguments
+    ---------
         xdata:      Independent variable data
                         2D numpy array [nDim, nPoints]
                         [[<--------- x1 ------------->]
@@ -30,25 +31,27 @@ def fit(xdata, ydata, K, ftype="ISMA", alpha0=10):
                         1D numpy array [nPoints,]
                         [<---------- y ------------->]
 
+        K:          Number of terms
+
+        ftype:      Function class ["MA", "SMA", "ISMA"]
+
+        alpha0:     Initial guess for smoothing parameter alpha
+
+    Returns
+    -------
+        FitConstraintSet
+
         rms_error:  float
             Root mean square error between original (not log transformed)
             data and function fit.
     """
 
-    # Check data is in correct form
     if ydata.ndim > 1:
         raise ValueError("Dependent data should be a 1D numpy array")
 
-    # Transform to column vector
     xdata = xdata.reshape(xdata.size, 1) if xdata.ndim == 1 else xdata.T
-
-    # Dimension of function (number of independent variables)
-    d = int(xdata.shape[1])
-
-    # begin fit data dict to generate fit constraint set
+    d = int(xdata.shape[1])  # Number of dimensions
     fitdata = {"ftype": ftype, "K": K, "d": d}
-
-    # find bounds of x
     if d == 1:
         fitdata["lb0"] = exp(min(xdata.reshape(xdata.size)))
         fitdata["ub0"] = exp(max(xdata.reshape(xdata.size)))
@@ -107,22 +110,10 @@ def residual(params):
         raise ValueError("Fitted constraint contains too large a value...")
 
     def evaluate(xdata):
-        """
-        Evaluate the y of this fit over a range of xdata.
-
-        INPUTS
-            xdata:      Independent variable data
-                            2D numpy array [nDim, nPoints]
-                            [[<--------- x1 ------------->]
-                             [<--------- x2 ------------->]]
-
-        OUTPUT
-            ydata:      Dependent variable data in 1D numpy array
-        """
+        """Evaluate the fit over a range of xdata."""
         xdata = xdata.reshape(xdata.size, 1) if xdata.ndim == 1 else xdata.T
         return CLASSES[ftype](xdata, params)[0]
 
-    # cstrt.evaluate = evaluate
     fitdata["rms_err"] = sqrt(mean(square(evaluate(xdata.T) - ydata)))
     fitdata["max_err"] = sqrt(max(square(evaluate(xdata.T) - ydata)))
 

gpfit/initialize.py

@@ -11,29 +11,28 @@ def get_initial_parameters(x, y, K):
     ensures that initialization has at least K+1 points per partition (i.e.
     per affine function)
 
-    INPUTS:
+    Arguments:
+    ----------
         x:      Independent variable data
                     2D column vector [nPoints x nDims]
 
         y:      Dependent variable data
                     2D column vector [nPoints x 1]
 
-    OUTPUTS:
+    Returns:
+    --------
         ba:     Initial b and a parameters
                     2D array [(dimx+1) x k]
 
     """
 
     npt, dimx = x.shape
-
-    X = hstack((ones((npt, 1)), x))
-    b = zeros((dimx + 1, K))
-
     if K * (dimx + 1) > npt:
         raise ValueError("Not enough data points")
 
-    # Choose K unique indices
-    randinds = randperm(npt)[0:K]
+    X = hstack((ones((npt, 1)), x))
+    b = zeros((dimx + 1, K))
+    randinds = randperm(npt)[0:K]  # Choose K unique indices
 
     # partition based on distances
     sqdists = zeros((npt, K))

gpfit/least_squares.py

@@ -21,8 +21,8 @@ def levenberg_marquardt(
     Levenberg-Marquardt alogrithm
     Minimizes sum of squared error of residual function
 
-    INPUTS
-    ------
+    Arguments
+    ---------
     residfun: function
         Mapping from parameters to residuals, of the form
         (r, drdp) = residfun(params)
@@ -44,7 +44,7 @@ def levenberg_marquardt(
     tolrms: float
         Tolerance on change in rms error per iteration
 
-    OUTPUTS
+    Returns
     -------
     params: np.array (1D)
         Parameter vector that locally minimizes norm(residfun, 2)

gpfit/logsumexp.py

@@ -11,14 +11,16 @@ def lse_implicit(x, alpha):
     - lse_implicit is a mapping R^n --> R, where n number of dimensions
     - implementation: newton raphson steps to find f(x,y) = 0
 
-    INPUTS:
+    Arguments:
+    ----------
             x:      independent variable data
                         2D numpy array [nPoints x nDim]
 
             alpha:  local softness parameter
                         1D array [K] (K=number of terms)
 
-    OUTPUTS:
+    Returns:
+    --------
             y:      ISMA approximation to log transformed data
                         1D numpy array [nPoints]
 
@@ -71,12 +73,9 @@ def lse_implicit(x, alpha):
             / sumexpo[i]
         )
         neval = neval + 1
-
-        # update inds that need to be evaluated
-        i[i] = abs(f[i]) > tol
+        i[i] = abs(f[i]) > tol  # update inds that need to be evaluated
 
     y = m + L
-
     dydx = alphaexpo / (tile(sumalphaexpo, (nx, 1))).T
     dydalpha = (h - Lmat) * expo / (tile(sumalphaexpo, (nx, 1))).T
 
@@ -89,14 +88,16 @@ def lse_scaled(x, alpha):
     - sums across the second dimension of x
     - note that lse_scaled is a mapping R^n --> R
 
-    INPUTS:
+    Arguments:
+    ----------
             x:      independent variable data
                         2D numpy array [nPoints x nDim]
 
             alpha:  local softness parameter
                         1D array [K] (K=number of terms)
 
-    OUTPUTS:
+    Returns:
+    --------
             y:      ISMA approximation to log transformed data
                         1D numpy array [nPoints]
 

gpfit/print_fit.py

@@ -4,64 +4,56 @@
 
 # pylint: disable=invalid-name
 def print_isma(A, B, alpha, d, K):
-    "prints ISMA fit from params"
+    """prints ISMA fit from params"""
+
     print("ISMA fit from params")
     string_list = [None] * K
-
     print_string = "1 = "
     for k in range(K):
         if k > 0:
             print(print_string)
             print_string = "    + "
-
         print_string += "({0:.6g}/w**{1:.6g})".format(exp(alpha[k] * B[k]), alpha[k])
-
         for i in range(d):
             print_string += " * (u_{0:d})**{1:.6g}".format(
                 i + 1, alpha[k] * A[d * k + i]
             )
-
         string_list[k] = print_string
-
     print(print_string)
+
     return string_list
 
 
 # pylint: disable=invalid-name
 def print_sma(A, B, alpha, d, K):
-    "prints SMA fit from params"
+    """prints SMA fit from params"""
+
     print("SMA fit from params")
     string_list = [None] * K
-
     print_string = "w**{0:.6g} = ".format(alpha)
     for k in range(K):
         if k > 0:
             print(print_string)
             print_string = "    + "
-
         print_string += "{0:.6g}".format(exp(alpha * B[k]))
-
         for i in range(d):
             print_string += " * (u_{0:d})**{1:.6g}".format(i + 1, alpha * A[d * k + i])
-
         string_list[k] = print_string
-
     print(print_string)
+
     return string_list
 
 
 # pylint: disable=invalid-name
 def print_ma(A, B, d, K):
-    "prints MA fit from params"
+    """prints MA fit from params"""
+
     print("MA fit from params")
     string_list = [None] * K
-
     for k in range(K):
         print_string = "w = {0:.6g}".format(exp(B[k]))
-
         for i in range(d):
             print_string += " * (u_{0:d})**{1:.6g}".format(i + 1, A[d * k + i])
-
         string_list[k] = print_string
         print(print_string)
 

gpfit/xfoil/wrapper.py

@@ -9,9 +9,10 @@
 
 def xfoil_comparison(airfoil, Cl, Re, Cd):
     """
-    comparison of XFOIL cd to input cd for given cl and re
-    INPUTS
-    ------
+    Comparison of XFOIL Cd to input Cd for given Cl and Re.
+
+    Arguments
+    ---------
     airfoil:                    airfoil name
                                     str - ("xxx.dat", "naca xxxx")
     Cl:                         lift coefficient
@@ -21,7 +22,7 @@ def xfoil_comparison(airfoil, Cl, Re, Cd):
     Cd:                         drag coefficient
                                     float or 1d-list or array
 
-    OUTPUTS
+    Returns
     -------
     err:                        error between Cd and XFOIL computed Cd
                                     1-d numpy array
@@ -61,9 +62,10 @@ def xfoil_comparison(airfoil, Cl, Re, Cd):
 
 def single_call(topline, cl, Re, M, max_iter=100, pathname="/usr/local/bin/xfoil"):
     """
-    single XFOIL call for given cl, re and mach number
-    INPUTS
-    ------
+    Single XFOIL call for given Cl, Re and Mach number.
+
+    Arguments
+    ---------
     topline:                        load airfoil call in XFOIL
                                         str - e.g. "load xxx.dat"
     cl:                             lift coefficient
@@ -76,7 +78,7 @@ def single_call(topline, cl, Re, M, max_iter=100, pathname="/usr/local/bin/xfoil
                                         int: default is 100
     pathname:                       system path to XFOIL
                                         str
-    OUTPUTS
+    Returns
     -------
     cd:                         XFOIL drag coefficient
                                     float