Implementing Ito Taylor Method to SDE

133 Views Asked by At

I am trying to apply the Ito Taylor method to the SDE. I have attached what I'm using for the Ito Taylor because I would like for it to be double checked please.enter image description here enter image description hereI'm supposed to get a strong speed of convergence of 1.5, but I am getting a negative number, and I don't understand why.

SDE: dX = Xdt + XdW, X(0)=1
Exact: X(t) = exp(0.5*t + W(t))

I believe the problem lies in either my X_app and/or X_exact. Here's my code:

% Ito Taylor method on linear SDE
% SDE is dX = X*dt + X*dW
% X(0) = 1
% For 0 <= t <= 1
% Time steps: 1/50, 1/100, 1/150, 1/200
% Exact Solution: X(t) = exp(0.5*t + W(t))
% K = 50, 100, 150, 200; % number of steps
% T = 1; % end of interval
% tau = T/K; % time step
% A for loop for the 4 timesteps (1/50, 1/100, 1/150, 1/200)
for f = 1:4
    K = 50*f; % So we can get 50, 100, 150, 200
    T = 1; 
    tau(f) = T/K;
    N_s = 2000; % number of trajectories
    W = zeros(N_s,K+1);  % setting up discrete Brownian path
    dW = zeros(N_s,K);  % setting up Brownian increments
    X_app = zeros(N_s,K+1); % setting up Approx solution
    X_exact = zeros(N_s,K+1); % setting up Exact solution
    Xerr = zeros(1,K+1); % setting up Global error
% Compute 2000 independent Brownian Motion trajectories
for q = 1:N_s 
    for p = 1:K 
        W(q,p+1) = W(q,p) + sqrt(tau(f))*normrnd(0,1);
    end
end
% Use these Brownian Motion trajectories to define the increments
for c = 1:N_s
    for n = 1:K
        dW(c,n) = W(c,n+1) - W(c,n);
    end
end
% Generate 2000 solution trajectories using the IT update & these increments
for j = 1:N_s
    for k = 1:K
        X_app(j,1) = 1;
        X_app(j,k+1) = X_app(j,k) + X_app(j,k)*tau(f) + X_app(j,k)*dW(j,k) + 
                       X_app(j,k)*((tau(f))^2)*0.5 + 
                       X_app(j,k)*tau(f)*dW(j,k) + X_app(j,k)*(0.5* 
                       ((dW(j,k)^2)-tau(f))) + X_app(j,k)* 
                       (((1/6)*dW(j,k)^3)-(tau(f)*dW(j,k)*0.5));
    end
 end
 % Estimate the global error of tau for the timestep using the exact 
 % solution
 for tk = 1:K+1
     for g = 1:N_s
         X_exact(g,tk) = exp(0.5.*tk + W(g,tk));
         Xerr(g,tk) = abs(X_exact(g,tk) - X_app(g,tk));
 end
 Err(tk) = (1/N_s)*sum(Xerr(:,tk));
 end
 GErr(f) = max(Err)
 % Make a log-log plot of the data through a scatter plot
 scatter(log(tau),log(GErr),'filled'), hold on
 end
 hold off
 % Finding the slop of the line is the strong order of convergence
 di = polyfit(log(tau),log(GErr),1);
 slope = di(1)