Limit eating all the memory b / c of unnecessary symbolic work

Question

The following code is a naive way to find the smallest number whose square has n divisors (the minimum should be its log, and x_i should be powers in its simple factorization). If I look at the case n = 2000 and use ten variables instead of twenty, this uses somewhere around 600 MB of memory. With a value of n, I'm actually trying to find the answer, I need about 20 variables so as not to miss the actual solution, and it quickly uses up all available memory, and then replaces the trash.

n=8*10^6; a = Table[N[Log[Prime[i]]], {i, 20}]; b = Table[Subscript[x, i], {i, 20}]; cond = Fold[And, Product[2 Subscript[x, i] + 1, {i, 20}] > n, Table[Subscript[x, i] >= 0, {i, 20}]] && b \[Element] Integers; NMinimize[{ab, cond}, b, MaxIterations -> 1000] 

It turns out that the problem is not at all related to whole programming, etc. (removing the limit on integers doesn't help).

My best guess is that the problem is that Mathematica is wasting all the memory expanding Product[2 Subscript[x, i] + 1, {i, 20}] . If I replaced the product with only Product[Subscript[x, i],{i,20}] and changed the restrictions to >= 1 , not 0 , I get results without fuss and without a kernel that uses more than 50 MB of memory. (Although this preserves the limitation of inequality and does not change the problem of minimizing the objective function, this gives rise to the requirement of integrality - I get even results that correspond to half-integers in the real problem.)

stack overflow.squite in their case, they had an objective function, which was evaluated symbolically with huge costs. They were able to fix this by forcing the function to accept only numerical input, effectively hiding it from Mathematica "I have an Expand [] hammer and everything looks like a nail trend. But you cannot hide the restriction of such a function (Mathematica will complain about the invalid restriction).

Any thoughts on how to fix this?

Edit: I know the correct answer - after my Mathematica code did not work, I used AMPL and the special MINLP solver to get it (pretty fast). I just want to know how I can ever hope that I can use Mathematica's built-in nonlinear optimization functions in the future, despite the crazy things that I think are related to my limitations when I introduced them in the only way I know .

Answer 1

May prevent this condition from doing anything if the inputs are not numeric, as shown below.

 n = 8*10^6; nvars = 20; a = Table[N[Log[Prime[i]]], {i, nvars}]; b = Table[Subscript[x, i], {i, nvars}]; c1[xx : {_?NumericQ ..}] := Times @@ (2 xx + 1); c2 = Map[# >= 0 &, b]; c3 = b \[Element] Integers; cond = Join[{c1[b] > n}, c2, {c3}]; In[29]:= Timing[NMinimize[{ab, cond}, b, MaxIterations -> 400]] Out[29]= {43.82100000000008, {36.77416664719056, {Subscript[x, 1] -> 3, Subscript[x, 2] -> 3, Subscript[x, 3] -> 2, Subscript[x, 4] -> 2, Subscript[x, 5] -> 1, Subscript[x, 6] -> 1, Subscript[x, 7] -> 1, Subscript[x, 8] -> 1, Subscript[x, 9] -> 1, Subscript[x, 10] -> 1, Subscript[x, 11] -> 1, Subscript[x, 12] -> 1, Subscript[x, 13] -> 0, Subscript[x, 14] -> 0, Subscript[x, 15] -> 0, Subscript[x, 16] -> 0, Subscript[x, 17] -> 0, Subscript[x, 18] -> 0, Subscript[x, 19] -> 0, Subscript[x, 20] -> 0}}} 

--- edit ---

I think I would like to point out that this can be posed as a problem of integer linear programming. We use variables 0-1 for all possible combinations of primes and degrees.

We can limit the number of primes using the fact that the solution cannot contain more primes than the minimum necessary, assuming that they are all raised to the first degree. The goal is then minimal if they start with 2.

We will assume that the maximum indicator is not more than 20. This is probably a convenient way to show this, but it has not yet occurred to me.

The purpose and limitations of this statement are given below. We get a completely linear problem by taking the logs equations of the sigma divisor.

 n = 8*10^6; nprimes = Ceiling[Log[2, n]]; maxexpon = 20; vars = Array[x, {maxexpon, nprimes}]; fvars = Flatten[vars]; c1 = Map[0 <= # <= 1 &, fvars]; c2 = {Element[fvars, Integers]}; c3 = Thread[Total[vars] <= 1]; c4 = {Total[N[Log[2*Range[maxexpon] + 1]].vars] >= N@Log [n]}; constraints = Join[c1, c2, c3, c4]; obj = Range[maxexpon].vars.N[Log[Prime[Range[nprimes]]]]; Timing[{min, vals} = NMinimize[{obj, constraints}, fvars];] Out[118]= {5.521999999999991, Null} Pick[fvars, fvars /. vals, 1] /. x[j_, k_] :> {Prime[k], j} Out[119]= {{11, 1}, {13, 1}, {17, 1}, {19, 1}, {23, 1}, {29, 1}, {31, 1}, {37, 1}, {5, 2}, {7, 2}, {2, 3}, {3, 3}} 

This approach processes n = 10 ^ 10 for about 23 seconds.

--- end of editing ---

Daniel Lichtblau

Answer 2

Instead, you can try:

 Catch[Do[If[DivisorSigma[0, k^2] > 2000, Throw[k]], {k, 1000000}]] 

which returns 180180.

---

Addition:

 Catch[Do[If[ Times@ @(2 FactorInteger[k][[All, 2]] + 1) > 2000, Throw[k]], {k, 1000000}]] 

It seems to be faster.

---

APPENDIX 2:

For this ultra-improved version (but not 100% proven):

 MinSquareWithDivisors[n_] := Min@Select [ Product[Prime[k]^#[[k]], {k, 1, Length[#]}] & /@ Flatten[IntegerPartitions /@ Range[Ceiling[Log[2, n]]], 1], DivisorSigma[0, #^2] > n &] 

MinSquareWithDivisors[2000000000] gives 2768774904222066200260800 in ~ 4 seconds

Explanation

First of all, it is necessary to prove that the sum of simple indicators in this minimum number does not exceed Log[2, n] . I have not yet found proof , but this may be due to the relationship between successive strokes.

Flatten[IntegerPartitions /@ Range[Ceiling[Log[2, n]]], 1] gives you all the lists with Total <= Log[2, n] , conveniently sorted from large to small.

Product[Prime[k]^#[[k]], {k, 1, Length[#]}] & use them as prime measures to create integers.

Min@Select [..., DivisorSigma[0, #^2] > n &] select the minimum one that is consistent with the original condition.