Functional style using lazy lists?

Question

Let's say I want to answer the question "what are the first 400 palindromic prime numbers?"

The first approach that comes to my mind from the set of languages that I know is to use some sort of lazy list materialization, a la IEnumerable (and yield) in C#, generators in Python, or sequence blocks in F#.

For example, in C#:

PrimesEnumerator().Where(n => n.ToString() == n.ToString().Reverse()).Take(400);

This would cause the PrimesGenerator to be pumped for primes long enough for the Where() clause to find enough numbers that meet the requirement for Take() to meet its quota.

The best I've come up with in Mathematica is something like:

i = 1; results = List[];
While[Length[results] != 400,
  If[IntegerDigits[Prime[i]] == Reverse[IntegerDigits[Prime[i]]],
    results = Append[results,Prime[i]]];
  i = i + 1]

It surprises me that I end up writing in such an imperative style in Mathematica. Am I missing something that would let me write this entirely functionally? Maybe even with lazy lists?

Update: I took inspiration from WReach's work of art answer, and made a package that took his ideas and expanded them into a broad, general solution for lazy data in Mathematica. I describe its usage in an answer below.

Answer 1

A "lazy list", "functional style" solution to this problem might look something like this:

sIntegers[] ~sMap~ Prime ~sFilter~ palindromicQ ~sTake~ 400 // sList

No such notation is built into Mathematica. However, creating such notations is Mathematica's strong suit. Let's do it.

First, we need to define the notion of a "stream". Streams are inherently lazy, so let's use HoldAll:

SetAttributes[stream, {HoldAll}]

A stream can be empty:

sEmptyQ[stream[]] := True

... or it can be non-empty, having two elements:

sEmptyQ[stream[_, _]] = False;

The first element of the stream is called the "head":

sHead[stream[h_, _]] := h

The remaining elements of the stream are called the "tail":

sTail[stream[_, t_]] := t

Armed with these definitions, we can now express an infinite stream of integers thus:

sIntegers[n_:1] :=
  With[{nn = n+1}, stream[n, sIntegers[nn]]]

sIntegers[] // sEmptyQ                 (* False *)
sIntegers[] // sHead                   (* 1 *)
sIntegers[] // sTail // sHead          (* 2 *)
sIntegers[] // sTail // sTail // sHead (* 3 *)

Infinite streams are difficult to display in a notebook. Let's introduce sTake which truncates a stream to a fixed length:

sTake[s_stream, 0] := stream[]
sTake[s_stream, n_] /; n > 0 :=
  With[{nn = n-1}, stream[sHead[s], sTake[sTail[s], nn]]]

Let's also introduce sList, which converts a (finite) stream into a list:

sList[s_stream] :=
  Module[{tag}
  , Reap[
      NestWhile[(Sow[sHead[#], tag]; sTail[#])&, s, !sEmptyQ[#]&]
    , tag
    ][[2]] /. {l_} :> l
  ]

Now we can inspect an integer stream directly:

sIntegers[] ~sTake~ 10 // sList
(* {1, 2, 3, 4, 5, 6, 7, 8, 9, 10} *)

sMap applies a function to every element of a stream:

sMap[stream[], _] := stream[]
sMap[s_stream, fn_] := stream[fn[sHead[s]], sMap[sTail[s], fn]]

sIntegers[] ~sMap~ Prime ~sTake~ 10 // sList
(* {2, 3, 5, 7, 11, 13, 17, 19, 23, 29} *)

sFilter selects elements from a stream that satisfy a given filter predicate:

sFilter[s_, pred_] :=
  NestWhile[sTail, s, (!sEmptyQ[#] && !pred[sHead[#]])&] /.
    stream[h_, t_] :> stream[h, sFilter[t, pred]]

sIntegers[] ~sFilter~ OddQ ~sTake~ 15 // sList
(* {1, 3, 5, 7, 9, 11, 13, 15, 17, 19, 21, 23, 25, 27, 29} *)

We now have almost all of the pieces in place to address the original problem. All that is missing is a predicate that detects palindromic numbers:

palindromicQ[n_] := IntegerDigits[n] /. d_ :> d === Reverse[d]

palindromicQ[123] (* False *)
palindromicQ[121] (* True *)

Now, we can solve the problem:

sIntegers[] ~sMap~ Prime ~sFilter~ palindromicQ ~sTake~ 400 // sList

(* {2,3,5,7,11,101, ... ,3528253,3541453,3553553,3558553,3563653,3569653} *)

The stream facility we have defined here is very basic. It lacks error checking, and further consideration should be given to optimization. However, it demonstrates the power of Mathematica's symbolic programming paradigm.

The following listing gives the complete set of definitions:

ClearAll[stream]
SetAttributes[stream, {HoldAll, Protected}]

sEmptyError[] := (Message[stream::empty]; Abort[])
stream::empty = "Attempt to access beyond the end of a stream.";

ClearAll[sEmptyQ, sHead, sTail, sTake, sList, sMap, sFilter, sIntegers]

sEmptyQ[stream[]] := True
sEmptyQ[stream[_, _]] = False;

sHead[stream[]] := sEmptyError[]
sHead[stream[h_, _]] := h

sTail[stream[]] := sEmptyError[]
sTail[stream[_, t_]] := t

sTake[s_stream, 0] := stream[]
sTake[s_stream, n_] /; n > 0 :=
  With[{nn = n-1}, stream[sHead[s], sTake[sTail[s], nn]]]

sList[s_stream] :=
  Module[{tag}
  , Reap[
      NestWhile[(Sow[sHead[#], tag]; sTail[#])&, s, !sEmptyQ[#]&]
    , tag
    ][[2]] /. {l_} :> l
  ]

sMap[stream[], _] := stream[]
sMap[s_stream, fn_] := stream[fn[sHead[s]], sMap[sTail[s], fn]]

sFilter[s_, pred_] :=
  NestWhile[sTail, s, (!sEmptyQ[#] && !pred[sHead[#]])&] /.
    stream[h_, t_] :> stream[h, sFilter[t, pred]]

sIntegers[n_:1] :=
  With[{nn = n+1}, stream[n, sIntegers[nn]]]



palindromicQ[n_] := IntegerDigits[n] /. d_ :> d === Reverse[d]

Answer 2

One way to get the lazy aspect is to use a closure, or the closest way for Mathematica to fake a closure.

This is the closures constructor:

makePalindromePrimeC[start_: 1] := Module[{p = Prime[start], r},
  ((r = NestWhile[NextPrime, p, 
       With[{d = IntegerDigits[#]}, d != Reverse[d]] &]); 
    p = NextPrime[r]; r) &]

This creates one:

palPrimeClosure = makePalindromePrimeC[]

Now you use it to generate some:

In[259]:= Table[palPrimeClosure[], {100}]

Out[259]= {2, 3, 5, 7, 11, 101, 131, 151, 181, 191, 313, 353, 373, \
383, 727, 757, 787, 797, 919, 929, 10301, 10501, 10601, 11311, 11411, \
12421, 12721, 12821, 13331, 13831, 13931, 14341, 14741, 15451, 15551, \
16061, 16361, 16561, 16661, 17471, 17971, 18181, 18481, 19391, 19891, \
19991, 30103, 30203, 30403, 30703, 30803, 31013, 31513, 32323, 32423, \
33533, 34543, 34843, 35053, 35153, 35353, 35753, 36263, 36563, 37273, \
37573, 38083, 38183, 38783, 39293, 70207, 70507, 70607, 71317, 71917, \
72227, 72727, 73037, 73237, 73637, 74047, 74747, 75557, 76367, 76667, \
77377, 77477, 77977, 78487, 78787, 78887, 79397, 79697, 79997, 90709, \
91019, 93139, 93239, 93739, 94049}

Generate some more:

In[260]:= Table[palPrimeClosure[], {50}]

Out[260]= {94349, 94649, 94849, 94949, 95959, 96269, 96469, 96769, \
97379, 97579, 97879, 98389, 98689, 1003001, 1008001, 1022201, \
1028201, 1035301, 1043401, 1055501, 1062601, 1065601, 1074701, \
1082801, 1085801, 1092901, 1093901, 1114111, 1117111, 1120211, \
1123211, 1126211, 1129211, 1134311, 1145411, 1150511, 1153511, \
1160611, 1163611, 1175711, 1177711, 1178711, 1180811, 1183811, \
1186811, 1190911, 1193911, 1196911, 1201021, 1208021}

Now create an entirely independent instance that starts searching at the 500th prime:

In[261]:= palPrimeClosure500 = makePalindromePrimeC[500]

Out[261]= (r$10054 = 
   NestWhile[NextPrime, p$10054, 
    With[{d = IntegerDigits[#1]}, d != Reverse[d]] &]; 
  p$10054 = NextPrime[r$10054]; r$10054) &

In[262]:= Table[palPrimeClosure500[], {30}]

Out[262]= {10301, 10501, 10601, 11311, 11411, 12421, 12721, 12821, \
13331, 13831, 13931, 14341, 14741, 15451, 15551, 16061, 16361, 16561, \
16661, 17471, 17971, 18181, 18481, 19391, 19891, 19991, 30103, 30203, \
30403, 30703}

Answer 3

I took inspiration from WReach's work of art answer, and made a package that took his ideas and expanded them into a broad, general solution for lazy data in Mathematica. You can find my implemenation on github.

To use the package to answer the original question from this post, you'd do something like:

palindromicQ[n_] := IntegerDigits[n] /. d_ :> d === Reverse[d]

Needs["Lazy`"]
Lazy[Primes] ~Select~ palindromicQ ~Take~ 400 // List

It can also do things like triangular numbers (from Project Euler #12):

divisorsLength[n_] := Apply[Times, #[[2]] + 1 & /@ FactorInteger[n]];

Needs["Lazy`"]
triangles = FoldList[Plus, 0, Lazy[Integers]];
triangles ~Select~ (divisorsLength[#] > 500 &) // First

Or some of the even kookier Project Euler questions:

Needs["Lazy`"]
Rest[Lazy[Integers]]~Take~ 9999  ~Select~
  ((Total@Most@Divisors@Total@Most@Divisors[#] === #) &) ~Select~
  ((Total@Most@Divisors[#] =!= #) &) // Total

Answer 4

Not the "lazy" but shortest:

Select[ToString /@ Prime[Range[10^4]], # == StringReverse[#] &]

Faster:

Select[Prime[Range[10^4]], IntegerDigits[#] == Reverse[IntegerDigits[#]] &]

Faster:

Pick[#, (# == Reverse[#]) & /@ IntegerDigits /@ #, True] &@ Prime[Range[10^4]]

Fastest: (twice faster than the rest, but slower than "lazy")

Select[Pick[#,(#==Reverse[#])&/@ IntegerDigits/@ #, True]&@ Range[10^6], PrimeQ]

Result:

{2, 3, 5, 7, 11, 101, 131, 151, 181, 191, 313, 353, 373, 383, 727, 757, 787, 797, 919, 929, 10301, 10501, 10601, 11311, 11411, 12421, 12721, 12821, 13331, 13831, 13931, 14341, 14741, 15451, 15551, 16061, 16361, 16561, 16661, 17471, 17971, 18181, 18481, 19391, 19891, 19991, 30103, 30203, 30403, 30703, 30803, 31013, 31513, 32323, 32423, 33533, 34543, 34843, 35053, 35153, 35353, 35753, 36263, 36563, 37273, 37573, 38083, 38183, 38783, 39293, 70207, 70507, 70607, 71317, 71917, 72227, 72727, 73037, 73237, 73637, 74047, 74747, 75557, 76367, 76667, 77377, 77477, 77977, 78487, 78787, 78887, 79397, 79697, 79997, 90709, 91019, 93139, 93239, 93739, 94049, 94349, 94649, 94849, 94949, 95959, 96269, 96469, 96769, 97379, 97579, 97879, 98389, 98689}

Answer 5

Here is an attempt to do this in a more "functional" way, without thought of efficiency:

pQ = # === StringReverse@# & @ IntegerString@# &;

ppf[x_?PrimeQ] /; pQ[x] := x
ppf[x_] := ppf @ NextPrime @ x

$IterationLimit = 1*^6;

NestList[ppf[# + 1] &, 2, 399]

Answer 6

I recently published a lazy list package on GitHub, so I might as well show it off here:

https://github.com/ssmit1986/lazyLists

Installation instructions can be found on the GitHub landing page above (or in the README.md file).

My code overloads the normal list processing functions like Map, MapIndexed, MapThread, FoldList, Pick, Cases, and Select. lazyRange[min, step] is the basic constructor for generating infinite ranges of integers. Take will return a lazyList in which the first argument is the extracted elements and the second is the tail of the list (which can be used to continue extracting more elements). I implemented it using ReplaceRepeated, since that seemed to be the fastest option I could find.

Example:

<< lazyLists`;
Take[Select[Map[Prime, lazyRange[]], PalindromeQ], 400]

Out[138]= lazyList[{<<primes>>}, <<tail>>]

Edit I'm still updating this repository with new functionalities and I'm open to suggestions. See the README file on the landing page for the update history.

Answer 7

A bit elaborate, but this works nicely for generating the first four hundred palindromic primes:

NestList[
  Function[p, NestWhile[NextPrime, p,
               (IntegerDigits[#] != Reverse[IntegerDigits[#]]) &, {2, 1}]],
         2, 399]

---

An alternative route:

NestList[
  Function[p, FixedPoint[NextPrime, p,
               SameTest -> (IntegerDigits[#2] == Reverse[IntegerDigits[#2]] &)]],
         2, 399]

This alternate is a bit faster than the previous snippet on my system, but you should do your own tests.

Answer 8

There are many possible ways. Here is one in a functional style:

PalindromicPrimeQ[k_Integer] := 
    IntegerDigits[k] == Reverse[IntegerDigits[k]] && PrimeQ[k]
Take[Select[Range[4 10^6], PalindromicPrimeQ], 400]

returns a list of length 400 like this:

{2,3,...101,131,...3558553, 3563653, 3569653}

Updating after rcollyer's suggestions I would rather follow this way :

PalindromicQ[k_Integer] := 
    IntegerDigits[k] == Reverse[IntegerDigits[k]]
Take[Select[Prime[Range[300000]], PalindromicQ], 400]; 

Edit

It is worty to point out that Spartacus' test pQ seems to be a bit faster than PalindromicQ:

pQ = # === StringReverse@# &@IntegerString@# &;
Take[Select[Prime[Range[300000]], pQ], 400]; // Timing

Its timing : 2.152 versus 2.465 in case of PalindromicQ

Answer 9

I like WReach's answer because it shows how to compose an expression that is similar to the example in the question. I like Sal Mangano's answer because it is concise and shows how to make something that behaves like a C# enumerator. I hope to provide a little of each.

Enumerator[state_:0, increment_:(# + 1 &)] := Module[
  {s = state},
  (s = increment[s]) &
];

Where[predicate_][enumerator_] := Function[
  NestWhile[enumerator, enumerator[], Not @* predicate]
];

TakeEnumerator[count_] := Table[#[], count] &;

PalindromicQ := PalindromeQ @* IntegerDigits;

Please note that Enumerator and Where both produce pure functions that behave like C# enumerators but TakeEnumerator does not.

---

With that the example in the question can be expressed as

Enumerator[0, NextPrime] // Where[PalindromicQ] // TakeEnumerator[400]

Result:

{2,3,...101,131,...3558553, 3563653, 3569653}

---

Like in Sal's answer, you can assign to something that behaves as an enumerator.

palindrome = Enumerator[0, NextPrime] // Where[PalindromicQ];
palindrome // TakeEnumerator[10]

Result:

{2, 3, 5, 7, 11, 101, 131, 151, 181, 191}

And repeat for the next 10:

palindrome // TakeEnumerator[10]

Result:

{313, 353, 373, 383, 727, 757, 787, 797, 919, 929}

---

And you can quickly compose enumerators:

Enumerator[] /* (#^2 &) // TakeEnumerator[10]

Result:

{1, 4, 9, 16, 25, 36, 49, 64, 81, 100}

or

Enumerator[] /* (#^2 &) /* Sqrt // TakeEnumerator[10]

Result:

{1, 2, 3, 4, 5, 6, 7, 8, 9, 10}