Michael Chuah's page for MAS.864 The Nature of Mathematical Modelling

mcx "at" mit "dot" edu

Home

Random Systems

Scanned work

Some MATLAB code used

6-4(d) Positon of a random walker in 1D
%% Meng Yee CHUAH
% Week 3 - Random Systems
close all; clear all; clc;

%% 6.2(c)
n = 1e5;
x = zeros(2,n);
y = zeros(2,n);
for i = 1:n
    x(:,i) = mod(8121*rand(2,1)+28411,134456)/134456;     % Initial seed generated by MATLAB (irony)
    y(:,i) = [sqrt(-2*log(x(1,i)))*sin(x(2,i)); sqrt(-2*log(x(1,i)))*cos(x(2,i))];
end

C1 = mean(y,2)
C2 = mean(y.^2,2)-mean(y,2).^2
C3 = mean(y.^3,2)-3.*C1.*C2-C1.^3



%% 6.4(d)
n = 1000;
x = 2^13-1; % Use 13 bits to store the information for an order 12 maximal LFSR
% dec2bin(x)
w = zeros(n,10);
for j = 1:10
    wn = 0; % Initial position of random walker
    for i = 1:n
        xn = mod(bitget(x,2)+bitget(x,5)+bitget(x,7)+bitget(x,13),2);
        x = bitshift(x,1,13);
        x = bitset(x,1,xn);
    %     dec2bin(x)
        if xn==1
            wn = wn+1;
        else
            wn = wn-1;
        end
        w(i,j) = wn;
    end
end

figure('Units','inches','Position',[0 0 6 4.5]); 
e = 1.5*std(w,1,2);
h = errorbar(zeros(1,n),e,'.','Color',0.9*ones(1,3));

hh = get(h,'children');    % Retrieve info from errorbar plot
xh = get(hh(2),'xdata');    % Get xdata from errorbar plot
wh = 2;
xh(4:9:end) = xh(1:9:end)-wh/2;  % Change xdata with respect to ratio
xh(7:9:end) = xh(1:9:end)-wh/2;
xh(5:9:end) = xh(1:9:end)+wh/2;
xh(8:9:end) = xh(1:9:end)+wh/2;
set(hh(2),'xdata',xh(:))  % Change error bars on the figure

hold on;
plot(w);
xlim([0 n]); title('6.4(d) Positon of a random walker in 1D');
xlabel('Time'); ylabel('Position');
saveas(gcf,'6-4(d) Positon of a random walker in 1D.png');
saveas(gcf,'6-4(d) Positon of a random walker in 1D.fig');

MATLAB code output

C1 =

    0.4034
    1.6376


C2 =

    0.0009
    0.0018


C3 =

   1.0e-05 *

    0.1339
    0.3859

>>