We develop and apply a numerical scheme for pricing options for the stochastic volatility model proposed by Barndorff-Nielsen and Shephard. This non-Gaussian Ornstein-Uhlenbeck type of volatility model gives rise to an incomplete market, and we consider the option prices under the minimal entropy martingale measure. To price numerically options with respect to this risk neutral measure, one needs to consider a Black & Scholes type of partial differential equation, with an integro-term arising from the volatility process. We suggest finite difference schemes to solve this parabolic integro-partial differential equation, and derive appropriate boundary conditions for the finite difference method. As an application of our algorithm, we consider price deviations from the Black & Scholes formula for call options, and the implications of the stochastic volatility on the shape of the volatility smile.