# How do I get this recurrent sequence to show?

What exactly am I doing wrong? The dsf function takes the digits of n, raises them to their own power, then adds them together.

Here is my dsf function.

dsf[n_] := Block[{m = n, t}
, t = IntegerDigits[m]; Sum[Max[1, t[[k]]]^t[[k]]
, {k, Length[t]}]
]


If I type dsf[256], it correctly outputs 49785.

If I type

RecurrenceTable[{c[n]==dsf[c[n-1]],c[1]==2},c,{n,1,4}]


It gives me the output {2,4,256, 32317006071311...), with the 4th term being a very big number. What's wrong?

• Why does it replace all 0 with 1? Do you want that behavior or iz 0 preferred? So dsf[100] is 1 or 3? Apr 25, 2017 at 21:02
• @volcanrb Your function dsf[] does NOT do what you describe. You have dsf[10] = 2 which is clearly wrong, since it mus be 1^2 = 1 Apr 26, 2017 at 19:29
• It does do what it describes. 10 -> 1^1+0^0 -> 1+1 -> 2. Remember that 0^0 =1 Apr 26, 2017 at 20:19
• @volcanrb 0^0 =1 is an additional assumption. You could as well set 0^0 = 0 because we are in the domain of integer numbers. See my modified answer. Apr 26, 2017 at 20:25

EDITS

29.04.17
- Cycle plot extended to base 20 using a code based on Tally[].
- Extend list of cycle lengths to b = 20

28.04.17

27.04.17
- Proof of periodicity

26.04.2017
- Modified function dsfc (assuming 0^0 = 0)
- Periodicity: Tools, results (1)

Introduction

The function dsf as defined in the OP is based on the assumption that a digit 0 raised to is own power gives 1, i.e. 0^0 -> 1.

Personally I would prefer that 0^x -> 0, and that x can be 0 as well. A good reason might be that this definition is now a "natural" extension to that of the sum of fixed powers of the digits.

Also in a closely related standard problem the assumption 0^0 = 1 is discarded (http://mathworld.wolfram.com/MuenchhausenNumber.html).

So I would prefer to see dsf[10] = 1^1 + 0 = 1 rather than dsf[10] = 2.

In the end it's a matter of taste, and I apologize for taking the freedom to provide a solution based on the modified function I call dsfc[], and also present more material on the interesting topic of periodicity.

Development using a modified basic function

Let us define the modified version of the basic function as

dsfc[n_] :=
Block[{t = Select[IntegerDigits[n], #1 > 0 &]}, Plus @@ (t^t)]


Here we have simplified (and modified) the calculation procedure by directy removing zero digits, thus avoiding the tests of Max[]. Furthermore, the sum of digits is calculated from the elegant t^t, avoiding the calculation of the Length[] of t.

Now we have

dsfc[10]

(* Out[3]= 1 *)


Interating the function gives the sequence

a[1] = s
a[k+1] = dsfc[ a[k] ], k = 1, 2, ...


In Mathematica we can generate a sequence conveniently using NestList[].

Example: the first 10 interations with a starting value s = 2 are

a = NestList[dsfc, 2, 10]

(* Out[256]= {2, 4, 256, 49785, 405024629, 387470792, 406668622, \
16964105, 387517185, 35207797, 389894275} *)


The sequence is interesting: although it rises very stepply it does not seem to diverge.

Making the list longer

a = NestList[dsfc, 2, 100];


and plotting it

Reveals ... periodicity!

Repetitions in this list can be found using Tally[]

at=Select[Tally[b], #[[2]] > 1 &]

(* Out[66]= {{2517290, 2}, {388247166, 2}} *)


The first repeated number is

r = at[[1, 1]]

(* Out[70]= 2517290 *)


This is located in a at

pos = Position[a, r] // Flatten

(* Out[46]= {29, 87} *)


The length of the period is then

p = pos[[2]] - pos[[1]]

(* Out[47]= 58 *)


Trying many different starting values I cold not find another period length, except, of course, for the trivial length 1 if s is a pure power of 10.

These results lead to some interesting questions

1. How does the starting value influence the periodicity properties?
2. Specifically, is there always periodicity?
3. Same questions for other number bases (2, 3, ...)

Periodic behaviour: tools

In oder to study the periodic properties of the iterated function we extend the basic function to an arbitrary number basis. We give it the same name but add the number base b as an optional second parameter

dsfc[n_, b_: 10] :=
Block[{t = Select[IntegerDigits[n, b], #1 > 0 &]}, Plus @@ (t^t)]


Instead of using Tally[] we shall employ the function NestWhileList[] and check for a reptition using the test with UnsameQ[].

Here is an illustrative example:

b = 10;
s = 2;
q0 = NestList[dsfc[#, b] &, s, 100]; q0
q = NestWhileList[dsfc[#, b] &, s, UnsameQ[##] &, All]; q
r = Last[q]
pos = Flatten[Position[q, r]]
p = #[[2]] - #[[1]] &[pos]
q1 = NestWhileList[dsfc[#, b] &, r, UnsameQ[##] &, All]; Short[q1]
p1 = Length[q1] - 1;
p == p1

(*

Out[74] q0 = {2, 4, 256, 49785, 405024629, 387470790, 406668617, ..., 34424686}

Out[75] q = {2, 4, 256, 49785, 405024629, 387470790, 406668617, ..., 17600772, 2517290}

Out[76] r = 2517290

Out[77] pos = {28, 86}

Out[78] p = 58

Out[79] q1 = {2517290, 388247166, 34471575, 1653876, 17697224, ..., 17600772, 2517290}

Out[80] = True

*)


Explanation

q0 is the list of just 100 elements of the iteration generated the starting value s = 2

q is the list from the starting value s = 2 to the first repetition of a number. That number r is the last element of q.

pos gives the positions of r in the list. There must be two of these, and their difference is the period p.

We could stop here, as we have all ingredients to carry out the study. But we have added an alternative method to calculate the period: q1 is the list with the starting value s = r. Its length minus one gives an alternative measure p1 of the period which is verified to be equal to p.

Periodic behaviour: results

The period lengths for the first few number bases in the format

{b, {list of period lengths}}


are

tpp = Table[
q2 = Table[
q = NestWhileList[dsfc[#, b] &, s, UnsameQ[##] &, All];
r = Last[q];
pos = Flatten[Position[q, r]];
p = #[[2]] - #[[1]] &[pos];
{s, p}, {s, 2, 500}];
pp = Union[#[[2]] & /@ q2];
{b, pp}, {d, 2, 15}]

{
{2, {1}},
{3, {1, 2}},
{4, {1}},
{5, {1, 2, 4}},
{6, {1, 5, 7}},
{7, {1, 4, 14, 17}},
{8, {1, 6}},
{9, {1, 2, 5, 7, 8, 11, 18, 29}},
{10, {1, 25, 52, 58}},
{11, {1, 2, 4, 55, 68, 113}},
{12, {1, 25, 48, 115, 494}},
{13, {1, 3, 66, 186, 195, 277, 299, 408}},
{14, {1, 19, 44, 204, 578}},
{15, {1, 9, 11, 90, 210, 2497}},
{16, {1, 410, 857, 2678}},
{17, {1, 23, 102, 760, 1612, 9115}},
{18, {1, 1034, 1695, 7856}},
{19, {1, 120, 326, 732, 3365, 3694, 3934, 23557, 28739}},
{20, {1, 688, 3103, 7340, 12083, 21401, 38745}}}
}


For b>=15 the results were obtained with a code used in the next paragraph.

Cycle parameters

A cycle has two parameters, the start index pos[[1]] and the end index pos[[2]], the difference being the length of the cycle. These depend on the number base b and the initial value s of the sequence.

The following code gives a graph of the dependence of the cycle parameters for number bases from 2 to 20 for a given initial value a1 of the sequence. The code can be used to experiment with varying initial values.

Remark: here I use the version with Tally[] which is much faster than the more elegant looking NestLWhileList[].

Here is the example for a[1] = 2

With[{s = 2, bmax = 20},
q2 = Table[
m = 1;
a = NestList[dsfc[#, b] &, s, m];
st = Select[Tally[a], #[[2]] > 1 &];
While[st == {},
m = 5*m;
a = NestList[dsfc[#, b] &, s, m];
st = Select[Tally[a], #[[2]] > 1 &]];
{b, pos = Flatten[Position[a, st[[1, 1]]]][[1 ;; 2]],
p = pos[[2]] - pos[[1]] }, {b, 2, bmax}];
q3 = #[[1 ;; 2]] & /@ q2;
ListLinePlot[Log[Transpose[#[[2]] & /@ (Join[{{1, {1, 1}}}, q3])]],
PlotLabel ->
"Interation of dsfc\nCycle data as a function of the number base b\
\nfor sequences with initial value a(1) = " <> ToString[s] <>
":" "\nblue curve = index of start of the cycle\nyellow curve = \
index of end of the cycle\nNotice: we have drawn curves rather than \
points for better visibility\n",
AxesLabel -> {"number base", "log(10,index k)"},
PlotRange -> {{2, bmax}, All}, ImageSize -> 500]]


And another one with a[1] = 5

Attempted proof of periodicity

The following proof is not strict but I believe the idea is correct. Critical comments are appreciated.

We shall sketch the proof of the claim that for any starting value the sequence of iterated application of dsfc eventually leads to a cycle, i.e. the sequence becomes periodic.

The elements of a sequence cannot grow indefinitely. This follows that there is only a finite set of numbers available for a sequence. Hence, by the drawer principle, at least two elements in a sequence must be equal. This leads to a cycle. QED.

Let us look at the first assertion in the case of base b = 10. The maximum value of an element with q decimal positions is m = 10^q - 1. The next element in the sequence is v = q 9^q. Now there is a number q1 of decimal positions for which v < u for q > q1. Numerically we find q1 = 34.

Notice that a trivial periodicity appears if some a[k] = b^k, especially the starting value. Then this number has only one digit 1, the rest are zero, and the period length is equal to 1.

• Possible exceptions to the proof could be numbers where dsfc[x] is equal to any multiple of 10, or for any number that has a multiple of 10 in it's following sequence. If you eventually got 10000 for example, it would immediately reduce to 1, and from there on it will then continue to be one. You may choose to define this also as periodicity, with a period of 1. Apr 27, 2017 at 23:09
• @volcanrb You are right. Thank you for pointing this out.I have modified the text to include this case explicitly. Apr 28, 2017 at 18:42

In:

dsf[n_] := With[{t = IntegerDigits[n]}, Sum[Max[1, k]^k, {k, t}]]
dsf[256]
RecurrenceTable[{c[n + 1] == dsf[c[n]], c[1] == 2}, c, {n, 1, 4}]


Out:

49785
{2, 4, 256, 49785}

• better use of Sum[], but why does that change how RecurrenceTable works? It should have produced a proper answer in the original poster's question too. Don't you think? Apr 25, 2017 at 22:28
• I've added use of your method to my answer to demonstrate how it can result in incredible function speedup! Thanks for that. Apr 25, 2017 at 22:40
• It's nothing wrong with RecurrenceTable. dsf is incorrect. If you use a symbol as input to your dsf function, such as dsf[n], the result would be Max[1, n]^n. When n = 256, Max[1,256]^256==323170060.... Apr 26, 2017 at 6:48
• I'm not following. @volcanrb is clearly setting t to IntegerDigits[]. Apr 26, 2017 at 6:56
• I know, but 256 is a number, but not a symbol. If you use a symbol but not a number, dsf will not work as expected. Apr 26, 2017 at 7:10

You're pointing out a very interesting oddity of RecurrenceTable[]. I do believe it's a bug. It shouldn't act this way.

As it has been pointed out, since you only want to obtain f[...f[f[x]]...], NestList will do just fine.

On another note, since you will be running this function many times, you should probably write it with a bit more optimization. You're using a VERY computationally expensive method for computing the Sum[].

Original method:

ClearAll[dsf];
dsf[n_Integer] := Block[{m = n, t}, t = IntegerDigits[m];
Sum[Max[1, t[[k]]]^t[[k]], {k, Length[t]}]]
dsf[#] & /@ Range@(10^5); // AbsoluteTiming


Out: {1.71049, Null}

Improving summation for the special case of a relatively short list of integers:

ClearAll[dsf];
dsf[n_Integer] := Total@Map[#^# &, IntegerDigits[n] /. {0 -> 1}]
dsf[#] & /@ Range@(10^5); // AbsoluteTiming


Out: {1.30601, Null}

Not bad, almost 25% improvement. But many numbers have same digits in a different order, and we can cache these results:

ClearAll[dsf];
dsf[n_List] := dsf[n] = Total@Map[#^# &, n];
dsf[n_Integer] := dsf[n] = dsf[Sort[IntegerDigits[n] /. {0 -> 1}]];
dsf[#] & /@ Range@(10^5); // AbsoluteTiming
dsf[#] & /@ Range@(10^5); // AbsoluteTiming


Out: {0.997551, Null} <--first pass is 42% faster {0.0793948, Null} <--next pass is LIGHTNING FAST!

Lastly, we can compile the function. This will kill our ability to cache at the function level (can be brought back with another function name), but the speed will go through the roof:

ClearAll[dsf];
dsf = Compile[{{t,_Integer}}, Total@Map[#^#&,IntegerDigits[t]/.{0->1}],
RuntimeOptions->"Speed", CompilationTarget->"C"];
dsf[#] & /@ Range@(10^5); // AbsoluteTiming
dsf[#] & /@ Range@(10^5); // AbsoluteTiming


Out: {0.494818, Null} <-- over 71% of execution time shaved off from original {0.495143, Null} <-- no caching means same speed every time

To bring back caching:

ClearAll[dsf, dsfInternal];
dsfInternal = Compile[{{n,_Integer,1}}, Total@Map[#^#&,n],
RuntimeOptions->"Speed", CompilationTarget->"C"];
dsf[n_Integer] := dsf[n] = dsfInternal[Sort[IntegerDigits[n]/.{0->1}]];
dsf[#] & /@ Range@(10^5); // AbsoluteTiming
dsf[#] & /@ Range@(10^5); // AbsoluteTiming


Out: {1.08752, Null} {0.0812414, Null}

First time is actually slower because of caching, but faster than without compilation. The next pass is the same speed as the first time we tried caching, because actual function is not executed anymore, it's just a hash table lookup.

EDIT: adding Sum[] version from UnchartedWorks

ClearAll[dsf];
dsf = Compile[{{n,_Integer}}, Sum[Max[1,k]^k, {k, IntegerDigits[n]}]],
RuntimeOptions->"Speed", CompilationTarget->"C"];
dsf[#] & /@ Range[10^5]; // AbsoluteTiming


Out: {0.0155958, Null}

Ultimately, we've achieved performance speed-up of 115X.

Enjoy!

Clear[dsf]

dsf[n_Integer] := Block[{m = IntegerDigits[n]},
Total[(Clip[#, {1, Infinity}] & /@ m)^m]]

seq = RecurrenceTable[{c[n] == dsf[c[n - 1]], c[1] == 2},
c, {n, 1, 11}]

(*  {2, 4, 256, 49785, 405024629, 387470792, 406668622, 16964105,
387517185, 35207797, 389894275}  *)