RSolve with integer and boolean variable

Question

I have a recurrence relation with two variables $i, b$ where $i \in \{0, \ldots, n-1\}$ for some parameter $n \in \mathbb{N}$; and $b\in \{0,1\}$. I tried to solve it using RSolve, but the function just outputs my command (RSolve[...]) without giving any results. What am I doing wrong?

The recurrence is:

eqn[i_, b_] = {
f[i, b] == (1 - b) ((i r f[-1 + i, b])/(1/n + i r) + ( i r f[i, 1 - b])/(1/n + i r)) + 
                 b (f[i, 1 - b]/(1 + r/(-1 + n)) + (r f[1 + i, b])/((-1 + n) (1 + r/(-1 + n)))), 
f[0, 1] == (r f[1, 1])/(-1 + n + r), 
f[n-1, 0] == ((n-1) r)/(1 + (n-1) r) + f[n-2, 0]/( 1 + (-1 + n) r)
}

where $r>1$ is just some parameter. I tried

RSolve[eqn[i, b], f[i, b], {i, b}]

Which outputs just

RSolve[{f[i, 
    b] == (1 - b) ((i r f[-1 + i, b])/(1/n + i r) + (
       i r f[i, 1 - b])/(1/n + i r)) + 
    b (f[i, 1 - b]/(1 + r/(-1 + n)) + (
       r f[1 + i, b])/((-1 + n) (1 + r/(-1 + n)))), 
  f[0, 1] == (r f[1, 1])/(-1 + n + r), 
  f[-1 + n, 0] == ((-1 + n) r)/(1 + (-1 + n) r) + f[-2 + n, 0]/(
    1 + (-1 + n) r)}, f[i, b], {i, b}]

i.e. the output is just the input we gave to mathematica.

Neither the comments here [1] nor the proposed solution here [2] helped me. Does anyone spot my the mistake I am making? In general, this should be a solvable recurrence.

Answer 1

With the boundary conditions omitted, as permitted by @black, RSolve can solve this problem.

nobc = Replace[f[i, b] == (1 - b) ((i r f[-1 + i, b])/(1/n + i r) + 
    (i r f[i, 1 - b])/(1/n + i r)) + b (f[i, 1 - b]/(1 + r/(-1 + n)) + 
    (r f[1 + i, b])/((-1 + n) (1 + r/(-1 + n)))), {{b -> 0}, {b -> 1}}, 
    Infinity] /. {f[z_, 0] -> h[z], f[z_, 1] -> g[z]}

(* {h[i] == (i r g[i])/(1/n + i r) + (i r h[-1 + i])/(1/n + i r), 
    g[i] == (r g[1 + i])/((-1 + n) (1 + r/(-1 + n))) + h[i]/(1 + r/(-1 + n))} *)

where I have used @Bill's substitution, {f[z_, 0] -> h[z], f[z_, 1] -> g[z]} for notational simplicity. RSolve immediately yields the solution

snobc = RSolve[nobc, {g[i], h[i]}, i] // Flatten

However, because this solution contains many instances of two distinct DifferenceRoots, displaying it as produced by Mathematica is quite cumbersome. Instead, do the substitution,

snobc /. Thread[(Cases[snobc, _DifferenceRoot, Infinity, Heads -> True] //
    Union) -> {dr1, dr2}]

(* {g[i] -> C[2] dr1[i] + C[1] dr2[i], 
    h[i] -> C[2] (dr1[i] - (r (-dr1[i] + dr1[1 + i]))/(-1 + n)) + 
        C[1] (dr2[i] - (r (-dr2[i] + dr2[1 + i]))/(-1 + n))} *)

where dr1 and dr2 represent respectively

{DifferenceRoot[Function[{y, m}, {-((1 + m)*n*r*(-1 + n + r)*y[m]) + 
     (-1 + n + r + 2*n*r^2 + 2*m*n*r^2)*y[1 + m] - 
     r*(1 + n*r + m*n*r)*y[2 + m] == 0, y[0] == 0, y[1] == 1}]], 
 DifferenceRoot[Function[{y, m}, {-((1 + m)*n*r*(-1 + n + r)*y[m]) + 
     (-1 + n + r + 2*n*r^2 + 2*m*n*r^2)*y[1 + m] - 
     r*(1 + n*r + m*n*r)*y[2 + m] == 0, y[0] == 1, y[1] == 0}]]}

The two boundary conditions can then be used to determine the constants C[1] and C[2], at least in principle.

Addendum: C[1] and C[2] computation

The i = 0 boundary condition can be applied without difficulty.

f[0, 1] == (r f[1, 1])/(-1 + n + r) /. f[z_, 1] -> g[z]
% /. g -> Function[i, C[2] (dr1[i] + (r dr2[i])/(-1 + n + r))] /. 
   Equal -> Rule;
snobc1 = Simplify[snobc /. %]

(* g[0] == (r g[1])/(-1 + n + r) *)
(* C[1] -> (r C[2])/(-1 + n + r) *)
(* {g[i] -> C[2] (dr1[i] + (r dr2[i])/(-1 + n + r)), 
    h[i] -> C[2] (dr1[i] + (r (dr1[i] - dr1[1 + i]))/(-1 + n) + 
    (r ((-1 + n + r) dr2[i] - r dr2[1 + i]))/((-1 + n) (-1 + n + r)))} *)

Note that I continue to use dr1 and dr2 to represent the two DifferenceRoots, which would replace dr1 and dr2 in the actual computations.

Determining C[2] is much more cumbersome.

f[n - 1, 0] == ((n - 1) r)/(1 + (n - 1) r) + 
    f[n - 2, 0]/(1 + (-1 + n) r) /. f[z_, 0] -> h[z]
% /. h -> Function[i, C[2] (dr1[i] + (r (dr1[i] - dr1[1 + i]))/(-1 + n) + 
         (r ((-1 + n + r) dr2[i] - r dr2[1 + i]))/((-1 + n) (-1 + n + r)))];
Simplify[Solve[%, C[2]][[1, 1]]]
snobc2 = Simplify[snobc1 /. %]

(* h[-1 + n] == ((-1 + n) r)/(1 + (-1 + n) r) + h[-2 + n]/(1 + (-1 + n) r) *)
(* C[2] -> ((-1 + n) r)/((1 + (-1 + n) r) (dr1[-1 + n] + (r (dr1[-1 + n] - 
   dr1[n]))/(-1 + n) - (dr1[-2 + n] + (r (dr1[-2 + n] - 
   dr1[-1 + n]))/(-1 + n) + (r ((-1 + n + r) dr2[-2 + n] - 
   r dr2[-1 + n]))/((-1 + n) (-1 + n + r)))/(1 + (-1 + n) r) + 
   (r ((-1 + n + r) dr2[-1 + n] - r dr2[n]))/((-1 + n) (-1 + n + r)))) *)
(* {g[i] -> ((-1 + n) r (dr1[i] + 
   (r dr2[i])/(-1 + n + r)))/((1 + (-1 + n) r) (dr1[-1 + n] + 
   (r (dr1[-1 + n] - dr1[n]))/(-1 + n) - (dr1[-2 + n] + 
   (r (dr1[-2 + n] - dr1[-1 + n]))/(-1 + n) + 
   (r ((-1 + n + r) dr2[-2 + n] - 
   r dr2[-1 + n]))/((-1 + n) (-1 + n + r)))/(1 + (-1 + n) r) + 
   (r ((-1 + n + r) dr2[-1 + n] - r dr2[n]))/((-1 + n) (-1 + n + r)))), 
    h[i] -> ((-1 + n) r (dr1[i] + (r (dr1[i] - dr1[1 + i]))/(-1 + n) + 
   (r ((-1 + n + r) dr2[i] - r dr2[1 + i]))/
   ((-1 + n) (-1 + n + r))))/((1 + (-1 + n) r) (dr1[-1 + n] + 
   (r (dr1[-1 + n] - dr1[n]))/(-1 + n) - (dr1[-2 + n] + 
   (r (dr1[-2 + n] - dr1[-1 + n]))/(-1 + n) + 
   (r ((-1 + n + r) dr2[-2 + n] - r dr2[-1 + n]))/
   ((-1 + n) (-1 + n + r)))/(1 + (-1 + n) r) + 
   (r ((-1 + n + r) dr2[-1 + n] - r dr2[n]))/((-1 + n) (-1 + n + r))))} *)

Here is a numerical result for specific parameters.

Table[Simplify[snobc2 /. {n -> 5, r -> 2} /. i -> ii], {ii, 0, 4}]

(* {{g[0] -> -(391468/3850125), h[0] -> 0}, 
    {g[1] -> -(391468/1283375), h[1] -> -(71176/256675)}, 
    {g[2] -> -(462644/1283375), h[2] -> -(467728/770025)}, 
    {g[3] -> 513484/3850125, h[3] -> -(117752/256675)}, 
    {g[4] -> 1691004/1283375, h[4] -> 215072/256675}} *)