# 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;
}
```