## Multithreading a Monte Carlo simulation

One of the more nebulous areas of computer science is concurrency, with multi-threaded programming being particularly difficult to understand if you haven’t had the benefit of lots of experience. Of course, the precise degree to which your multi-threaded code is “difficult” depends hugely on the problem that the code is intended to solve (multi-threaded graphical user interfaces, for instance, are notoriously difficult to get working). There are, however, a number of problem domains that are eminently suited to multi-threading.

One of these is the problem of Monte Carlo simulation. I had a discussion about this yesterday with someone who is in the process of learning mathematical finance. He loves the subject but is frustrated by the fact that he’s essentially using only one of the four cores on his machine whenever he writes code to, say, price an option using a Monte Carlo algorithm. Surely, he said, there must be an easy way to get a multi-threaded Monte Carlo routing going?

Given that he’s been reading Mark Joshi’s excellent books on mathematical finance and computational finance, I thought that I’d run through the simplest possible example of how multi-threading his routine would result in faster run times. Here’s the code I knocked together in a few minutes based on a trivial Monte Carlo pricer for a call option. (Note: this was for demonstration purposes only, intended to illustrate the ideas behind threading and the Boost Thread library in particular; writing code such as this in a bank will get you fired at your first code review!)

#include <iostream> using namespace std; #include <boost/thread.hpp> #include <boost/timer.hpp> #if defined(_USE_BOOST_MATH) #include <boost/random/mersenne_twister.hpp> #include <boost/random/normal_distribution.hpp> #include <boost/random/variate_generator.hpp> typedef boost::mt19937 T_base_prng; typedef boost::normal_distribution<> T_norm_dist; typedef boost::variate_generator<T_base_prng&, T_norm_dist> T_norm_varg; #endif #if !defined(_MSC_VER) #define EXIT_SUCCESS 0 #define EXIT_FAILURE 1 #endif double getOneGaussianByBoxMuller(); double SimpleMC(double expiry, double strike, double spot, double rate, double vol, unsigned long numPaths, unsigned long base_seed = 42); void SimpleMC_threaded(double expiry, double strike, double spot, double rate, double vol, unsigned long numPaths, double* target, unsigned long base_seed = 42); int main(int argc, char* argv[]) { try { double expiry = 1; double strike = 100; double spot = 100; double rate = 0.05; double vol = 0.20; unsigned long numPaths = 1e8; boost::timer timer; double result_unthreaded = SimpleMC(expiry, strike, spot, rate, vol, numPaths); double time_unthreaded = timer.elapsed(); cout << "Single-threaded run\n" << "-------------------\n" << "Call: " << result_unthreaded << endl << "Time: " << time_unthreaded << "\n\n"; double result1; double result2; double result3; double result4; boost::timer timer_threaded; boost::thread worker1(SimpleMC_threaded, expiry, strike, spot, rate, vol, numPaths / 4, &result1, 41); boost::thread worker2(SimpleMC_threaded, expiry, strike, spot, rate, vol, numPaths / 4, &result2, 42); boost::thread worker3(SimpleMC_threaded, expiry, strike, spot, rate, vol, numPaths / 4, &result3, 43); boost::thread worker4(SimpleMC_threaded, expiry, strike, spot, rate, vol, numPaths / 4, &result4, 44); worker1.join(); worker2.join(); worker3.join(); worker4.join(); double result_threaded = (result1 + result2 + result3 + result4) / 4; double time_threaded = timer_threaded.elapsed(); cout << "Multi-threaded run\n" << "------------------\n" << "Call: " << result_threaded << endl << "Time: " << time_threaded << "\n\n"; #if !defined(_MSC_VER) cout << "Press Enter to continue."; cin.get(); #else system("PAUSE"); #endif return (EXIT_SUCCESS); } catch (exception& e) { cerr << e.what() << endl; return (EXIT_FAILURE); } catch (...) { cerr << "Non-standard exception caught... " << endl; return (EXIT_FAILURE); } } double SimpleMC(double expiry, double strike, double spot, double rate, double vol, unsigned long numPaths, unsigned long base_seed) { double variance = vol * vol * expiry; double rootVariance = sqrt(variance); double itoCorrection = -0.5 * variance; double movedSpot = spot * exp(rate * expiry + itoCorrection); double runningSum = 0; double thisSpot; #if !defined(_USE_BOOST_MATH) for (unsigned long i = 0; i < numPaths; ++i) { double thisGaussian = getOneGaussianByBoxMuller(); thisSpot = movedSpot * exp(rootVariance * thisGaussian); double thisPayoff = thisSpot - strike; thisPayoff = thisPayoff > 0 ? thisPayoff : 0; runningSum += thisPayoff; } #else T_base_prng base_prng(base_seed); T_norm_dist norm_dist(0, 1); T_norm_varg norm_varg(base_prng, norm_dist); for (unsigned long i = 0; i < numPaths; ++i) { thisSpot = movedSpot * exp(rootVariance * norm_varg()); double thisPayoff = thisSpot - strike; thisPayoff = thisPayoff > 0 ? thisPayoff : 0; runningSum += thisPayoff; } #endif double mean = runningSum / numPaths; mean *= exp(-rate * expiry); return mean; } void SimpleMC_threaded(double expiry, double strike, double spot, double rate, double vol, unsigned long numPaths, double* target, unsigned long base_seed) { double variance = vol * vol * expiry; double rootVariance = sqrt(variance); double itoCorrection = -0.5 * variance; double movedSpot = spot * exp(rate * expiry + itoCorrection); double runningSum = 0; double thisSpot; #if !defined(_USE_BOOST_MATH) for (unsigned long i = 0; i < numPaths; ++i) { double thisGaussian = getOneGaussianByBoxMuller(); thisSpot = movedSpot * exp(rootVariance * thisGaussian); double thisPayoff = thisSpot - strike; thisPayoff = thisPayoff > 0 ? thisPayoff : 0; runningSum += thisPayoff; } #else T_base_prng base_prng(base_seed); T_norm_dist norm_dist(0, 1); T_norm_varg norm_varg(base_prng, norm_dist); for (unsigned long i = 0; i < numPaths; ++i) { thisSpot = movedSpot * exp(rootVariance * norm_varg()); double thisPayoff = thisSpot - strike; thisPayoff = thisPayoff > 0 ? thisPayoff : 0; runningSum += thisPayoff; } #endif double mean = runningSum / numPaths; mean *= exp(-rate * expiry); *target = mean; } double getOneGaussianByBoxMuller() { double result; double x; double y; double sizeSquared; do { x = 2.0 * rand() / static_cast<double>(RAND_MAX) - 1; y = 2.0 * rand() / static_cast<double>(RAND_MAX) - 1; sizeSquared = x * x + y * y; } while (sizeSquared >= 1.0); result = x * sqrt(-2 * log(sizeSquared) / sizeSquared); return result; }

## Short straws

Phew! Just my luck: I finally decide to commit to starting a blog when all hell breaks loose with work and I have to pull eighteen-hour days just to keep on top of things. What’s even more annoying is that I actually have a ton of posts in preparation but I’ve discovered that – just like when I’m trying to write a paper – I’m never happy with the results and keep trying to hone and tweak them to an inevitably unattainable level of perfection.

The moral of all this is of course that chasing perfection is a great idea when you’re trying to get a paper published, but it’s an absolutely terrible idea when you’re trying to get a blog off the ground. I shall have to pay attention to this in future.

Anyway, while we’re waiting for the good stuff to come along, I thought it might be worthwhile mentioning a conversation I had at lunch yesterday. I got into a discussion with a colleague about finance, a discussion which took rather a surprising turn. He’s someone whose knowledge on all things physics related genuinely impresses me but he hit me with somewhat of a bombshell half way through our conversation. We were talking about how large – sometimes stupidly large – amounts of money are made each day on stock exchanges by trading in shares. He understood how someone could make money by buying a stock and holding it in the expectation that its price could rise in value. But – and this is a big but – he had no idea how short selling works.

I’ve talked to some other people since then and almost all of them seem to have the same difficulty. I admit, short-selling stocks is a weird idea when you first encounter it, but it really couldn’t be simpler. Here’s how I explained it.

Suppose your next-door neighbour loves video games. In particular, he’s very fond of his Playstation 3 and tries to keep up with all of the latest game releases. Several months ago, he went to Amazon and pre-ordered a copy of an eagerly awaited game so that he’d be sure of having it delivered on its day of release. He waits and waits until, eventually, the release date arrives, his credit card is charged £40, and the postman pops a copy of it through his letterbox. However, being a busy man, he’s due to go on a business trip later that day and he won’t be able to play the game until he gets back in, say, a month’s time.

You’ve been watching his growing excitement about the release of the game with interest. You think to yourself, “Hey, if he’s this excited about the game surely others are too. And if others are excited about it, perhaps there’s a chance to make a few bob here.” (You say this to yourself in your best Del Boy voice.) Knowing that he’s going away for a while – and knowing that the game is in such high demand that there are actually shortages in stores – you ask him if he would mind you borrowing it while he’s gone, promising to give it back to him the moment he arrives back from his trip. In return for allowing you to borrow the game, you might decide to thank him by promising to keep an eye on his house while he’s gone.

Notice that you promised only to give it back to him when he comes home again; in particular, you didn’t actually tell him what you were going to do with the game! Instead of sitting down to play it, you decide to put it to work for you. You go to ebay, say, and place the copy of the game for sale there. You find to your delight that such is the enormous demand for the game, people are willing to pay a huge premium in order to get their hands on it. The ebay auction actually ends with one poor sucker bidder buying the game from you for £80 – *double what your friend paid!* You smile as your PayPal account grows by £80 and pop the game in the post to the buyer.

A month passes and your neighbour is due back in a few hours. Remembering that you promised to give him back a copy of the game, you stroll down to your local video games store and find – surprise! – that plenty of copies are now available. Even better, they’re available for only £30 now that the initial excitement has worn off. You pick up a copy and on the way home you buy a couple of six packs of beer, one of which is to thank your friend for letting you have the game.

That evening, you smile to yourself as you enjoy your beer.

Contrived though this example is, it captures essentially everything you need to know about short selling stocks. Broadly speaking, a trader will short sell a stock if he or she believes that the stock is over valued or that the stock price will drop over a certain time frame. Although different regulatory requirements mean that the details can sometimes change, the basic procedure is always the same: the trader seeks out someone who already has the stock and asks to borrow it from them for a fixed period of time in return for a small fee. He then sells the stock on the open market and puts the cash in his account. When the time comes for him to return the stock to the person who lent it to him, the trader buys the stock from the market, hoping that its price has gone down in the meantime. If it has, he makes a profit; if it hasn’t, he has to bear the loss involved with buying the stock and giving it back to its original owner.

There are, of course, plenty of reasons why shorting stocks in practice isn’t as simple as this (regulatory requirements, temporary freezes on shorts, managing the dividend payments on stocks you’ve shorted, stock splits that happen during a short, and so forth). Shorting also has rather a bad image among the public. But it’s really quite a straightforward – and ultimately beneficial – practice.

And it’s all done so that you can enjoy that beer when your short comes home!

## Physics, but from a hedgerow

Welcome one and all to physics from a hedgerow. Sure, the name is unwieldy but it has a significance I might explain one day.

For the moment I’m still trying to get used to wordpress so don’t expect a great deal of activity over the coming week. Beyond that, however, we’ll be taking some — hopefully interesting — trips into physics, artificial intelligence, and money. Mostly physics of course but rather a lot about the others too.