%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%  Brownian Motion Drift Density
%
%  This code plots the exact density of Brownian Motion with Drift
%
%  Inputs:
%  1. n=number of discretization points
%  2. T=maximum time
%  3. mu= mean of the process
%  4. sigma= standard deviation of each jump
%  5. L= length of spatial domain
%
%  Ouputs:
%   1. Contour plot of the time varying density
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% Inputs
n=100;
T=1;
mu=3;
sigma=1;
L=3;

%% Discretization
t=linspace(.01,T,n);
x=linspace(-L,L,n);

[X,T]=meshgrid(x,t);
%% Movie of Process

figure('units','normalized','outerposition',[0 0 1 1]);
xlabel('$$x$$','interpreter','latex')
ylabel('$$t$$','interpreter','latex')
ax=gca;
ax.FontSize = 40;

hold on;

for i=1:n
    clf;
    density=1./sqrt(2*pi*sigma^2*t(i)).*exp(-(x-mu*t(i)).^2./(2*sigma^2*t(i)));
    plot(x,density,'linewidth',2);
    xlabel('$$x$$','interpreter','latex')
    ylabel('$$X_t$$','interpreter','latex')
    ax=gca;
    ax.FontSize = 40;
    axis([-L,L,0,2]);
    drawnow;

    
end

%% Contour Plot
figure('units','normalized','outerposition',[0 0 1 1]);
hold on;
xlabel('$$x$$','interpreter','latex')
ylabel('$$t$$','interpreter','latex')
ax=gca;
ax.FontSize = 40;

density=1./sqrt(2*pi*sigma^2*T).*exp(-(X-mu*T).^2./(2*sigma^2*T));

contours=linspace(0,1,100);
contourf(X,T,density,contours,'edgecolor','none');
caxis([0,1]);




