How to iteratively build a list?
Let me see if I get this right. So you have
max = 10^12
and you are iterating over the same list of primes with a nested loop. Your inner loop iterates over q
and we should first calculate the largest prime we need in the list to ensure that we are able to hit max
.
Your very first iteration, where p=2
is critical because p^2*q^3
will be the smallest. So how about we find out for which prime q
we make it above the max
mark? Just for clarity, I'm making the code a bit verbose
With[{p = 2},
MinimalBy[Range[2000],
Abs[p^2*Prime[#]^3 - max] &
]
]
(* {819} *)
This tells us that we are passing the max
mark for the 820th prime
2^2*Prime[820]^3
(* 1000664355604 *)
For collecting the values, I suggest you use Reap
and Sow
instead of AppendTo
as it is faster. For the iteration, there are many choices. Let me use a simple Do
in this case
result = Reap[With[
{
primes = Prime[Range[820]]
},
Do[
If[p =!= q && p^2*q^3 < max,
Sow[p^2*q^3]
],
{p, primes},
{q, primes}
]
]][[2, 1]];
That last [[2,1]]
looks a bit clumsy, but when you read how Reap
works, you'll understand its meaning. I got about 21k results and they all have the required form
result[[500]] // FactorInteger
(* {{2, 2}, {3581, 3}} *)
The GeneralUtilities``
package contains undocumented iterator functionality that resembles the generator-based iterators from C++. Perhaps these will make it into the core language some day. Until then, the following iterator-based answer is mainly out of academic interest rather than a practical recommendation...
We start by defining sqube
as in the question:
sqube[{x_,y_}] := x*x*y*y*y
We then define an iterator over all the primes (RangeIterator
does not support Infinity
so we use a googol instead):
primes[] := RangeIterator[1*^100] // MapIterator[Prime]
So then:
primes[] // TakeIterator[10] // Normal
(* {2, 3, 5, 7, 11, 13, 17, 19, 23, 29} *)
We will also define a useful helper iterator (which really ought to be built-in):
takeWhileIterator[crit_] := MapIterator[If[crit[#], #, IteratorExhausted]&]
Armed with these definitions, we can define our results:
results[maxi_] :=
(* for (int p : primes) *)
primes[] //
(* not present in the C++, removes the need to know in advance how many primes *)
takeWhileIterator[sqube[{2, #}] <= maxi &] //
(* for (int q : primes) *)
JoinMapIterator[
( primes[] //
(* if (p != q) *)
SelectIterator[Curry[#2 != # &][#]] //
(* compute p*p*q*q*q *)
MapIterator[Curry[sqube[{#2, #}]&][#]] //
(* if (p*p*q*q*q > maxi) break *)
takeWhileIterator[# <= maxi &]
) &
]
So then:
results[500] // Normal
(* {108, 500, 72, 200} *)
results[5000] // Normal
(* {108, 500, 1372, 72, 1125, 3087, 200, 675, 392, 1323} *)
results[10^12] // Normal // Short
(* { 108, 500, 1372, 5324, 8788, 19652, 27436
, <<20952>>
, 52810620731, 87171249997, 194935071113
, 272147293459, 482754937967, 967692132989
}
*)
% // Last // FactorInteger
(* {{29, 3}, {6299, 2}} *)