This page shows recipe 15.10, page 583, which I contributed for Sal Mangano's "Mathematica Cookbook", published by O'Reilly.
The method is very fast as it is 100% explicit. It solves the linear complementarity problem of the open boundary value problem with a PSOR (projective successive overrelaxation) method using three nodes for every time step, being fully Compile[]d. The basic steps are transformation, iteration, and backtransformation.
What makes this rectangular grid method so powerful is the fact that although it is faster than most tree-based implementations, it computes the option prices for the whole interval, not just for ONE particular price of the underlying, which is a limitation all tree-based methods possess.

The following example uses time=1 year, vol=0.4, rate=0.05, and strikes from 50 through 100 in steps of 5. a is 5, nn is 100, and mm is (initially) 20. This corresponds to an alpha of 0.4. For mm = 16 alpha is exactly 0.5, that is the stability limit. For lower values of mm alpha is larger then 0.5, resulting in complete blow-up (instability).

Discussion:
The essential items are already discussed in the book, page 584. For a 100% explicit method, it is crucial to keep s/h^2 ≤ 1/2. This is the stability limit. That means, if the spatial step size is reduced by a factor of f, the temporal step size has to be reduced by a factor of f^2. Thus, if mm, the number of temporal steps, is lowered below 16 a value of alpha is used that is larger than 1/2, violating the stability limit, and the instability becomes quite apparent. On the other hand, increasing the value of mm results in more time steps to be used (smaller step sizes), causing the method to run longer, as can be seen by dragging the slider and observing the number of milliseconds for the running time in the plot label. Therefore, the proper choice of mm is best taken as the lowest number that gives accurate results. Of course, nn and a also factor into this (nn can be increased from the default value of 100, but that also increases the number of spatial steps used, slowing down the method). As the purpose of this webMathematica example is to illustrate the behavior, nn is fixed at 100, a is fixed at 5, and only mm can be changed.
The precision and speed results are not only remarkable because of the high precision results being computed in ~2 ms per strike (~6 ms for mm=50). This in itself is much faster than several binomial tree implementations in C++ I have seen over the years computing just ONE option and delivering inferior precision! They are particularly remarkable because the method returns the option prices for ALL 201 (2 nn +1) stock prices used as nodes in the spatial interpolation, unlike traditional tree methods, which compute option prices only "strike by strike". The total running time shown above is for 11 x 201 = 2211 spot values and 11 x mm time values. It is the fastest 100% explicit trinomial method I am aware of (I welcome competing algorithms that deliver the same precision quality with faster speed).

Literature
This is a direct implementation of a derivation presented in:
Prof. Ansgar Jüngel: "Modellierung und Numerik von Finanzderivaten", Vorlesungsmanuskript, Johannes-Gutenberg-Universität Mainz, 2002.

This website is powered by webMathematica.
Please visit www.webMathematica.com to find out how you can empower your website with Mathematica.

Copyright (c) 2010 Andreas Lauschke Consulting. All Rights Reserved.