I want to price exotic options under the exponential VG model and Merton's model to compare both models.
To price exotics under Merton's model, I have written the code below. The output is the price of a Call option, Asian, Digital and Up and in Barrier Call option. However, the use of loops leads to a very slow computation. Is there a clever way to not use loops here? In case of the VG model, I can do it but in this case, I do not see it.
function [Call,Asian,Digital,UIBP] = ExoticPricingMerton(S0,K,mu,delta,lambda,sigma,r,q,Maturity,H) ht = 1/252; %trading days P = 10^3; %Number of simulations grid = (0:ht:Maturity); N = length(grid); omega = r-q-((1/2)*sigma^2+lambda*(exp(mu+(1/2)*delta^2)-1)); S = zeros(P,N); S(:,1) = S0; for i=1:P for j=2:length(grid) N = poissrnd(lambda*ht); J = cumsum([0, normrnd(mu,delta,1,N)]); Z = normrnd(0,1); S(i,j) = S(i,j-1)*exp(omega*ht + sigma*sqrt(ht)*Z + J(end)); end end %European Call option A = max(S(:,end)-K,0); Call = exp(-r*Maturity)*(1/P)*sum(A); %Asian option A = max(mean(S,2) - K,0); Asian = exp(-r*Maturity)*(1/P)*sum(A); %Digital price A = max(S(:,end) - K, 0)./(S(:,end)-K); Digital = exp(-r*Maturity)*(1/P)*sum(A); %Up-and-in out Barrier A1 = (max(S,,2)-H)./abs(max(S,,2)-H); A2 = max(A1,0); A = (max(S(:,end)-K,0)).*A2; UIBP = exp(-r*Maturity)*(1/P)*sum(A); end