a = {{5, 3}, {1, 2}} b = {{1, 1}, {1, 0}} c = {{1, 1}, {1, 1}} some functions to check the criteria: c1[result_] := And @@ Table[ c[[i, 1]] >= Total@result[[All, i]] >= c[[i, 2]] , {i, Length@result[[1]]}] c2[result_] := And @@ Table[ b[[i, 1]] >= Total@result[[i]] >= b[[i, 2]] , {i, Length@result}] and calculate the score: score[result_] := Total@MapThread[ Dot, { a, result }] now for this small problem `NMaximise` gets the solution.. NMaximize[{score[result], c1[result] && c2[result] && And @@ (0 <= # <= 1 & /@ Flatten[result])} , {v1, v2, v3, v4}, Integers ] > {7., {v1 -> 1, v2 -> 0, v3 -> 0, v4 -> 1}} result /. %[[2]] > {{1, 0}, {0, 1}}