# How to plot a system of recurrence equations

I have the following code working well ... but I would then like to plot M[k] and L[k] on the same axes with k along the horizontal and M[k] and L[k] values along the vertical. The graph I'm looking to automatically generate from the code is below. I've tried various constructs of DiscretePlot and ListPlot but to no avail ... I'd appreciate suggestions!!

Clear[a, b, k, M, L, M0, L0, NN]
NN = 25; a = 2/10; b = 1/10; M0 = 30; L0 = 50;
N[TableForm[
Prepend, {RecurrenceTable[{M[k + 1] == (1 - a)*M[k] + a*L[k],
L[k + 1] == (1 - b)*L[k] + b*M[k], M[0] == M0, L[0] == L0}, {M,
L}, {k, 0, NN}], Range[0, NN]}],
TableHeadings -> {{}, {"k", "M[k]", "L[k]"}}]]


Straightforward method:

With[{NN = 25, a = 2/10, b = 1/10, M0 = 30, L0 = 50},
ListLinePlot[Transpose[RecurrenceTable[{M[k + 1] == (1 - a) M[k] + a L[k],
L[k + 1] == (1 - b) L[k] + b M[k],
M[0] == M0, L[0] == L0},
{M, L}, {k, 0, NN}]],
DataRange -> {0, NN}, PlotRange -> All]]


Slick method:

With[{NN = 25, a = 2/10, b = 1/10, M0 = 30, L0 = 50},
ListLinePlot[Transpose[NestList[{{1 - a, a}, {b, 1 - b}}.# &, {M0, L0}, NN]],
DataRange -> {0, NN}, PlotRange -> All]]


Even slicker method:

With[{NN = 25, a = 2/10, b = 1/10, M0 = 30, L0 = 50},
DiscretePlot[MatrixPower[{{1 - a, a}, {b, 1 - b}}, k, {M0, L0}] // Evaluate,
{k, 0, NN}, Filling -> None, Joined -> True,
PlotRange -> All]]


All three versions produce the following figure:

• wonderfully creative, J. M. ... can DiscretePlot also be used in some way to show the points along the curve?
– PRG
Mar 29 '20 at 15:25
• Sure, just add the option setting PlotMarkers -> {"Point", Medium} to DiscretePlot[]. Mar 29 '20 at 15:28
• Thank you very much, JM, very clever coding!!!!
– PRG
Mar 29 '20 at 15:36
Clear["Global*"]

eqns = {
M[k + 1] == (1 - a)*M[k] + a*L[k],
L[k + 1] == (1 - b)*L[k] + b*M[k],
M[0] == M0, L[0] == L0};


RSolve provides the exact solution to the recurrence equations

sol = RSolve[eqns, {L, M}, k][[1]]

{* {L -> Function[{k}, (a L0 + (1 - a - b)^k b L0 + b M0 - (1 - a - b)^k b M0)/(
a + b)], M ->
Function[{k}, -((-a L0 + a (1 - a - b)^k L0 - a (1 - a - b)^k M0 - b M0)/(
a + b))]} *}


Verifying,

eqns /. sol // Simplify

{* {True, True, True, True} *}


solEx[k_] = {L[k], M[k]} /. sol /.
{a -> 2/10, b -> 1/10, M0 -> 30, L0 -> 50} // Simplify

{* {1/3 10^(1 - k) (2 7^k + 13 10^k), 1/3 10^(1 - k) (-4 7^k + 13 10^k)} *}


The functions share a common limit

lim = Limit[solEx[k], k -> Infinity]

(* {130/3, 130/3} *)


Plotting

With[{NN = 25},
Plot[Evaluate@solEx[k], {k, 0, NN},
PlotRange -> All,
PlotLegends -> Placed[{"L", "M"}, {0.5, 0.3}],
Prolog -> {Gray, Dashed,
Line[{{0, lim[[1]]}, {NN, lim[[1]]}}]}]]


EDIT: If RSolve is unable to solve:

Clear["Global*"]

NN = 25; a = 2/10; b = 1/10; M0 = 30; L0 = 50;
rt = RecurrenceTable[{M[k + 1] == (1 - a)*M[k] + a*L[k],
L[k + 1] == (1 - b)*L[k] + b*M[k], M[0] == M0, L[0] == L0}, {M, L}, {k, 0,
NN}];

M[k_Integer] := rt[[k + 1, 1]] /; 0 <= k <= NN

L[k_Integer] := rt[[k + 1, 2]] /; 0 <= k <= NN

DiscretePlot[{L[k], M[k]}, {k, 0, NN},
PlotRange -> All, Filling -> None,
PlotLegends -> Placed["Expressions", {0.5, 0.3}]]


• Hi Bob: This is very nice and clever. The convergence line is a nice feature; when RSolve can't solve a system, I resort to the RecurrenceTable. How would sol = work when storing values from a recurrence table?
– PRG
Mar 29 '20 at 18:23
• sol[k_] := table[[k]] Mar 29 '20 at 18:31
• but, what if I started with the RecurrenceTable instead of RSolve, as RSolve doesn't always produce a solution?
– PRG
Mar 29 '20 at 18:47
• Bob: Many thanks for the edit -- very helpful and I've learned much. Thanks for sharing your talents.
– PRG
Mar 30 '20 at 1:26
• Corrected index in definitions of L[k] and M[k] Mar 30 '20 at 2:04