% Simple OU  %%%%%%%%%%%%%%%%%%%
% Approximate the solution to
% dx = 0.25 dt + \sqrt(x) dW
% x(0) = 1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear, close all
Nmax = 1000; % number of ensemble members
N = 1e3;% length of time series
x = zeros(N+1,Nmax); % time series
dt = 0.001; %time step
Tf = N*dt;
tt = linspace(0,Tf, N+1); % time vector
for jj = 1:Nmax
    dW = sqrt(dt)*randn(N,1);
    x(1,jj) = 1;
    for ii = 1:N
       x(ii+1,jj) = x(ii,jj) + 0.25*dt + x(ii,jj)*dW(ii,1);
    end
end
avgx=mean(x,2);

plot(tt,avgx)

