Add MTG for Generalized Assignment Problem (#28)

* Add MTG

* Keep original fixed size Fortran arrays

* Make error codes private

.flake8

@@ -1,3 +1,5 @@
 [flake8]
 max-line-length = 79
-exclude = .svn,CVS,.bzr,.hg,.git,__pycache__,.eggs,*.egg,*/_vendor/*,node_modules
+exclude = .svn,CVS,.bzr,.hg,.git,__pycache__,.eggs,*.egg,*/_vendor/*,node_modules
+per-file-ignores =
+     tests/test__generalized_assignment.py: E127, E128, E131, E201

README.md

@@ -13,6 +13,7 @@ Solving knapsack problems with Python using algorithms by [Martello and Toth](ht
 * Multiple 0-1 knapsack problem: MTM, MTHM
 * Change-making problem: MTC2
 * Bounded change-making problem: MTCB
+* Generalized assignment problem: MTG
 
 Documentation is available [here](https://mknapsack.readthedocs.io).
 
@@ -28,7 +29,7 @@ Documentation is available [here](https://mknapsack.readthedocs.io).
      * `conda install -c conda-forge m2w64-toolchain_win-64`, or
      * [Install MSYS2](https://www.msys2.org) and `pacman -S --needed base-devel mingw-w64-x86_64-toolchain`
 
-2. `pip install -U mknapsack`
+2. `pip install mknapsack`
 
 
 ## Example usage
@@ -135,6 +136,27 @@ capacity = 190
 res = solve_bounded_change_making(weights, n_items, capacity)
 ```
 
+### Generalized Assignment Problem
+
+```python
+from mknapsack import solve_generalized_assignment
+
+# Given seven item types with the following knapsack dependent weights:
+weights = [[4, 1, 2, 1, 4, 3, 8],
+           [9, 9, 8, 1, 3, 8, 7]]
+
+# ...and the following knapsack dependent weights:
+profits = [[6, 9, 4, 2, 10, 3, 6],
+           [4, 8, 9, 1, 7, 5, 4]]
+
+# ...and two knapsack with the following capacities:
+capacities = [11, 22]
+
+# Assign items into the knapsacks while maximizing profits
+res = solve_generalized_assignment(weights, profits, capacities)
+```
+
+
 ## References
 
 * [Knapsack problems: algorithms and computer implementations](https://dl.acm.org/doi/book/10.5555/98124) by S. Martello and P. Toth, 1990

mknapsack/__init__.py

@@ -5,9 +5,13 @@
     'solve_bounded_knapsack',
     'solve_bounded_change_making',
     'solve_change_making',
+    'solve_generalized_assignment',
     'solve_multiple_knapsack',
     'solve_single_knapsack',
-    'solve_unbounded_knapsack'
+    'solve_unbounded_knapsack',
+    'FortranInputCheckError',
+    'NoSolutionError',
+    'ProblemSizeError'
 ]
 
 import os
@@ -21,9 +25,13 @@
     os.add_dll_directory(extra_dll_dir)
 
 
+from mknapsack._exceptions import FortranInputCheckError, NoSolutionError, \
+    ProblemSizeError  # noqa: E402
+
 from mknapsack._bounded import solve_bounded_knapsack  # noqa: E402
 from mknapsack._bounded_change_making import solve_bounded_change_making  # noqa: E402, E501
 from mknapsack._change_making import solve_change_making  # noqa: E402
+from mknapsack._generalized_assignment import solve_generalized_assignment  # noqa: E402, E501
 from mknapsack._multiple import solve_multiple_knapsack  # noqa: E402
 from mknapsack._single import solve_single_knapsack  # noqa: E402
 from mknapsack._unbounded import solve_unbounded_knapsack  # noqa: E402

mknapsack/_algos.f

@@ -7179,6 +7179,8 @@ subroutine mtg ( n, m, p, w, c, minmax, z, xstar, back, jck, jb )
 c y(l,j) =   packed solution of the relaxed problem for item  j  at
 c            level  l  of the branch-decision tree.
 c
+cf2py intent(in) n, m, p, w, c, minmax, back, jck
+cf2py intent(out) z, xstar, jb
       integer       p(10,100),w(10,100),c(10),xstar(100),z,back
       integer       h,s,r,u,su,vc,sb,t,qh,zbound,vjjub,penalt
       integer       dd(10),ud(10),q(10),pakl(10),ip(10),ir(10),il(10),

mknapsack/_exceptions.py

@@ -5,13 +5,17 @@ class NoSolutionError(Exception):
     """Error when no solution to a problem exists."""
 
 
+class ProblemSizeError(Exception):
+    """Error when problem size is too large to be solved."""
+
+
 class FortranError(Exception):
     """Error when running Fortran code."""
 
 
 class FortranInputCheckError(Exception):
     """Error in Fortran source code input validation."""
-    error_codes = {
+    _error_codes = {
         'mtm': {
             -1: 'Number of items/knapsacks is either too small or large',
             -2: 'Profit, weight or capacity is <= 0',
@@ -79,6 +83,15 @@ class FortranInputCheckError(Exception):
                 'greater than knapsack capacity',
             -5: 'Total weight of all items is smaller than or equal to '
                 'knapsack capacity'
+        },
+        'mtg': {
+            -1: 'Number of knapsacks is less than 2',
+            -2: 'Number of items is less than 2',
+            -3: 'Number of knapsacks cannot be larger than number of bits - 2',
+            -4: 'Profit, weight, or capacity is <= 0',
+            -5: 'One or more of weights is greater than knapsack capacity',
+            -6: 'One or more knapsacks cannot fit any items',
+            -7: 'Number of branching trees is too small for the problem size'
         }
     }
 
@@ -90,9 +103,9 @@ def __init__(self, method=None, z=None):
     def __str__(self):
         if self.method is None or self.z is None:
             return 'Generic exception from Fortran side'
-        elif (self.method not in self.error_codes or
-              self.z not in self.error_codes[self.method]):
+        elif (self.method not in self._error_codes or
+              self.z not in self._error_codes[self.method]):
             return ('Generic exception from Fortran side, could not resolve '
                     f'error code {self.z}')
         else:
-            return self.error_codes[self.method][self.z]
+            return self._error_codes[self.method][self.z]

mknapsack/_generalized_assignment.py

@@ -0,0 +1,172 @@
+"""Module for solving generalized assignment problem."""
+
+
+import logging
+
+from typing import List, Optional
+
+import numpy as np
+
+from mknapsack._algos import mtg
+from mknapsack._exceptions import FortranInputCheckError, NoSolutionError, \
+    ProblemSizeError
+from mknapsack._utils import preprocess_array, pad_array
+
+
+logger = logging.getLogger(__name__)
+
+
+def solve_generalized_assignment(
+    profits: List[List[int]],
+    weights: List[List[int]],
+    capacities: List[int],
+    maximize: bool = True,
+    method: str = 'mtg',
+    method_kwargs: Optional[dict] = None,
+    verbose: bool = False
+) -> np.ndarray:
+    """Solves the generalized assignment problem.
+
+    Given a set of items with knapsack-dependent weights and profits, assign
+    items exactly to one knapsack while maximizing or minimizing profit.
+
+    Args:
+        profits: Profit of item when assigned to a knapsack, where
+            knapsacks are in the rows and items in the columns.
+        weights: Weight of each item when assigned to a knapsack, where
+            knapsacks are in the rows and items in the columns.
+        capacities: Capacity of each knapsack, which should match the number of
+            rows in `profits` and `weights`.
+        maximize: Whether to maximize or minimize profits. Defaults to True.
+        method:
+            Algorithm to use for solving, should be one of
+
+                - 'mtg' - provides a fast heuristical solution that might not
+                  be the global optimum, but is suitable for larger problems,
+                  or an exact solution if required
+
+                  The problem size is limited as follows:
+                    * Number of knapsacks <= 10
+                    * Number of items <= 100
+
+            Defaults to 'mtg'.
+        method_kwargs:
+            Keyword arguments to pass to a given `method`.
+
+                - 'mtg'
+                    * **require_exact** (int, optional) - Whether to require an
+                      exact solution or not (0=no, 1=yes). Defaults to 0.
+                    * **max_backtracks** (int, optional) - The maximum number
+                      of backtracks to perform when ``require_exact=0``.
+                      Defaults to 100000.
+                    * **check_inputs** (int, optional) - Whether to check
+                      inputs or not (0=no, 1=yes). Defaults to 1.
+
+            Defaults to None.
+
+    Returns:
+        np.ndarray: Assigned knapsack for each item.
+
+    Raises:
+        NoSolutionError: No feasible solution found.
+        ProblemSizeError: When problem is too large to be solved with the given
+            method.
+        FortranInputCheckError: Something is wrong with the inputs when
+            validated in the original Fortran source code side.
+        ValueError: Something is wrong with the given inputs.
+
+    Example:
+        .. code-block:: python
+
+            from mknapsack import solve_generalized_assignment
+
+            res = solve_generalized_assignment(
+                weights=[[4, 1, 2, 1, 4, 3, 8],
+                         [9, 9, 8, 1, 3, 8, 7]],
+                profits=[[6, 9, 4, 2, 10, 3, 6],
+                         [4, 8, 9, 1, 7, 5, 4]],
+                capacities=[11, 22],
+                maximize=True
+            )
+
+    References:
+        * Silvano Martello, Paolo Toth, Knapsack Problems: Algorithms and
+          Computer Implementations, Wiley, 1990, ISBN: 0-471-92420-2,
+          LC: QA267.7.M37.
+
+        * `Original Fortran77 source code by Martello and Toth\
+          <https://people.sc.fsu.edu/~jburkardt/f77_src/knapsack/knapsack.f>`_
+    """
+    profits = preprocess_array(profits)
+    weights = preprocess_array(weights)
+    capacities = preprocess_array(capacities)
+
+    m_profits = profits.shape[0]
+    m_weights = weights.shape[0]
+    n_profits = profits.shape[1]
+    n_weights = weights.shape[1]
+    n = n_profits
+    m = len(capacities)
+
+    if m_profits != m or m_weights != m:
+        raise ValueError(
+            f'Number of rows in profits ({m_profits}) or weights '
+            f'({m_weights}) is not equal to the length of capacities')
+
+    if n_profits != n_weights:
+        raise ValueError('Profits length must be equal to weights '
+                         f'({n_profits != n_weights}')
+
+    method = method.lower()
+    method_kwargs = method_kwargs or {}
+    if method == 'mtg':
+        maxm = 10
+        maxn = 100
+
+        # The problem size is limited on Fortran side
+        error_msg = (
+            'The number of {0} {1} exceeds {2}, and there is no guarantee '
+            'that the algorithm will work as intended or provides a solution '
+            'at all')
+        if m > maxm:
+            raise ProblemSizeError(error_msg.format('knapsacks', m, maxm))
+        if n > maxn:
+            raise ProblemSizeError(error_msg.format('items', n, maxn))
+
+        p = pad_array(profits, (maxm, maxn))
+        w = pad_array(weights, (maxm, maxn))
+        c = pad_array(capacities, maxm)
+
+        if method_kwargs.get('require_exact', 0):
+            back = -1
+        else:
+            back = method_kwargs.get('max_backtracks', 100_000)
+
+        z, x, jb = mtg(
+            n=n,
+            m=m,
+            p=p,
+            w=w,
+            c=c,
+            minmax=2 if maximize else 1,
+            back=back,
+            jck=method_kwargs.get('check_inputs', 1)
+        )
+        print(z)
+        print(x)
+        print(jb)
+        if z == 0:
+            raise NoSolutionError('No feasible solution found')
+        elif z < 0:
+            raise FortranInputCheckError(method=method, z=z)
+
+        if verbose:
+            logger.info(f'Method: "{method}"')
+            logger.info(f'Total profit: {z}')
+            logger.info(f'Solution vector: {x}')
+            bound_name = 'Upper' if maximize else 'Lower'
+            logger.info(f'{bound_name} bound: {jb}')
+    else:
+        raise ValueError(f'Given method "{method}" not known')
+
+    return np.array(x)[:n]

mknapsack/_utils.py

@@ -10,12 +10,12 @@ def preprocess_array(ar, dtype='int32'):
 
 
 def pad_array(ar, width):
-    """Pad array with zeros to given length."""
-    assert ar.ndim == 1
-    new_ar = np.zeros((width, ), dtype=ar.dtype, order='F')
-    n = len(ar)
-    new_ar[:n] = ar
-    return new_ar
+    """Pad one or two-dimensional array with zeros."""
+    if ar.ndim == 1:
+        pad_width = (0, width - ar.shape[0])
+    elif ar.ndim == 2:
+        pad_width = ((0, width[0] - ar.shape[0]), (0, width[1] - ar.shape[1]))
+    return np.pad(ar, pad_width, mode='constant')
 
 
 def check_all_int(ar):