Message Boards Message Boards

Solving Black-Scholes PDE with NDSolve

Sandu Ursu
Posted 11 years ago

Hi, there

I've tried to replicate the procedure described in the Black-Scholes PDE numerically solved in Mathematica on Youtube, and encountered a problem.

So.. I do exactly as described in the video, but as shown below, I get a different result and I don't know why.

T = 1;
sigma = 0.3;
r = 0.05;
K = 10;

sol = NDSolve[{D[V[S, t], t] + r*S*D[V[S, t], S] + 1/2 sigma^2 S^2 D[V[S, t], {S, 2}] - r V[S, t] == 0, V[S, T] == Max[S - K, 0], V[1000, t] == 1000 - K, V[0, t] == 0},  V, {S, 0.1, 1000}, {t, 0, T}]

When I evaluate the resulting function for a price of 10 at time 0:

 V[10, 0] /. sol

I get {5.96577}, which is different from the price of {1.42313} given by the formula:

FinancialDerivative[{"European", "Call"}, {"StrikePrice" -> 10.00, "Expiration" -> 1},  {"InterestRate" -> 0.05, "Volatility" -> sigma, "CurrentPrice" -> 10}]

Is there a problem with NDSolve[], or I'm just ignoring some recent change implemented in Mathematica 10?
POSTED BY: Sandu Ursu
2 Replies
Sort By:

The difference is caused by precision: the interval for S is quite long, so step size will be large. You could decrease it to max 100, or you could set PrecisionGoal.
POSTED BY: Sjoerd de Vries

Dear Sandu,

it's not MMA 10. MMA9 gives 5.96577, too. I cannot get the 1.42269 from that Youtube video in either MMA9 nor MMA10. I also do not get the warning shown in the Youtube video. I also cannot find a mistake in the equation.

The 3D plot 

Plot3D[V[S, t] /. sol, {S, 0.1, 25}, {t, 0, 1}]

also looks quite different from what it should look. More confusingly, in order to get the 3D plot of the video, I have to solve this differential equation:

sol = NDSolve[{D[V[S, t], t] + r*S*D[V[S, t], S]^2 + 
     0.5 sigma*sigma* S ^2*D[V[S, t], {S, 2}] - r *V[S, t] == 0, 
   V[S, T] == Max[S - K, 0], V[1000, t] == 1000 - K, V[0, t] == 0}, 
  V, {S, 0.1, 1000}, {t, 0, T}]

which has the derivative of the second term squared. I must be overlooking something, but I do not know what.

If I set the exponent of the S slightly larger than 1 it seems to work as well:

sol = NDSolve[{D[V[S, t], t] + r*S^1.04953 D[V[S, t], S] + 
     0.5 sigma*sigma* S ^2*D[V[S, t], {S, 2}] - r *V[S, t] == 0, 
   V[S, T] == Max[S - K, 0], V[1000, t] == 1000 - K, V[0, t] == 0}, 
  V, {S, 0.1, 1000}, {t, 0, T}]

There is a very sharp transition somewhere between 1.04952 and 1.04953. Okay. So this looks funny. I suppose that it has something to do with numerics. Let's try:

sol = NDSolve[{D[V[S, t], t] + r*S D[V[S, t], S] + 
     0.5 sigma*sigma* S ^2*D[V[S, t], {S, 2}] - r *V[S, t] == 0, 
   V[S, T] == Max[S - K, 0], V[1000, t] == 1000 - K, V[0, t] == 0}, 
  V, {S, 0.1, 1000}, {t, 0, T}, PrecisionGoal -> 10]

Now

V[S0, 0] /. sol

gives 1.42291. Success! The graph looks better, too:

enter image description here

Cheers, Marco
POSTED BY: Marco Thiel
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract