analytic solution of Black-Scholes PDE

Here we presentan analyticalsolution for the Black-Scholes partial differential equation,

$rf={\displaystyle \frac{\partial f}{\partial t}}+rx{\displaystyle \frac{\partial f}{\partial x}}+{\displaystyle \frac{1}{2}}{\sigma }^2{x}^2{\displaystyle \frac{{\partial }^{2}f}{\partial x^2}},f=f(t,x),$(1)

over the domain ,with terminal condition $f(T,x)=\psi (x)$,by reducing this parabolic PDE tothe heat equation of physics.

We begin by making the substitution:

$$u=e^{-rt}f,$$

which is motivated by the fact that it is the portfolio valuediscounted by the interest rate $r$ (see the derivation of theBlack-Scholes formula)that is a martingale.Using the product rule on $f=e^{rt}u$,we derive the PDE that the function $u$ must satisfy:

$$rf=r{e}^{rt}u=r{e}^{rt}u+e^{rt}\frac{\partial u}{\partial t}+rx{e}^{rt}\frac{\partial u}{\partial x}+\frac{1}{2}{\sigma }^2{x}^2{e}^{rt}\frac{{\partial }^{2}u}{\partial x^2};$$

or simply,

$0={\displaystyle \frac{\partial u}{\partial t}}+rx{\displaystyle \frac{\partial u}{\partial x}}+{\displaystyle \frac{1}{2}}{\sigma }^2{x}^2{\displaystyle \frac{{\partial }^{2}u}{\partial x^2}}.$(2)

Next, we make the substitutions:

$$y=\log x,s=T-t.$$

These changes of variables can be motivatedby observing that:

• The underlying process described by the variable $x$is a geometric Brownian motion(as explained in the derivation of the Black-Scholes formula itself),so that $\log x$ describes a Brownian motion, possiblywith a drift.Then $\log x$ should satisfy some sort of diffusion equation(well-known in physics).

• The evolution of the systemis backwards from the terminal state of the system. Indeed,the boundary condition is given as a terminal state,and the coefficient of $\partial u/\partial t$ is positive in equation (2).(Compare with the standard heat equation,$0=-\partial u/\partial t+\partial u/\partial x$, which describes atemperature evolving forwards in time.)So to get to the heat equation,we have to use a substitution to reverse time.

Since

$$\frac{\partial u}{\partial s}=-\frac{\partial u}{\partial t},\frac{\partial u}{\partial x}=\frac{\partial u}{\partial y}\frac{dy}{dx}=\frac{1}{x}\frac{\partial u}{\partial y},$$

and

$$\frac{{\partial }^{2}u}{\partial x^2}=\frac{\partial }{\partial x}\left(\frac{1}{x}\frac{\partial u}{\partial y}\right)=-\frac{1}{x^2}\frac{\partial u}{\partial y}+\frac{1}{x^2}\frac{{\partial }^{2}u}{\partial y^2},$$

substituting in equation (2),we find:

$0=-{\displaystyle \frac{\partial u}{\partial s}}+(r-\frac{1}{2}{\sigma }^2){\displaystyle \frac{\partial u}{\partial y}}+{\displaystyle \frac{1}{2}}{\sigma }^2{\displaystyle \frac{{\partial }^{2}u}{\partial y^2}}.$(3)

The first partial derivative with respect to $y$does not cancel (unless $r=\frac{1}{2}{\sigma }^2$)because we have not take into accountthe drift of the Brownian motion.To cancel the drift (which is linear in time),we make the substitutions:

$$z=y+(r-\frac{1}{2}{\sigma }^2)\tau ,\tau =s.$$

Under the new coordinate system $(z,\tau )$, we have the relationsamongst vector fields:

$$\frac{\partial }{\partial z}=\frac{\partial }{\partial y},\frac{\partial }{\partial \tau }=-(r-\frac{1}{2}{\sigma }^2)\frac{\partial }{\partial y}+\frac{\partial }{\partial s},$$

leading to the following of equation (3):

$0=-{\displaystyle \frac{\partial u}{\partial \tau }}-(r-\frac{1}{2}{\sigma }^2){\displaystyle \frac{\partial u}{\partial z}}+(r-\frac{1}{2}{\sigma }^2){\displaystyle \frac{\partial u}{\partial z}}+{\displaystyle \frac{1}{2}}{\sigma }^2{\displaystyle \frac{{\partial }^{2}u}{\partial z^2}};$

or:

${\displaystyle \frac{\partial u}{\partial \tau }}={\displaystyle \frac{1}{2}}{\sigma }^2{\displaystyle \frac{{\partial }^{2}u}{\partial z^2}},u=u(\tau ,z),$(4)

which is one form of the diffusion equation.The domain is onand $0\le \tau \le T$;the initial condition is to be:

$u(0,z)=e^{-rT}\psi (e^z):=u_0(z).$

The original function $f$ can be recoveredby

$$f(t,x)=e^{rt}u(T-t,\log x+(r-\frac{1}{2}{\sigma }^2)\tau ).$$

The fundamental solution of the PDE (4)is known to be:

$$G_{\tau }(z)=\frac{1}{\sqrt{2\pi {\sigma }^2\tau }}\exp \left(-\frac{z}{2{\sigma }^2\tau }\right)$$

(derived using the Fourier transform);and the solution $u$ with initial condition $u_0$is given by the convolution:

$$u(\tau ,z)=u_0*G_{\tau }(z)=\frac{e^{-rT}}{\sqrt{2\pi {\sigma }^2\tau }}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\psi (e^{\zeta })\exp \left(-\frac{{(z-\zeta )}^2}{2{\sigma }^2\tau }\right)𝑑\zeta .$$

In terms of the original function $f$:

$$f(t,x)=\frac{e^{-r\tau }}{\sqrt{2\pi {\sigma }^2\tau }}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\psi (e^{\zeta })\exp \left(-\frac{{\left(\log x+(r-\frac{1}{2}{\sigma }^2)\tau -\zeta \right)}^2}{2{\sigma }^2\tau }\right)𝑑\zeta ,$$

($\tau =T-t$) which agrees with the result derived using probabilistic methods (http://planetmath.org/BlackScholesFormula).