Find the set of integers for which two linear equalities are valid

Question

Which algorithm can I use to find the set of all positive integer values n1, n2, ... ,n7 for which the following inequalities hold.

 97n1 + 89n2 + 42n3 + 20n4 + 16n5 + 11n6 + 2n7 - 185 > 0 -98n1 - 90n2 - 43n3 - 21n4 - 17n5 - 12n6 - 3n7 + 205 > 0 n1 >= 0, n2 >= 0, n3 >=0. n4 >=0, n5 >=0, n6 >=0, n7 >= 0 

For example, one set n1= 2, n2 = n3 = ... = n7 =0 makes the inequality true. How to find out all the other values? A similar question was published in M.SE.

ADDED :: I need to generalize the solution for n variables (maybe big). What procedure can I apply? For another special case, n=8

 97n1 + 89n2 + 42n3 + 20n4 + 16n5 + 11n6 + 6n7 + 2n8 - 185 > 0 -98n1 - 90n2 - 43n3 - 21n4 - 17n5 - 12n6 - 7 - 3n8 + 205 > 0 n1 >= 0, n2 >= 0, n3 >=0. n4 >=0, n5 >=0, n6 >=0, n7 >= 0, n8 >= 0 

Python takes forever. Wolfram Mathematica shows that there are 4015 solutions in less than a minute.

 Length[Solve[{97 n1 + 89 n2 + 42 n3 + 20 n4 + 16 n5 + 11 n6 + 6 n7 + 2 n8 - 185 > 0, -98 n1 - 90 n2 - 43 n3 - 21 n4 - 17 n5 - 12 n6 - 7 n7 - 3 n8 + 205 > 0, n1 >= 0, n2 >= 0, n3 >= 0, n4 >= 0, n5 >= 0, n6 >= 0, n7 >= 0, n8 >= 0}, {n1, n2, n3, n4, n5, n6, n7, n8}, Integers]] 

Answer 1

Reti43 has the right idea, but there is a quick recursive solution that works with less strict assumptions about your inequalities.

 def solve(smin, smax, coef1, coef2): """ Return a list of lists of non-negative integers `n` that satisfy the inequalities, sum([coef1[i] * n[i] for i in range(len(coef1)]) > smin sum([coef2[i] * n[i] for i in range(len(coef1)]) < smax where coef1 and coef2 are equal-length lists of positive integers. """ if smax < 0: return [] n_max = ((smax-1) // coef2[0]) solutions = [] if len(coef1) > 1: for n0 in range(n_max + 1): for solution in solve(smin - n0 * coef1[0], smax - n0 * coef2[0], coef1[1:], coef2[1:]): solutions.append([n0] + solution) else: n_min = max(0, (smin // coef1[0]) + 1) for n0 in range(n_min, n_max + 1): if n0 * coef1[0] > smin and n0 * coef2[0] < smax: solutions.append([n0]) return solutions 

You would apply this to your original problem, for example,

 smin, coef1 = 185, (97, 89, 42, 20, 16, 11, 2) smax, coef2 = 205, (98, 90, 43, 21, 17, 12, 3) solns7 = solve(smin, smax, coef1, coef2) len(solns7) 1013 

and to a longer problem like this,

 smin, coef1 = 185, (97, 89, 42, 20, 16, 11, 6, 2) smax, coef2 = 205, (98, 90, 43, 21, 17, 12, 7, 3) solns8 = solve(smin, smax, coef1, coef2) len(solns8) 4015 

On my Macbook, both of these cases are resolved in milliseconds. This should scale well enough for a few larger problems, but it’s basically O (2 ^ N) in terms of the number of factors N. How well it actually scales depends on how large the additional coefficients are - the more coefficients (compared to smax-smin), the fewer possible solutions, the faster it will work.

Updated . From the discussion of the related message by M.SE, I see that the connection between the two inequalities here is part of the structure of the problem. Given this, a slightly simpler solution can be given. The code below also contains a couple of additional optimizations that speed up the solution for the case with 8 variables from 88 milliseconds to 34 milliseconds on my laptop. I tried it with examples with 22 variables and got the results in less than a minute, but it will not be practical for hundreds of variables.

 def solve(smin, smax, coef): """ Return a list of lists of non-negative integers `n` that satisfy the inequalities, sum([coef[i] * n[i] for i in range(len(coef)]) > smin sum([(coef[i]+1) * n[i] for i in range(len(coef)]) < smax where coef is a list of positive integer coefficients, ordered from highest to lowest. """ if smax <= smin: return [] if smin < 0 and smax <= coef[-1]+1: return [[0] * len(coef)] c0 = coef[0] c1 = c0 + 1 n_max = ((smax-1) // c1) solutions = [] if len(coef) > 1: for n0 in range(n_max + 1): for solution in solve(smin - n0 * c0, smax - n0 * c1, coef[1:]): solutions.append([n0] + solution) else: n_min = max(0, (smin // c0) + 1) for n0 in range(n_min, n_max + 1): solutions.append([n0]) return solutions 

You would apply this for example with 8 variables, for example,

 solutions = solve(185, 205, (97, 89, 42, 20, 16, 11, 6, 2)) len(solutions) 4015 

This solution directly lists the lattice points in a restricted area. Since you want all these solutions, the time required to obtain them will be proportional (at best) to the number of connected lattice points, which grows exponentially with the number of dimensions (variables).

Answer 2

In both examples you gave, I noticed the same pattern. For example, for the first case, if you add two equations together, you will get

 -n1 - n2 - n3 - n4 - n5 - n6 - n7 + 20 > 0 

which can be rearranged to

 n1 + n2 + n3 + n4 + n5 + n6 + n7 < 20 

This is a good limited equation that you can use to iterate over. In particular, you can iterate for n1 from 0 to 19, for n2 from 0 to 19-n1, etc. A possible solution to this may be (0, 0, 0, 0, 0, 0, 0) but we notice that this does not satisfy our original equation. So, just generate all possible values ​​for (n1, n2, ..., n7) and save only those that satisfy your equation. Hardcoding all this leads to

 def find_solutions(N): sols = [] for n1 in xrange(N): for n2 in xrange(N-n1): for n3 in xrange(N-n1-n2): for n4 in xrange(N-n1-n2-n3): for n5 in xrange(N-n1-n2-n3-n4): for n6 in xrange(N-n1-n2-n3-n4-n5): for n7 in xrange(N-n1-n2-n3-n4-n5-n6): if (97*n1 + 89*n2 + 42*n3 + 20*n4 + 16*n5 + 11*n6 + 2*n7 - 185 > 0 and -98*n1 - 90*n2 - 43*n3 - 21*n4 - 17*n5 - 12*n6 - 3*n7 + 205 > 0): sols.append((n1, n2, n3, n4, n5, n6, n7)) return sols 

find_solutions(20) finds all 1013 solutions in 0.6 seconds. Similarly, for the second case, he finds all 4015 solutions in 2.3 seconds. Now this is not at all easy to generalize, but it shows that with a smart approach, Python or any other language should not be slow.

On the other hand, recursion allows us to generalize this to any number of nested loops, but at a slightly slower price.

 def find_solutions(N, coeffs, depth=0, variables=None, subtotal=None, solutions=None): if variables is None: solutions = [] subtotal = [0 for _ in xrange(len(coeffs[0]))] variables = [0 for _ in xrange(len(coeffs[0])-1)] if depth == len(coeffs[0])-2: for v in xrange(N-sum(variables[:depth])): conditions = all( subtotal[i]+coeffs[i][depth]*v > coeffs[i][-1] for i in xrange(len(coeffs)) ) if conditions: variables[depth] = v solutions.append(tuple(variables)) else: for v in xrange(N-sum(variables[:depth])): variables[depth] = v total = [subtotal[i]+coeffs[i][depth]*v for i in xrange(len(coeffs))] find_solutions(N, coeffs, depth+1, variables, total, solutions) if depth == 0: return solutions 

To run this, generate coefficients for each equation and pass them to the function. Keep in mind that the sign of constants is inverted!

 coeffs = [ [97, 89, 42, 20, 16, 11, 2, 185], [-98, -90, -43, -21, -17, -12, -3, -205] ] solutions = find_solutions(20, coeffs) print(len(solutions)) 

This completes the case of n = 7 in 1.6 seconds and the case of n = 8 in 5.8. I will consider any possible optimizations if you expect your n to become very large, but so far it looks satisfactory.

The question remains whether the sum of your equations will always be simplified to something like n1 + n2 + ... nn < N . There is a simple solution around this if it is not, but I decided not to prematurely exceed the code beyond the examples that you provided to us.

Finally, I assume that the same approach could be implemented in Java or C #, and that would probably be even faster. I would not mind doing this if it is expected that your general cases will require much more time.

Answer 3

There are 1013 solutions, but I do not know how the most effective way to solve this problem.

Looking at the second inequality, 17 * n5 cannot be greater than 205 (otherwise the entire left side cannot be positive). This leads to n5 <= 12 . By calculating a similar estimate for each of the other variables, you can reduce the problem to one that can be quickly solved with nested loops.

This Java code displays all the solutions.

 for (int n1 = 0; n1 <= 2; n1++) { for (int n2 = 0; n2 <= 2; n2++) { for (int n3 = 0; n3 <= 4; n3++) { for (int n4 = 0; n4 <= 9; n4++) { for (int n5 = 0; n5 <= 12; n5++) { for (int n6 = 0; n6 <= 17; n6++) { for (int n7 = 0; n7 <= 68; n7++) { if (97 * n1 + 89 * n2 + 42 * n3 + 20 * n4 + 16 * n5 + 11 * n6 + 2 * n7 - 185 > 0 && -98 * n1 - 90 * n2 - 43 * n3 - 21 * n4 - 17 * n5 - 12 * n6 - 3 * n7 + 205 > 0) { System.out.println(Arrays.asList(n1, n2, n3, n4, n5, n6, n7)); } } } } } } } } 

You can make this more efficient by stopping each cycle as soon as the sum of the first few members

 98n1 + 90n2 + 43n3 + 21n4 + 17n5 + 12n6 + 3n7 

reaches 205.

For example, cycle n4 can be replaced by

 for (int n4 = 0; 98 * n1 + 90 * n2 + 43 * n3 + 21 * n4 < 205; n4++) 

If you change all 7 cycles in this way, all 1013 solutions can be found very quickly.

Answer 4

It seems like a bit of linear programming with whole solutions. I think some algorithms are already implemented. Try the Linear Programming Tool / Library for Java .

Answer 5

(up to 13 integers) Here is an ugly brute force with gpgpu (opencl) on 1600 cores, finding 1013 (7 ints) solutions in 9.3 milliseconds, including the time it takes to load an array from gpu into the CPU:

Edit: Fixed n1, n2, n3, as they were 1,20,400 instead of 20,20,20.

 __kernel void solver1(__global __write_only float * subCount){ int threadId=get_global_id(0); int ctr=0; int n1=threadId/400; int n2=(threadId/20)%20; int n3=threadId%20; for(int n4=0;n4<=20;n4++) for(int n5=0;n5<=20;n5++) for(int n6=0;n6<=20;n6++) for(int n7=0;n7<=20;n7++) if ( (97*n1 + 89*n2 + 42*n3 + 20*n4 + 16*n5 + 11*n6 + 2*n7 - 185 > 0) && (-98*n1 - 90*n2 - 43*n3 - 21*n4 - 17*n5 - 12*n6 - 3*n7 + 205 > 0)) {ctr++;} subCount[threadId]=ctr; } 

then something like subCount.Sum () gives the number of solutions. (in microseconds)

global working size = 8000 (threads) (getting n4 in this will make it 160,000 and increase performance)

local work size = 160 (inefficient on my machine, which makes its power-2 optimal)

This is just a string that you need to send to gpu to compile. You simply add additional loops (or just threadId relationships) to the string, such as n8, n9, and change the body of the "if" to sum them up.

Edit: Adding 1 extra integer to the equations increases the solution time to 101 milliseconds even with greater optimization (finds 4015 solutions).

  __kernel void solver2(__global __write_only float * subCount){ int threadId=get_global_id(0); int ctr=0; int n1=threadId/160000; int c1n1=97*n1; int c2n1=-98*n1; int n2=(threadId/8000)%20; int c1n2=89*n2; int c2n2=- 90*n2 ; int n3=(threadId/400)%20; int c1n3=42*n3; int c2n3=- 43*n3 ; int n4=(threadId/20)%20; int c1n4=20*n4 ;int c2n4=- 21*n4 ; int n5=threadId%20; int c1n5=16*n5 ;int c2n5=- 17*n5 ; int t1=c1n1+c1n2+c1n3+c1n4+c1n5; int t2=c2n1+c2n2+c2n3+c2n4+c2n5; for(int n6=0;n6<=20;n6++) for(int n7=0;n7<=20;n7++) for(int n8=0;n8<=20;n8++) if(t1+ 11*n6 + 2*n7+6*n8 > 185 && t2 - 12*n6 - 3*n7-7*n8 > -205) ctr++; subCount[threadId]=ctr; } 

workgroup size = 256 global work size = 3200000 threads

Editing: 9-integer version with just one solution found for 14 cycles for c1 = 3 and c2 = -4 in 1581 milliseconds. Changing the upper bound from 20 to 19 for n6-n9 gave the same result in 1310 milliseconds.

Edit: Adding some of the answers from Reti43 and Molecule , 9 integer versions took 14 milliseconds (but still inefficient for the first 5 integers) using 3.2M streams.

 __kernel void solver3(__global __write_only float * subCount){ int threadId=get_global_id(0); int ctr=0; int n1=threadId/160000; int c1n1=97*n1; int c2n1=-98*n1; int n2=(threadId/8000)%20; int c1n2=89*n2; int c2n2=- 90*n2 ; int n3=(threadId/400)%20; int c1n3=42*n3; int c2n3=- 43*n3 ; int n4=(threadId/20)%20; int c1n4=20*n4 ;int c2n4=- 21*n4 ; int n5=threadId%20; int c1n5=16*n5 ;int c2n5=- 17*n5 ; int t1=c1n1+c1n2+c1n3+c1n4+c1n5; int t2=c2n1+c2n2+c2n3+c2n4+c2n5; int m=max( max( max( max(n1,n2),n3),n4),n5); for(int n6=0;n6<=20-m;n6++) for(int n7=0;n7<=20-m-n6;n7++) for(int n8=0;n8<=20-m-n6-n7;n8++) for(int n9=0;n9<=20-m-n6-n7-n8;n9++) if(t1+ 11*n6 + 2*n7+6*n8 +3*n9> 185 && t2 - 12*n6 - 3*n7-7*n8-4*n9 > -205) ctr++; subCount[threadId]=ctr; } 

• 10 ints: solutions 46k, 35 ms (c1 = 7, c2 = -8)
• 11 ints: 140k solutions, 103ms (c1 = 1, c2 = -2)
• 12 ints: 383k solutions, 274ms (c1 = 5, c2 = -6)

• probably too slow at 15 tf due to hardware limitations. Could not even compile for 20 ints.

  Note: m = max (max (max (max (n1, n2), n3), n4), n5) is not optimal, but m = n1 + n2 + n3 + n4 + n5.

I will try the Monte Carlo method to overcome bottlenecks in hardware resources. Maybe there is a convergence of solution numbers somehow.

Answer 6

pseudo code:

 for each inequation: find all real roots of the equivalent equation, ie the zero-crossings for each interval between two adjacent roots: pick any number strictly inside the interval evaluate the polynomial in that point if the evaluated polimonial is positive: add every integer in the interval to the list of solutions to that inequation (treat the open-ended intervals outside the extreme roots as special cases, they may contain infinite solutions) find the integers that are in all the lists of solutions to the individual equations