I am trying to simulate a CIR process using Matlab and my GPU for effeciency. At the moment i run into some implementation problems due to the recursive nature of the discretization.
The sheme I currently use is the simple Euler
$$V_t \approx \max \{V_{t-1} + \kappa_V(\bar{V}-V_{t-1})\Delta h + \sigma_V\sqrt{V_{t-1}}\sqrt{\Delta h}\mathcal{N}(0,1),0\}. $$
To proper simulate on the GPU I need some help function, like discribed in this Matlab guide
function [ V ] = simulation_fun( V, W, N, h, kappa_V, V_bar, sigma_V)
t=1;
while t < N % a path has N steps
V(:,t+1)=max(V(:,t)+kappa_V*(V_bar-V(:,t))*h + sigma_V*sqrt(V(:,t))*sqrt(h)*W(:,t),0);
t=t+1;
end
end
Than I define the model parameters
N=1000;
M=1000; % Number of MC simulations
h=0.01;
kappa_V=1;
V_bar=1;
sigma_V=0.2;
V=gpuArray(ones(M,N));
W=gpuArray.randn(M,N);
After that I call the simulation via
V = arrayfun(@simulation_fun,V, W, N, h, kappa_V, V_bar, sigma_V);
and get the follwing error
Error using gpuArray/arrayfun
Indexing is not supported. error at line: 3 .
However, if I use
V = simulation_fun(V, W, N, h, kappa_V, V_bar, sigma_V);
it works, but the simulation is roughly 2.5 times slower than only using the CPU.
Does anyone know how to correctly implement the simulation?