<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def solve_mle_rev(C, tol=1e-10, maxiter=100, show_progress=False, full_output=False, return_statdist=True, **kwargs):
 """Number of states"""

    M = C.shape[0]

    """Initial guess for primal-point"""
    z0 = np.zeros(2*M)
    z0[0:M] = 1.0

    """Inequality constraints"""
    # G = np.zeros((M, 2*M))
    # G[np.arange(M), np.arange(M)] = -1.0
    G = -1.0*scipy.sparse.eye(M, n=2*M, k=0)
    h = np.zeros(M)

    """Equality constraints"""
    A = np.zeros((1, 2*M))
    A[0, M] = 1.0
    b = np.array([0.0])

    """Scaling"""
    c0 = C.max()
    C = C/c0

    """Symmetric part"""
    Cs = C + C.T

    """Column sum"""
    c = C.sum(axis=0)

    if scipy.sparse.issparse(C):
        Cs = Cs.tocsr()
        c = c.A1
        A = scipy.sparse.csr_matrix(A)
        F = objective_sparse.F
        DF = objective_sparse.DFsym
        convert_solution = objective_sparse.convert_solution
    else:
        F = objective_dense.F
        DF = objective_dense.DF
        convert_solution = objective_dense.convert_solution

    """PDIP iteration"""
    res = primal_dual_solve(F, z0, DF, A, b, G, h,
                            args=(Cs, c),
                            maxiter=maxiter, tol=tol,
                            show_progress=show_progress,
                            full_output=full_output)
    if full_output:
        z, info = res
    else:
        z = res
    pi, P = convert_solution(z, Cs)

    result = [P]
    if return_statdist:
        result.append(pi)
    if full_output:
        result.append(info)

    return tuple(result) if len(result) > 1 else result[0]