In [14]:
import math
import time

import numpy as np

from numba import cuda


@cuda.jit(inline=True, device=True)
def cnd_cuda(d):
    A1 = 0.31938153
    A2 = -0.356563782
    A3 = 1.781477937
    A4 = -1.821255978
    A5 = 1.330274429
    RSQRT2PI = 0.39894228040143267793994605993438
    K = 1.0 / (1.0 + 0.2316419 * math.fabs(d))
    ret_val = (RSQRT2PI * math.exp(-0.5 * d * d) *
               (K * (A1 + K * (A2 + K * (A3 + K * (A4 + K * A5))))))
    if d > 0:
        ret_val = 1.0 - ret_val
    return ret_val

#https://numba.pydata.org/numba-examples/examples/finance/blackscholes/results.html

@cuda.jit
def black_scholes_cuda_kernel(callResult, X, T, R, V, S):
    i = cuda.threadIdx.x + cuda.blockIdx.x * cuda.blockDim.x
    if i >= X.shape[0]:
        return
    sqrtT = math.sqrt(T[i])
    d1 = (math.log(S / X[i]) + (R + 0.5 * V * V) * T[i]) / (V * sqrtT)
    d2 = d1 - V * sqrtT
    cndd1 = cnd_cuda(d1)
    cndd2 = cnd_cuda(d2)

    expRT = math.exp((-1. * R) * T[i])
    callResult[i] = (S * cndd1 - X[i] * expRT * cndd2)


def black_scholes_cuda(optionStrike, optionYears, Riskfree, Volatility, StockPrice):

    blockdim = 512, 1
    griddim = int(math.ceil(float(len(optionStrike))/blockdim[0])), 1

    stream = cuda.stream()

    d_callResult = cuda.device_array_like(optionStrike, stream)
    d_optionStrike = cuda.to_device(optionStrike, stream)
    d_optionYears = cuda.to_device(optionYears, stream)

    black_scholes_cuda_kernel[griddim, blockdim, stream](
            d_callResult, d_optionStrike,
            d_optionYears, Riskfree, Volatility, StockPrice)
    callResult = d_callResult.copy_to_host(stream=stream)
    stream.synchronize()

    return callResult

In [15]:
optionStrike = np.linspace(95, 100, 1000000).astype(np.float32)

In [16]:
optionYears = np.linspace(0.01, 2, 1000000).astype(np.float32)

In [17]:
StockPrice = 101.1

In [18]:
Riskfree = 0.02

In [19]:
Volatility = 0.1

In [20]:
a = black_scholes_cuda(optionStrike, optionYears, Riskfree, Volatility, StockPrice)

In [21]:
a.shape


Out[21]:
(1000000,)

In [22]:
#b.shape

In [25]:
%timeit black_scholes_cuda(optionStrike, optionYears, Riskfree, Volatility, StockPrice)


100 loops, best of 3: 9.41 ms per loop

In [26]:
1/0.00799  #million per second


Out[26]:
125.15644555694617

In [ ]: