kill(all)$

maxima : casino  ( Why play or not: random walk),version 1.0,part 3
peter.vlasschaert@gmail.com, 03/02/2020

see :youtube :https://www.youtube.com/watch?v=9FVlkg2bLlI
                         =>   " Lecture 25 Random Walks "

load packages 
============

load(lrats)$
load(solve_rec)$ /*use solving Difference equations */

see :yt:15min28sec

claim : Mathematical Model (  = second order system )
=====================
We want get  T == n ( = dollars to start) + m( = dollars to win)
            /   0 if  n =0 ( Broke )               (1)     BC
X[n]  = |   1 if  n = T ( We get target)   (2)    BC
             \  p*X[n+1] +q*X[n-1]
constraint : p+q  = 1 => q = 1-p  'loose' , probability
                                         p            'win   '  , probability
X[n] = p*X[n+1] +(1-p)*X[n-1]              (3)   Difference Equations
-------------------------------------------------------------------------------------
W=X[n+1] = to win     one dollar with the game @ every bet
L = X[n-1]  = to loose one dollar with the game @ every bet
BC = Boundary Condition
-------------------------------------------------------------------------------------
Play : gives sequences 
WLWWLLLW.........

general : model ' Homogene Second Order Difference equation '
=================================================
1e)  solution  'r' : r1 ≠ r2
2e)  solution  'r' : r1 = r2 =r

pa:a2*a[n+2]+a1*a[n+1]+a0*a[n] = 0;

see  (3), model mutual exclusive.' see'

pb: X[n] = p*X[n+1]+(1-p)*X[n-1];

1e)  solution  'r' : r1 ≠ r2 
---------------------------------

p1:rhs(pb) -lhs(pb) = 0;

to find : form 'pa'
-----------------------
=>  replace n with n+1 into p1

p1a:subst([n=n+1],p1);

  =>  a2,a1,a0

p1b1:a2 = ratcoef(lhs(p1a),X[n+2]);
p1b2:a1 = ratcoef(lhs(p1a),X[n+1]);
p1b3:a0 = ratcoef(lhs(p1a),X[n]);

find solution : pa
-----------------------

p1c:[p1b1,p1b2,p1b3];
p1c1:subst(p1c,pa);

p1c2:solve_rec(p1c1,a[n]);

If  p ≠ 1/2  then 
===========

p1d: subst([a[n] =X[n]],p1c2);

use : BC  (1),(2)

p1d1:subst([n=0],rhs(p1d))=0;
p1d2:subst([n=T],rhs(p1d)) =1;

how to find %k[1],%k[2]

p1d3:solve([p1d1,p1d2],[%k[1],%k[2]]);

p1d4:flatten(p1d3);

find the solution : X[n]
use :p1d

p1d5:subst([p1d4],p1d);

rem : bet dollar a time
-------------------------------

p1d6:ratsimp(p1d5);

roulette: bet on black or bet on black
36 numbers  = 18 black + 18 white
1    number   0 :   green zero  
1     number  00 : green double  zero
============================
total : 38

result :probability red  'or'  black
=======================

p1d7a: p= 18/38,float;

rem :almost fair game ( = fair game p=1/2)

We want play the game : n = 1000$
                 when win           m = 100$ 
                                            T = n+m  we go home 'Happy'

p1d7b:n=1000;

p1d7c:T =1000+100;

p1d7d:m = 100;

p1d7:[p1d7a,p1d7b,p1d7c];

solve : p1d6

p1d8:lhs(p1d6) = subst(p1d7,rhs(p1d6));

see :video 

value :

1/37648,float;

limit calculation  ( n => ∞)
============

p1e:subst([T=n+m],p1d6);

assume(p<1/2,-p<0);

declare(m,integer)$

p1e1:ev(limit(p1e,n,inf),nouns);

p1e2:subst([p1d7a,p1d7d],p1e1);

If  p  = 1/2  then   (fair game)
======================
use : p1a, p=1/2

p2:subst([p=1/2],p1a);

solution : r1=r2 =r =1

p2a:solve_rec(p2,X[n]);

use : BC  (1),(2)

p2a1:subst([n=0],rhs(p2a))=0;
p2a2:subst([n=T],rhs(p2a)) =1;

how to find %k[1],%k[2]

p2a3:solve([p2a1,p2a2],[%k[1],%k[2]]);
p2a4:flatten(p2a3);

find the solution : X[n]
use :p2a

p2a5:subst(p2a4,p2a);

result :

p2a6:lhs(p2a5) =subst([p1d7b,p1d7c],rhs(p2a5));

p2a7:float(p2a6);

rem : r1 ≠ r2  ,different root  ( not fair game),Difference  equation
          r1 = r2  ,double   root  ( fair game       ),Difference equation

other topics : our in the video.