# RSolve does not evaluate this recursion with two boundary conditions

I am using RSolve to solve for a function defined recursively, with two boundary conditions:

• First boundary condition describes the relationship between $$f(1)$$ and $$f(0)$$
• Second boundary condition describes the relationship between $$f(S)$$ and $$f(S-1)$$

This way the function should be defined between $$0$$ and $$S$$, and for all $$0 it should follow the recursive relationship. Here is my code:

Clear[f, g, s, \[Sigma], \[Theta], \[Delta], \[Beta], S]
eq1 = f[s - 1]*\[Sigma]*\[Theta] +
f[s + 1]*\[Delta]*(s +
1) - (\[Sigma]*\[Theta] + \[Delta]*s + \[Beta])  f[s] == 0;
eq2 = f[1]*\[Delta] - f[0]*\[Sigma]*\[Theta] + \[Beta]*(1 - f[0]) == 0;
eq3 = f[S - 1]*\[Sigma]*\[Theta] - (\[Delta]*S + \[Beta]) f[S] == 0;
sols = RSolveValue[{eq1, eq2, eq3}, {f[s]}, s]
\$RecursionLimit = 100000;
newfn[s_] := Evaluate[sols];
c1 = Solve[Sum[newfn[s], {s, 0, Infinity}] == 1, {C[1]}][[1]];
newfn[s] /. c1


After solving for the function I ideally want to impose the condition that $$\sum_{s\geq 0} f(s) = 1$$ to pin down the constant $$c_1$$. However, when I run the code above

sols = RSolveValue[{eq1, eq2, eq3}, {f[s]}, s]


part does not run: outputting me the evaluation for the input.

Any help or tips would be greatly appreciated! Thanks.

This recursion can be solved as follows.

r = RSolve[{eq1, eq2}, f, s] // Flatten
(* {f -> DifferenceRoot[Function[{y, n},
{θ*σ*y[n] + (-β - δ - n*δ - θ*σ)*y[1 + n] + (2 + n)*δ*y[2 + n] == 0,
y[0] == C[1], y[1] == -(β/δ) - ((-β - θ*σ)*C[1])/δ}]]} *)


To determine the unknown parameter, C[1], specify S. For instance, with S = 5,

c = Solve[eq3 /. S -> 5 /. r, C[1]] // Flatten // Simplify
(* {C[1] -> (β^5 + 120*δ^5 + β^4*(15*δ + 4*θ*σ) + β^3*(85*δ^2 +
36*δ*θ*σ + 6*θ^2*σ^2) + β^2*(225*δ^3 + 104*δ^2*θ*σ + 27*δ*θ^2*σ^2 +
4*θ^3*σ^3) + β*(274*δ^4 + 96*δ^3*θ*σ + 27*δ^2*θ^2*σ^2 + 6*δ*θ^3*σ^3 + θ^4*σ^4))/
(β^5 + 120*δ^5 + 120*δ^4*θ*σ + 60*δ^3*θ^2*σ^2 + 20*δ^2*θ^3*σ^3 + 5*δ*θ^4*σ^4 +
θ^5*σ^5 + 5*β^4*(3*δ + θ*σ) + 5*β^3*(17*δ^2 + 10*δ*θ*σ + 2*θ^2*σ^2) +
5*β^2*(45*δ^3 + 35*δ^2*θ*σ + 12*δ*θ^2*σ^2 + 2*θ^3*σ^3) +
β*(274*δ^4 + 250*δ^3*θ*σ + 110*δ^2*θ^2*σ^2 + 30*δ*θ^3*σ^3 + 5*θ^4*σ^4))} *)


Finally, recompute f with C[1] evaluated.

sols = RSolve[{eq1, eq2, f[0] == C[1] /. c}, f, s] // Flatten


The result is too long to reproduce here in InputForm but looks like

Incidentally, renormalizing f[s], as suggested in the question, is not possible, because C[1] is determined by eq3 and thus is not a free parameter.

Some general tips:

1. RSolveValue is not managing to solve this equation in a meaningful way, so resolving that is your first problem.
2. I would try progressively to simplify the problem e.g. reduce the number of parameters in the equation (e.g. by merging parameters), give numerical values for these parameters. If Mathematica can solve a simplified version of your problem, you can perhaps build up the solution from there. If the simplified version is insoluble, you will need to find a different approach.
3. Recurrence relations are generally easier to solve if you specify initial conditions, rather than final conditions. Can you specify initial conditions, then find the initial conditions that give the desired final conditions?
4. It might be useful to attempt a numerical solution, to understand the behaviour. For example, you might discover that you have made a mistake in the formulation that is worth correcting before proceeding further.

Good luck