MathGroup Archive 2008

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Re: How do little quickest ?


Artur wrote:
> Now the quickest procedure is:
> 
> aa = {}; Do[Print[x]; rmin = 10^10; k = 2^x; w = Floor[(k - 1)/2];
>  Do[m = FactorInteger[k n (k - n)]; rad = 1;
>     Do[rad = rad m[[s]][[1]], {s, 1, Length[m]}];
>    If[rad < rmin, rmin = rad], {n, 1, w, 2}];
>  AppendTo[aa, rmin], {x, 2, 30}]; aa
> 
> because all possible numbers was odd <{n, 1, w, 2} inspite 
> {n, 1, w}> and GCD checking was delated because GCD[2^x-(2n-1),(2n-1)]=1 
> 
> What to do yet because time is still expotential and for x=27 my computer need whole day to count.
> 
> Best wishes
> Artur
> 
> Mathematical problem is following:
> 
> I need help with finding rule on follwing problem
> 
> Sequences:
> 
> B={1, 3, 7, 15, 27, 63, 125, 243, 343, 999, 1805, 3721, 8181, 16335,
> 
> 32761, 65533, 112847, 190269, 519375, 1046875, 1953125}
> 
> C={2, 4, 8, 16, 32, 64, 128, 256, 512, 1024, 2048, 4096, 8192, 16384,
> 
> 32768, 65536, 131072, 262144, 524288, 1048576, 2097152}
> 
> A={1, 1, 1, 1, 5, 1, 3, 13, 169, 25, 243, 375, 11, 49, 7, 3, 18225,
> 
> 71875, 4913, 1701, 144027}
> 
> C(n) = 2^n
> 
> A(n) = C(n) - B(n)
> 
> A(n) < B(n) < C(n)
> 
> is obtained by algorhitm such
> 
> that GCD[A(n),B(n),C(n)] = 1 and product of distinct prime divisors of
> 
> A(n)*B(n)*C(n) have minimal value.
> 
> I'm looking for formula for B(n).
> 
> Who can help?
> [...]

There is no way to avoid the exponential aspect unless you have some 
theory to use that allows you to skip large swatches of the potential 
products. Otherwise you need to evaluate twice as many for each 
increment in size. And of course it gets worse if the complexity of 
computing those factorizations grows. (This will not play a role until 
one tries to handle 2^32 or so, as that is where Mathematica will no 
longer always succeed in factoring via trial division.)

The code below is around 3.5 times faster than what I last sent, for the 
size range in question. It avoids direct factorization, instead using a 
standard sieving approach. I opted to work with sums of (approximated) 
logs rather than products of factors, just to be certain I could use 
packed arrays to sufficient size. The disadvantage is that this is 
memory intensive and will require a 64 bit machine in order to have any 
hope for handling the full range you have in mind (2^30).

flatfax = Compile[{{n,_Integer}}, Module[
   {fx, len=2^(n-1), pr=3, start, logpr},
   fx = Table[0.,{len}];
   While[pr<=2*len,
     start = Round[(pr+1)/2.];
     logpr = Log[N[pr]];
     Do [fx[[j]] += logpr, {j,start,len,pr}];
     pr = NextPrime[pr];
     ];
   fx
   ]]

minprods = Compile[{{fx,_Real,1}},
   Module[{len=Length[fx], loglen, n, parts},
     loglen = Floor[Log[2,N[len]]];
     Table[
       n = 2^(j-1);
       parts = Take[fx,2*n];
       parts = Take[parts,n] + Reverse[Drop[parts,n]];
       Min[parts]
       , {j,1,loglen}]
   ]]

In[78]:= Timing[Round[2*Exp[minprods[flatfax[26]]]]]

Out[78]= {168.059, {6, 14, 30, 30, 42, 30, 78, 182, 1110, 570,
   1830, 6666, 2310, 2534, 5538, 9870, 20010, 141270, 14070,
   480090, 155490, 334110, 1794858, 2463270, 2132130}}

In a new session tried to go to 28, and ran out of memory. As I said, 
this will require a large machine and even then I'm not sure one will 
squeak past the memory limits. Alternatively, try working with products 
of prime factors rather than sums of logs of said factors, see if you 
can stay within integer packed array size limits. I think this will 
break down at the minprods step, because some of the products will 
exceed machine integer size even though the minimal ones do not.

One thing I can say is, if you seek further speed improvement, the place 
to focus is flatfax. particularly in the loops. Possibly creating a 
SparseArray vector with Log[pr] in nonzero positions, then adding to the 
fx table (I have no idea whether this will make things faster or slower)


Daniel Lichtblau
Wolfram Research


  • Prev by Date: PDE problem with Wentzell-Robin boundary condition: wrong solution without any warning
  • Next by Date: Visulizing the content of a variable
  • Previous by thread: Re: How do little quickest ?
  • Next by thread: Re: Re: Re: How do little quickest