# How to simulate a CIR process using GPU and Matlab

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?

You should write some kernel functions in CUDA (Nvidia language) for your matlab code. Arrayfun is quite restrictive and not appropriate. Look at this link http://fr.mathworks.com/help/distcomp/run-cuda-or-ptx-code-on-gpu.html for more details about matlab and parallel computing.

There are some restrictions to using arrayfun. You can read the restrictions here.

Judging from the error, you cannot use indexes the way you are.

You probably have to create separate GPU arrays for $V_{t+1}$ and $V_t$.

I suggest that you find similar examples in Matlab's website and try to replicate its functionality. Here is an article with information on how to program the second order wave equation on a GPU with Matlab.