How do I define a function in mathemathica that uses a complicated recurrence relation?

Question

I am trying to write a little script to compute an integer sequence. The function I'm trying to write in code is the one int the blackboard, a(n). The problem is that I was expecting the function h(n) I defined in the script to give a number as a result, but it is giving something else: for h(2) it gives ArgMax[{p, HarmonicNumber[p] <= 1}, p, Integers] How can I correct that? (You have to understand that I am by no means a programmer, nor know much about mathematica. Thanks in advnace. The script I wrote is this one:

    h[n_] := (ArgMax[{p, 
      Sum[1/s, {s, 1 + Sum[h[k], {k, 1, (n - 1)}], p}] <= 1}, p, 
     Integers]) - Sum[h[k], {k, 1, (n - 1)}]; h[1] = 1;

a(n)=(maximum p such that the sum from s equals r to p is less or equal than one)-r+1, where r=1+the sum from k=1 to (n-1) of a(k), and a(1)=1

PD: Those letters that look like v's are r's. Sorry.

I think there's a Mathematica-specific Stackexchange forum so you might get more and better answers there. That said, you really need to try to state the problem so that others can understand it; a link to a messy image doesn't help others understand. Help others help you by writing out the problem statement. — Robert Dodier, Feb 12 '20 at 21:10
Okay. I've seen a major flaw in the script so I'll just rewrite the question. I haven't found that forum though. — b0nes, Feb 12 '20 at 21:31

score 1 · Answer 1 · answered Feb 13 '20 at 06:56

The main reason you are receiving this answer is that the Sum command in Mathematica makes sophisticated symbolic computations to rewrite the terms and expressions into classical mathematical functions, trying to make the least possible assumptions on everything.

Of course, if you sum inverses of integers it tries to go back to the harmonic number.

It is very often useful but could be at cross-purposes.

To avoid this, just replace

Sum[ , {var, start, end}]

by

Total[Table[ , {var, start, end}]]

It should give you what you want, without the added time of symbolic interpretation of the summand.

Anyway, there are major troubles in your code.

You should use memoization for h
You want to find the maximum p such that the sum of inverses starting at r is inferior or equal to 1. Your code is not going to do that. You cannot simply put inequalities into range arguments. You have to loop around your sum with control or use inverse functions.

It's interesting that `Maximize[{Sum[1/s <= 1, {s, 4, p}]}, Element[p, Integers]]` doesn't process. I'm intrigued by how an inverse function could be used. Possibly not this particular case. — Chris Degnen, Feb 13 '20 at 18:49

Chris Degnen · Accepted Answer · 2020-02-13T18:51:45.827

a[1] = 1;

a[n_] := Module[{sum = 0},
  r = 1 + Sum[a[k], {k, n - 1}];
  x = r;
  While[sum <= 1, sum += 1/x++];
  p = x - 2;
  p - r + 1]

Table[a[n], {n, 6}]

{1, 2, 6, 16, 43, 117}

The result for a[4] is 16 not 14.

To illustrate, when n = 4

r = 1 + Sum[a[k], {k, 4 - 1}]

  = 1 + a[1] + a[2] + a[3]       (* refer to established results for a[n] *)

  = 1 +  1  +  2  +  6  =  10

sum = 0;
x = r;
While[sum <= 1, sum += 1/x++];
p = x - 2;
p - r + 1

16

Or in another form

Total[Table[1/s, {s, 10, 25}]] <= 1   (* True *)

p - r + 1 = 25 - 10 + 1 = 16

Using memoisation, as mentioned by ogerard

Clear[a]

a[1] = 1;

a[n_] := a[n] = Module[{sum = 0},
   r = 1 + Sum[a[k], {k, n - 1}];
   x = r;
   While[sum <= 1, sum += 1/x++];
   p = x - 2;
   p - r + 1]

only reduces the following run time by 9 seconds

Timing[Table[a[n], {n, 14}]]

{40.8906, {1, 2, 6, 16, 43, 117, 318, 865, 2351, 6391, 17372, 47222, 128363, 348927}}

How do I define a function in mathemathica that uses a complicated recurrence relation?

2 Answers2