0w1

Classic Linear Inequality ( O( 2 ^ ( N + M ) * ( M ^ 3 + M * N ) ) )

I have given up using Simplex Method....
Given an augmented matrix, where each row, say { a, b, c, .. z }, represents a * x1 + b * x2 + .. ≤ z, we will find the maximum value for some vector of coefficient for these variables. A brute force solution is to take advantage of the corollary that maximal / minimal value is always in one of the vertices that all the linear equations forms. That is, we can enumerate each subset of all the vectors, find a solution, and update our answer with that solution.
N : number of variables
M : number of constraints

const double EPS = 1e-9;
const double DINF = 1e9;

struct Gauss_Solve_Matrix{
    vvd mat;
    vd ans;
    Gauss_Solve_Matrix( const vvd &m ){
        mat = m;
        ans.resize( m[ 0 ].size() - 1 );
    }
    void Gauss_elimination(){
        for( int i = 0; i < min( mat.size(), mat[ 0 ].size() - 1 ); ++i ){
            if( mat[ i ][ i ] == 0 )
                for( int j = i + 1; j < mat.size(); ++j )
                    if( mat[ j ][ i ] ){
                        swap( mat[ i ], mat[ j ] );
                        break;
                    }
            if( mat[ i ][ i ] == 0 ) continue;

            for( int j = mat[ i ].size() - 1; j >= i; --j )
                mat[ i ][ j ] /= mat[ i ][ i ];

            for( int j = i + 1; j < mat.size(); ++j )
                for( int k = mat[ j ].size() - 1; k >= i; --k )
                    mat[ j ][ k ] -= mat[ j ][ i ] * mat[ i ][ k ];
        }
    }
    void back_substitution(){
        for( int i = mat.size() - 1; i >= 0; --i ){
            int t; // first non zero coefficient
            for( t = 0; t + 1 < mat[ i ].size(); ++t )
                if( mat[ i ][ t ] )
                    break;
            if( t >= ( int ) mat[ i ].size() - 1 ) continue;
            double s = 0;
            for( int j = t + 1; j + 1 < mat[ i ].size(); ++j )
                s += mat[ i ][ j ] * ans[ j ];
            ans[ t ] = ( mat[ i ].back() - s ) / mat[ i ][ t ];
        }
    }
};

bool valid( const vvd &mat, const vd &sol ){
    for( int i = 0; i < mat.size(); ++i ){
        double s = 0;
        for( int j = 0; j + 1 < mat[ i ].size(); ++j )
            s += mat[ i ][ j ] * sol[ j ];
        if( s > EPS + mat[ i ].back() ) return false;
    }
    return true;
}

void solve(){
    int M, N; cin >> M >> N; // M constrictions, N variables
    vvd mat( N + M + 1, vd( N + 1 ) );
    for( int i = 0; i < M; ++i )
        for( int j = 0; j < N; ++j )
            cin >> mat[ i ][ j ];
    for( int i = 0; i < M; ++i )
        cin >> mat[ i ][ N ];
    for( int i = M; i < M + N; ++i )
        mat[ i ][ i - M ] = -1; // -a[ 0 ~ N - 1 ] <= 0 -> a[ 0 ~ N - 1 ] >= 0
    vd qry( N );
    for( int i = 0; i < N; ++i )
        cin >> qry[ i ];
    for( int i = 0; i < N + 1; ++i )
        mat[ M + N ][ i ] = ( i != N ? qry[ i ] : DINF );

    vd sol;
    double maxv = -DINF;
    for( int S = 1; S < ( 1 << N + M + 1 ); ++S ){
        vvd smat;
        for( int i = 0; i < N + M + 1; ++i )
            if( S >> i & 1 )
                smat.push_back( mat[ i ] );
        Gauss_Solve_Matrix *gsm = new Gauss_Solve_Matrix( smat );
        gsm->Gauss_elimination();
        gsm->back_substitution();
        if( not valid( mat, gsm->ans ) )
            continue;
        double value = 0;
        for( int i = 0; i < gsm->ans.size(); ++i )
            value += gsm->ans[ i ] * qry[ i ];
        if( upmax( maxv, value ) )
            sol = gsm->ans;
    }

    if( fabs( maxv - ( -DINF ) ) < EPS ){
        cout << "No solution\n";
        return;
    }

    if( fabs( maxv - DINF ) < EPS ){
        cout << "Infinity\n";
        return;
    }

    cout << "Bounded solution\n";
    cout << fixed << setprecision( 16 );
    for( int i = 0; i < sol.size(); ++i )
        cout << sol[ i ] << " \n"[ i + 1 == sol.size() ];
}

int main(){
    ios::sync_with_stdio( false );
    solve();
    return 0;
}