auction.pyx

#!python
#cython: boundscheck=False
#cython: wraparound=False
#cython: cdivision=True

cimport cython
import sys
from collections import namedtuple

import numpy as np
from cython.parallel import parallel, prange
from libcpp cimport bool
cimport numpy as np
from libc.stdlib cimport malloc, free

cdef float minimum = sys.float_info.min

cdef struct SearchResult:
    float price
    int item
    float payoff

cdef SearchResult new_search_results() nogil:
    cdef SearchResult sr
    sr.price = minimum
    sr.item = -1
    sr.payoff = 0.
    return sr

# TODO rename, linear_sum_assignment
# - package -> auction_matching
cdef void search_bid(float [:, :] X,
               int i,
               int start,
               int end,
               SearchResult best,
               SearchResult second,
               float *_prices) nogil:

    cdef SearchResult _best = new_search_results()
    cdef SearchResult _second = new_search_results()
    cdef int j
    cdef float price, net_payoff

    for j in range(start, end):
        price = X[i,j]
        net_payoff = price - _prices[i]
        if net_payoff > best.price:
            _second = _best
            _best.price = net_payoff
            _best.item = i
            _best.payoff = price
        elif net_payoff > second.price:
            _second.price = net_payoff
            _second.item = i
            _second.payoff = price

    #with search_lck:
    if _best.price > best.price:
        if best.price > _second.price:
            second = best
        else:
            second = _second
            best = _best
    elif _best.price > second.price:
        second = _best

cdef void bid(float [:, :] X,
              int n,
              int bidder,
              int n_blocks,
              float eps,
              int *_unassigned,
              int *_p_unassigned_loc,
              float *_payoffs,
              float *_prices,
              int *_assignments,
              int *_belongs_to) nogil:

    cdef int partition_size = int((n-1) / n_blocks + 1)

    cdef SearchResult best = new_search_results()
    cdef SearchResult second = new_search_results()

    cdef int start = 0
    cdef int end = partition_size

    cdef float new_bid = 0.
    cdef int worker

    for worker in prange(n_blocks):
        start = worker * partition_size
        end = start + partition_size
        if end > n:
            end = n

        search_bid(
            X,
            bidder,
            start,
            end,
            best,
            second,
            _prices,
        )

    #with _lck:
    with gil:
        new_bid = _prices[best.item] \
                    + (best.payoff - _prices[best.item]) \
                    - (second.payoff - _prices[second.item]) \
                    + eps
        if new_bid > _prices[best.item]:
            previous_owner = _belongs_to[best.item]
            if previous_owner != -1:
                _assignments[previous_owner] = -1
                _p_unassigned_loc[0] = _p_unassigned_loc[0] + 1
                _unassigned[_p_unassigned_loc[0]] = previous_owner
            _assignments[bidder] = best.item
            _belongs_to[best.item] = bidder
            _payoffs[bidder] = best.payoff
            _prices[best.item] = new_bid
        else:
            _p_unassigned_loc[0] = _p_unassigned_loc[0] + 1
            _unassigned[_p_unassigned_loc[0]] = bidder

def match(np.ndarray X, int n_bidders=1, int _n_blocks=1):
    cdef int m, n
    cdef bool reverse = False

    m = X.shape[0]
    n = X.shape[1]

    if m > n:
        # TODO work with transpose
        # note if width of matrix is smaller, then the algo doesn't work, we need to transpose, and
        # record the transposition so that we can reverse the results
        reverse = True

    cdef double eps = 1 / (m+1)

    cdef int *_unassigned = <int *>malloc(m * sizeof(int))
    cdef int _unassigned_loc = m
    cdef int *_p_unassigned_loc = &_unassigned_loc

    cdef int *_assignments = <int *>malloc(m * sizeof(int))
    cdef int *_belongs_to = <int *>malloc(n * sizeof(int))

    cdef float *_payoffs = <float *>malloc(m * sizeof(float))
    cdef float *_prices = <float *>malloc(n * sizeof(float))

    cdef int n_iter = 0
    cdef int bidder, i
    cdef int n_blocks = _n_blocks

    cdef np.ndarray[int, ndim=1] col_ind
    cdef np.ndarray[int, ndim=1] row_ind

    cdef float [:, :] X_view = X

    for i in range(m):
        _unassigned[i] = 0
        _assignments[i] = -1
        _payoffs[i] = 0.

    for i in range(n):
        _belongs_to[i] = -1
        _prices[i] = 0.

    while _unassigned_loc > 0:
        n_iter += 1

        for i in prange(n_bidders, nogil=True):
            #with _lck:
            with gil:
                if _p_unassigned_loc[0] == 0:
                    break
                bidder = _unassigned[_p_unassigned_loc[0]]
                _p_unassigned_loc = _p_unassigned_loc[0] - 1

            # BID
            bid(X_view,
                n,
                bidder,
                n_blocks,
                eps,
                _unassigned,
                _p_unassigned_loc,
                _payoffs,
                _prices,
                _assignments,
                _belongs_to)

    if reverse:
        pass # TODO needs to return opposite row & col ind since we operated on the transpose

    row_ind = np.asarray(range(m))
    col_ind = np.zeros(m, dtype=int)
    for i in range(m):
        col_ind[i] = _assignments[i]

    free(_unassigned)
    free(_assignments)
    free(_belongs_to)
    free(_payoffs)
    free(_prices)

    return row_ind, col_ind