## In this set of demo's we are going to study
# random walks on the integers.
## Uniform random walks
# In this problem, we can move around the integers.
# At each step, we move -1 or +1 with prob 1/2 each.
# We start at the point 0 and stop when
# we get to +6 or - 4.
# How long does this take to stop?
# In this problem, we'll show that the answer is 24.
# A simulation is easy to code in Julia
function simlength()
S = 0
len = 0
while S != -4 && S != 6
len += 1
S += rand(-1:2:1)
end
return len
end
## Now we need to run this a large number of times.
# in-class: write the code to simulate this
# in-class: typo on function to show scoping rules
N = 10^6
avg = 0
for i=1:N
avg += simlength()/N
end
avg
# this fails because of Julia's weird scoping rules with Global :(
##
avg = 5
for i=1:10
avg = 6
end
## You get the expected behavior in a function
function mysim(N)
avg = 0
for i=1:N
avg += simlength()/N
end
return avg
end
mysim(10^6)
## Here's a fancy Julia way.
# run simlength 10^6 times
N = 10^6
expectedlen = mapreduce(x -> simlength()/N, +, 1:N)
## Now, let's derive the following system.
using LinearAlgebra
A = collect(SymTridiagonal(ones(11), -0.5*ones(10)))
##
A[1,2] = 0
A[end,end-1] = 0
b = [0;ones(9);0]
##
x = A\b
## Convert A as a sparse matrix.
using SparseArrays
B = sparse(A)
## Create A as a sparse matrix.
# we list all the non-zero entries in A.
function build_A()
nz = 9*3 + 2
I = zeros(Int, nz)
J = zeros(Int, nz)
V = zeros(nz)
# the arrays I,J,V satisfy:
# A[I[i],J[i]] = V[i]
index = 1
I[index] = 1
J[index] = 1
V[index] = 1.0
index += 1
for i=1:9
I[index] = i+1
J[index] = i
V[index] = -0.5
index += 1
I[index] = i+1
J[index] = i+1
V[index] = 1.0
index += 1
I[index] = i+1
J[index] = i+2
V[index] = -0.5
index += 1
end
I[index] = 11
J[index] = 11
V[index] = 1.0
index += 1
return sparse(I,J,V,11,11)
end
C = build_A()