# RSolve fails on a two-equation non-linear system

I have a system of just two recurrence equations, I need to solve them and so I use RSolve:

RSolve[{H[n + 1] == (1 + (7 h[n])/10) H[n],
h[n + 1] == 1 - (1 + (7 h[n])/10) H[n]},
{H[n], h[n]}, n]


But this just returns my code, without any errors:

RSolve[{H[1 + n] == (1 + (7 h[n])/10) H[n],
h[1 + n] == 1 - (1 + (7 h[n])/10) H[n]},
{H[n], h[n]}, n]


Where am I going wrong? Why am I not even seeing any errors?

• I am modelling the frequency of two variants, H and h, through discrete time (n)
• Total frequency in the population: h + H = 1
• The two variant's frequencies depend on each other as you can see in the above code
• rewriting your equation shows that there are no Solutions. RSolve[{H[n + 1] == (1 + (7 h[n])/10) H[n], h[n + 1] + H[n + 1] == 1}, {H[n], h[n]}, n] gives {} – Max1 Apr 17 '14 at 20:09
• Thanks @Max1 hmm... perhaps I could just use RecurrenceTable and get some numerical "solution" – hello_there_andy Apr 17 '14 at 20:41

I know it doesn't answer your question directly, but with the definition above, your function H[n] asymptotically approaches 1 pretty fast regardless of the starting value H[0]:

(* Generates a table for the first 20 values of H given H[0] == alpha *)
f[alpha_] := Module[{H},
H[0] = alpha;
H[n_] := H[n] = H[n - 1] (1 + .7 (1 - H[n - 1]));
Table[H[n], {n, 20}]]

ListLinePlot[Table[f[x], {x, 0, 1, .05}]]


produces

You could also just program iteratively:

f[{x_, y_}] := With[{ch = {1, 0.7 x}.{x, y}}, {ch, 1 - ch}]
hdt[p_, n_] := Transpose@NestList[f, {p, 1 - p}, n]


If you just want $\{H(n),h(n)\}$ for starting values $\{H(0),h(0)\}=\{p,1-p\}$:

hd[p_, n_] := Nest[f, {p, 1 - p}, n];


Visualizing:

lp[p_] :=
ListPlot[hdt[p, 10], Joined -> True, PlotMarkers -> {Automatic, 8},
PlotLegends -> (Placed[
SwatchLegend[Automatic, {"H(t)", "h(t)"},
LegendLayout -> (Framed@Grid[#] &)], Scaled[{0.5, 0.5}]])];
Manipulate[lp[p], {p, 0.1, 0.9, Appearance -> "Labeled"}]


You could also do something like:

GraphicsRow[ListLinePlot@Transpose@RecurrenceTable[{
H[n + 1] == N@(1 + (7 h[n])/10) H[n],
h[n + 1] == N@1 - (1 + (7 h[n])/10) H[n], H[0] == #/100,
h[0] == 1/2},
{H[n], h[n]}, {n, 1, 15}] & /@ {179, 180}]


• Note that H + h == 1 as per the question, but you have H[0]+h[0] != 1, and also since these are frequencies, H <= 1. – Victor K. Apr 18 '14 at 3:59

In order to make my answer more comprehensible I decoupled the set of equations first:

(*1*) H[n+1]==7/10*H[n](17/7-H[n])
(*2*) h[n+1]==7/10*h[n](3/7+h[n])


with the constraints:

(h[n]+H[n]==1) && (0<=H[n]<=1) && (0<=h[n]<=1)


These can be obtained by using h[n]+H[n]==1. The recurence relations can now be solved independently. In this form it becomes clear, that the recurrence relations 1 and 2 are very similar to the logistic map.RSolve can find solutions of the logistic map for cases where an analytic solution exists (r=2 and r=4).

RSolve[a[n+1]==r(1-a[n])a[n]/.r->2,a[n],n]


Out:

{{a[n] -> 1/2 - 1/2 E^(2^n C[1])}}


Since there is no analytical solution available for aribtary r, Mathematica just gives back the expression (can be tested with r->3 for example). The case of H[n+1]==7/10(17/7-H[n])H[n] is even more general, therefore there is probably no known analytic solution which could be given. Thus mathematica just returns the expression.

Regarding the fixed points of the system: If one views the recurrence relations as functions and plots them,

Plot[{7/10 x (3/7 + x), 7/10 x (17/7 - x), x}, {x, 0, 1},PlotLegends->"Expressions"]


one can find the fixed points as the intersections with f[x_]=x. For any 1>H[0]>0 H[n] will converge towards 1 and h[n] will converge towards 0. However H[n]==0 and h[n]==1 are also intersecting points obeying the constraints and thus a solution. With sligthly different parameters there will probably occure more complex behavior similar to the logistic map.

In case anyone spots an error please comment.