Converging Sums of a Fractal Sequence

Question

Background

A fractal sequence is an integer sequences where you can remove the first occurrence of every integer and end up with the same sequence as before.

A very simple such sequence is called Kimberling's paraphrases. You start with the positive natural numbers:

1, 2, 3, 4, 5, 6, 7, 8, 9, ...

Then you riffle in some blanks:

1, _, 2, _, 3, _, 4, _, 5, _, 6, _, 7, _, 8, _, 9, ...

And then you repeatedly fill in the blanks with the sequence itself (including the blanks):

1, 1, 2, _, 3, 2, 4, _, 5, 3, 6, _, 7, 4, 8, _, 9, ...
1, 1, 2, 1, 3, 2, 4, _, 5, 3, 6, 2, 7, 4, 8, _, 9, ...
1, 1, 2, 1, 3, 2, 4, 1, 5, 3, 6, 2, 7, 4, 8, _, 9, ...
1, 1, 2, 1, 3, 2, 4, 1, 5, 3, 6, 2, 7, 4, 8, 1, 9, ...

That's our fractal sequence! Now let's take the partial sums:

1, 2, 4, 5, 8, 10, 14, 15, 20, 23, 29, 31, 38, 42, 50, 51, 60, ...

But what if we iterate this process? "Fractalise" the new sequence (i.e. the partial sums obtained from the steps above):

1, _, 2, _, 4, _, 5, _, 8, _, 10, _, 14, _, 15, _, 20, _, 23, ...
1, 1, 2, _, 4, 2, 5, _, 8, 4, 10, _, 14, 5, 15, _, 20, 8, 23, ...
1, 1, 2, 1, 4, 2, 5, _, 8, 4, 10, 2, 14, 5, 15, _, 20, 8, 23, ...
1, 1, 2, 1, 4, 2, 5, 1, 8, 4, 10, 2, 14, 5, 15, _, 20, 8, 23, ...
1, 1, 2, 1, 4, 2, 5, 1, 8, 4, 10, 2, 14, 5, 15, 1, 20, 8, 23, ...

And take the partial sums again:

1, 2, 4, 5, 9, 11, 16, 17, 25, 29, 39, 41, 55, 60, 75, 76, 96, ...

Rinse, repeat. It turns out that this process converges. Every time you repeat this process, a larger prefix of the sequence will remain fixed. After an infinite amount of iterations, you would end up with OEIS A085765.

Fun fact: This process would converge to the same sequence even if we didn't start from the natural numbers as long as the original sequence starts with 1. If the original sequence starts with any other x, we'd obtain x*A085765 instead.

The Challenge

Given a positive integer N, output the Nth element of the converged sequence.

You may write a program or function, taking input via STDIN (or closest alternative), command-line argument or function argument and outputting the result via STDOUT (or closest alternative), function return value or function (out) parameter.

You may choose whether the index N is 0- or 1-based.

Test Cases

The sequence starts with:

1, 2, 4, 5, 9, 11, 16, 17, 26, 30, 41, 43, 59, 64, 81, 82, 108, 117, 147, 151, 192, 203, 246, 248, 307, 323, 387, 392, 473, 490, 572, 573, 681, 707, 824, 833, 980, 1010, 1161, 1165, 1357, 1398, 1601, 1612, 1858, 1901, 2149, 2151, 2458, 2517

So input 5 should result in output 9.

Here is a naive CJam reference implementation which generates the first N numbers (given N on STDIN). Note that your code should only return the Nth element, not the entire prefix.

Answer 1

CJam (23 22 bytes)

The partial sums are given at the even indexes of the fractal sequence, which is A086450. The recurrence given there as the definition of A086450 is the basis for these implementations.

Using an explicit "stack" (in scare quotes because it's not LIFO):

{),){2md~)\),>+$)}h+,}

Online demo

Dissection

{         e# Anonymous function body; for clarify, pretend it's f(x)
          e# We use a stack [x_0 ... x_i] with invariant: the result is sum_j f(x_j)
  ),      e# Initialise the stack to [0 ... x]
  )       e# Uncons x, because our loop wants one value outside the stack
  {       e# Loop. Stack holds [x_0 ... x_{i-1}] x_i
    2md   e# Split x_i into (x_i)/2 and (x_i)%2
    ~)\   e# Negate (x_i)%2 and flip under (x_i)/2
    ),>   e# If x_i was even, stack now holds [x_0 ... x_{i-1}] [0 1 ... (x_i)/2]
          e# If x_i was odd, stack now holds [x_0 ... x_{i-1}] [(x_i)/2]
    +     e# Append the two arrays
    $     e# Sort to get the new stack
    )     e# Uncons the greatest element in the new stack
  }h      e# If it is non-zero, loop
          e# We now have a stack of zeroes and a loose zero
  +,      e# Count the total number of zeroes, which is equivalent to sum_j f(0)
}

---

At 23 bytes there's a much more efficient approach, with memoisation:

{2*1a{2md~)\){j}%>:+}j}

Online demo

Answer 2

Python 2, 55 49 42

I have no idea what's going on, but it seems hard to beat the Maple formula from the OEIS page. This uses 0-based indexing.

f=lambda n,t=0:n<1or f(n/2,n%2)-~-t*f(n-1)

Thanks to @PeterTaylor for -6 bytes.

Answer 3

Haskell, 65

s l=[0..]>>=(\i->[l!!i,s l!!i])
r=1:(tail$scanl1(+)$s r)
f n=r!!n

Answer 4

Templates Considered Harmful, 124

Fun<If<A<1>,Add<Ap<Fun<Ap<If<Sub<A<1>,Mul<I<2>,Div<A<1>,I<2>>>>,A<0>,A<0,1>>,Div<A<1>,I<2>>>>,A<1>>,Ap<A<0>,Sub<A<1>,T>>>,T>>

This is an anonymous function. It's more or less the same as my Python answer the Maple formula on the OEIS page, except that I didn't implement modulus, so I had to use n-n/2*2 instead of n%2.

Expanded:

Fun<If<
    A<1>,
    Add<
        Ap<
            Fun<Ap<
                If<
                    Sub<
                        A<1>,
                        Mul<
                            I<2>,
                            Div<A<1>,I<2> >
                        >
                    >,
                    A<0>,
                    A<0,1>
                >,
                Div<A<1>,I<2>>
            >>,
            A<1>
        >,
        Ap<
            A<0>,
            Sub<A<1>, T>
        >
    >,
    T
>> 

Answer 5

Mathematica, 47 44 bytes

If[#<1,1,#0[Floor@#/2]+(1-2#~Mod~1)#0[#-1]]&

Answer 6

Matlab 108 103

I am using the fact that the desired series is the partial sum of https://oeis.org/A086450

But the computation complexity of my implementation is far from optimal, even for this simple recurrence.

n=input('')+1;
z=zeros(1,n);z(1)=1;
for k=1:n;
z(2*k)=z(k);
z(2*k+1)=sum(z(1:k+1));
end;
disp(sum(z(1:n)))