The main scope of this thesis is to implement a structured numerical analysis to check the exactness and applicability of the famous Kirk formula (1995) and the newer Bjerksund-Stensland formula (2011) widely used by energy markets practitioners while pricing and hedging (bivariate and trivariate) spread options when the strike price is different from zero. This research found that by varying volatilities, drifts, correlations, strikes, exercise times, heating rates and initial price of emission-certificates, these two analytical approximations have limitations for pricing and hedging spread options. Notably the more recent Bjerksund-Stensland formula, which is supposed to be an improvement on the Kirk formula, is not better to provide reliable result in three-dimensional trading markets. This is important, as energy markets often are three-dimentional. It will be shown mathematically with numerical experiments that both approximations provide acceptable results for pricing bivariate spread options with respect to positive strike prices. But their performances are unsatisfactory for negative strike prices. Furthermore, neither of them performed well to price trivariate spread options. And both performed poorly in hedging trivariate spread options. Although using a closed-form formula is very attractive for practitioners, this research proposes that it is safer to keep using the slower Monte Carlo numerical method, until future researches perfect existing closed-form formulas or discover a new one.