This is the eightieth seventh part of the ILP series. For your convenience you can find other parts in the table of contents in Part 1 – Boolean algebra
We want to create a sequence of numbers of increasing length where each number is a square of some integer and has same digits as previous number in the sequence (plus one additional digit, obviously). Example sequence is:
25,256,5625,65025,105625,9150625
Let’s solve it with ILP:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 |
int firstNumberLength = 2; int sequenceLength = 6; List<IVariable> numbers = new List<IVariable>(); List<IVariable> lastDigitsHash = new List<IVariable>(); for(int i=0;i<sequenceLength;++i){ int numberLength = firstNumberLength + i; IVariable[] digits = Enumerable.Range(0, numberLength).Select(x => solver.CreateAnonymous(Domain.PositiveOrZeroInteger).Set<LessOrEqual>(solver.FromConstant(9))).ToArray(); digits[0].Set<GreaterOrEqual>(solver.FromConstant(1)); IVariable fullNumber = solver.FromConstant(0); long power = 1; for(int j=numberLength-1;j>=0;--j){ fullNumber = fullNumber.Operation<Addition>(digits[j].Operation<Multiplication>(solver.FromConstant(power))); power *= 10; } numbers.Add(fullNumber); IVariable root = solver.CreateAnonymous(Domain.PositiveOrZeroInteger); root.Operation<Multiplication>(root).Set<Equal>(fullNumber); List<IVariable> digitsHash = new List<IVariable>(); for(int j=0;j<=9;++j){ IVariable count = solver.FromConstant(0); for(int k=0;k<numberLength;++k){ count = count.Operation<Addition>(digits[k].Operation<IsEqual>(solver.FromConstant(j))); } digitsHash.Add(count); } if(i > 0){ for(int j=0;j<=9;++j){ digitsHash[j].Operation<Subtraction>(lastDigitsHash[j]).Set<LessOrEqual>(solver.FromConstant(1)).Set<GreaterOrEqual>(solver.FromConstant(0)); } } lastDigitsHash = digitsHash; } solver.AddGoal("g", solver.FromConstant(0)); solver.Solve(); foreach(var n in numbers){ Console.WriteLine(solver.GetValue(n)); } |
We define sequence parameters and the iterate through numbers (line 6).
For each number we first create digits (line 8) and then construct the full number (line 13).
Next, we make the number a root of some integer (line 21).
Next, we count digits (lines 24-32). If it’s another element of the sequence then we compare digits with previous number and make sure only one digit differs (lines 34-40).
Finally, we print the solution. Output:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 |
Tried aggregator 2 times. MIP Presolve eliminated 819 rows and 148 columns. MIP Presolve modified 4802 coefficients. Aggregator did 1276 substitutions. Reduced MIP has 3216 rows, 1765 columns, and 9845 nonzeros. Reduced MIP has 1529 binaries, 236 generals, 0 SOSs, and 0 indicators. Presolve time = 0.01 sec. (19.66 ticks) Probing fixed 610 vars, tightened 521 bounds. Probing changed sense of 47 constraints. Probing time = 0.02 sec. (10.98 ticks) Cover probing fixed 0 vars, tightened 298 bounds. Tried aggregator 2 times. MIP Presolve eliminated 1714 rows and 1024 columns. MIP Presolve modified 607 coefficients. Aggregator did 120 substitutions. Reduced MIP has 1382 rows, 621 columns, and 4249 nonzeros. Reduced MIP has 494 binaries, 127 generals, 0 SOSs, and 0 indicators. Presolve time = 0.02 sec. (11.69 ticks) Probing fixed 1 vars, tightened 60 bounds. Probing time = 0.00 sec. (2.73 ticks) Cover probing fixed 0 vars, tightened 1 bounds. Tried aggregator 1 time. MIP Presolve eliminated 4 rows and 4 columns. MIP Presolve modified 82 coefficients. Reduced MIP has 1378 rows, 617 columns, and 4237 nonzeros. Reduced MIP has 490 binaries, 127 generals, 0 SOSs, and 0 indicators. Presolve time = 0.00 sec. (3.41 ticks) Probing time = 0.01 sec. (2.18 ticks) Cover probing fixed 0 vars, tightened 37 bounds. Clique table members: 2786. MIP emphasis: balance optimality and feasibility. MIP search method: dynamic search. Parallel mode: deterministic, using up to 8 threads. Root relaxation solution time = 0.00 sec. (11.29 ticks) Nodes Cuts/ Node Left Objective IInf Best Integer Best Bound ItCnt Gap 0 0 0.0000 95 0.0000 543 0 0 0.0000 179 Cuts: 244 820 0 0 0.0000 107 Cuts: 51 891 0 0 0.0000 146 Cuts: 252 1058 0 2 0.0000 105 0.0000 1058 Elapsed time = 0.47 sec. (328.37 ticks, tree = 0.01 MB, solutions = 0) 22 22 0.0000 141 0.0000 9845 85 65 0.0000 184 0.0000 34909 316 166 0.0000 106 0.0000 75479 982 267 0.0000 72 0.0000 107310 1450 350 0.0000 61 0.0000 143017 1940 410 0.0000 51 0.0000 178585 2411 459 0.0000 34 0.0000 216550 2810 505 0.0000 39 0.0000 247172 3303 520 0.0000 31 0.0000 286418 5107 548 0.0000 89 0.0000 418904 Elapsed time = 5.47 sec. (3571.19 ticks, tree = 0.12 MB, solutions = 0) 6677 514 0.0000 53 0.0000 546687 6904 506 0.0000 179 -0.0000 568526 7670 266 0.0000 82 -0.0000 616168 8937 320 0.0000 110 -0.0000 703318 10053 687 0.0000 85 -0.0000 777066 11072 233 0.0000 39 -0.0000 808490 13952 351 0.0000 121 -0.0000 879987 16596 613 0.0000 96 -0.0000 970191 19124 995 0.0000 45 -0.0000 1071157 21397 1402 0.0000 56 -0.0000 1178302 Elapsed time = 23.98 sec. (14652.33 ticks, tree = 0.26 MB, solutions = 0) 25537 1576 0.0000 69 -0.0000 1293558 27855 2005 0.0000 62 -0.0000 1397107 29683 2393 0.0000 84 -0.0000 1512387 31656 2629 infeasible -0.0000 1620551 33486 2930 infeasible -0.0000 1734357 35221 3248 0.0000 57 -0.0000 1843613 37566 3448 0.0000 103 -0.0000 1966558 40724 3750 infeasible -0.0000 2055454 * 40763 3536 integral 0 0.0000 -0.0000 2056606 0.00% Root node processing (before b&c): Real time = 0.47 sec. (327.00 ticks) Parallel b&c, 8 threads: Real time = 44.83 sec. (22037.17 ticks) Sync time (average) = 4.44 sec. Wait time (average) = 0.00 sec. ------------ Total (root+branch&cut) = 45.30 sec. (22364.16 ticks) 25 256 5625 65025 105625 9150625 |