For the recurrence relation:
$f(0)=1$
$f(1)=1$
$f(2)=2$
$f(2n)=f(n)+f(n+1)+n$
$f(2n+1)=f(n)+f(n−1)+1$
How to find square matrices $M_0, M_1$ and vectors $u, v$ such that if the base-2 expansion of $n$ is given by $e_1 e_2 \cdots e_j$, then $$f(n) = u M_{e_1} \cdots M_{e_j} v.$$ ??
This is how I would approach this question, and I guess my answer would not be unique.
For each $n = e_1e_2\cdots e_j$, there are $7$ numbers I should store (as a row vector) to calculate $f(2n)$ and $f(2n+1)$: $$uM_{e_1}M_{e_2}\cdots M_{e_j} = w(n) = \pmatrix {1&n&f(n-2)&f(n-1)&f(n)&f(n+1)&f(n+2)}.$$
To calculate $f(2n)$, the matrix $M_0$ is appended to the above product: $$w(2n) = \pmatrix{1\\2n\\f(2n-2)\\f(2n-1)\\f(2n)\\f(2n+1)\\f(2n+2)}^T =\pmatrix{1\\2n\\f(n-1)+f(n)+n-1\\f(n-1)+f(n-2)+1\\f(n)+f(n+1)+n\\f(n)+f(n-1)+1\\f(n+1)+f(n+2)+n+1}^T = w(n)M_0$$
To calculate $f(2n+1)$, the matrix $M_1$ is appended to the above product: $$ w(2n+1) = \pmatrix{1\\2n+1\\f(2n-1)\\f(2n)\\f(2n+1)\\f(2n+2)\\f(2n+3)}^T =\pmatrix{1\\2n+1\\f(n-1)+f(n-2)+1\\f(n)+f(n+1)+n\\f(n)+f(n-1)+1\\f(n+1)+f(n+2)+n+1\\f(n+1)+f(n)+1}^T = w(n)M_1$$
So $M_0$ can be written as
$$M_0 = \pmatrix{1&0&-1&1&0&1&1\\ 0&2&1&0&1&0&1\\ 0&0&0&1&0&0&0\\ 0&0&1&1&0&1&0\\ 0&0&1&0&1&1&0\\ 0&0&0&0&1&0&1\\ 0&0&0&0&0&0&1}$$
And similarly for $M_1$. I imagine the base case to be
$$w(2) = uM_1M_0 = \pmatrix{1&2&1&1&2&3&7}$$
For the rest, I have not tried out the calculation yet, but I suppose this involves solving the equation $$uM_1M_0 = w(2)$$ to solve for $uM_1$ and then $u$, and verify if that satisfies $w(0)$, $w(1)$ and $w(3)$. Lastly
$$v = \pmatrix{0\\0\\0\\0\\1\\0\\0}$$
Edit: Just found out $uM_1M_0 = w(2)$ is inconsistent and there is no $uM_1$ which satisfy the equation. So hope this answer could give you some idea on other ways to solve this.