4
$\begingroup$

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[
  MapThread[
   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]"}}]]

enter image description here

enter image description here

$\endgroup$

2 Answers 2

10
$\begingroup$

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:

figure

$\endgroup$
3
  • 1
    $\begingroup$ wonderfully creative, J. M. ... can DiscretePlot also be used in some way to show the points along the curve? $\endgroup$
    – user42700
    Mar 29, 2020 at 15:25
  • $\begingroup$ Sure, just add the option setting PlotMarkers -> {"Point", Medium} to DiscretePlot[]. $\endgroup$ Mar 29, 2020 at 15:28
  • $\begingroup$ Thank you very much, JM, very clever coding!!!! $\endgroup$
    – user42700
    Mar 29, 2020 at 15:36
4
$\begingroup$
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} *}

For your specific parameters,

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]]}}]}]]

enter image description here

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}]]

enter image description here

$\endgroup$
6
  • $\begingroup$ 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? $\endgroup$
    – user42700
    Mar 29, 2020 at 18:23
  • $\begingroup$ sol[k_] := table[[k]] $\endgroup$
    – Bob Hanlon
    Mar 29, 2020 at 18:31
  • $\begingroup$ but, what if I started with the RecurrenceTable instead of RSolve, as RSolve doesn't always produce a solution? $\endgroup$
    – user42700
    Mar 29, 2020 at 18:47
  • $\begingroup$ Bob: Many thanks for the edit -- very helpful and I've learned much. Thanks for sharing your talents. $\endgroup$
    – user42700
    Mar 30, 2020 at 1:26
  • $\begingroup$ Corrected index in definitions of L[k] and M[k] $\endgroup$
    – Bob Hanlon
    Mar 30, 2020 at 2:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.