ils.py 文件源码

def mils(A, B, y, p=1):
    # x_hat,z_hat = mils(A,B,y,p) produces p pairs of optimal solutions to
    #               the mixed integer least squares problem min_{x,z}||y-Ax-Bz||, 
    #               where x and z are real and integer vectors, respectively.
    #
    # Input arguments:
    #    A - m by k real matrix
    #    B - m by n real matrix
    #          [A,B] has full column rank
    #    y - m-dimensional real vector
    #    p - the number of optimal solutions
    #
    # Output arguments:
    #    x_hat - k by p real matrix
    #    z_hat - n by p integer matrix (in double precision). 
    #           The pair {x_hat(:,j),z_hat(:,j)} is the j-th optimal solution
    #           i.e., its residual is the j-th smallest, so
    #           ||y-A*x_hat(:,1)-B*z_hat(:,1)||<=...<=||y-A*x_hat(:,p)-B*z_hat(:,p)||

    m, k = A.shape
    m2, n = B.shape
    if m != m2 or m != len(y) or len(y[1]) != 1:
        raise ValueError("Input arguments have a matrix dimension error!")

    if rank(A) + rank(B) < k + n:
        raise ValueError("hmmm...")

    Q, R = qr(A, mode='complete')
    Q_A = Q[:, :k]
    Q_Abar = Q[:, k:]
    R_A = R[:k, :]

    # Compute the p optimal integer least squares solutions
    z_hat = ils(dot(Q_Abar.T, B), dot(Q_Abar.T, y), p)

    # Compute the corresponding real least squares solutions
    x_hat = lstsq(R_A, dot(Q_A.T, (dot(y, ones((1, p))) - dot(B, z_hat))))

    return x_hat, z_hat