Optimal Oil Production and Taxation in Presence of Global Disruptions

TL;DR

This paper develops a mathematical model for optimal oil extraction and taxation considering global market fluctuations, deriving equilibrium policies for producers and governments using stochastic processes and differential game theory.

Contribution

It introduces a regime-switching jump diffusion model for oil prices and characterizes the Nash equilibrium for extraction and taxation policies.

Findings

Existence of a Nash equilibrium in the model

Optimal extraction and tax policies derived

Numerical example illustrating the theoretical results

Abstract

This paper studies the optimal extraction policy of an oil field as well as the efficient taxation of the revenues generated. Taking into account the fact that the oil price in worldwide commodity markets fluctuates randomly following global and seasonal macroeconomic parameters, we model the evolution of the oil price as a mean reverting regime-switching jump diffusion process. Given that oil producing countries rely on oil sale revenues as well as taxes levied on oil companies for a good portion of the revenue side of their budgets, we formulate this problem as a differential game where the two players are the mining company whose aim is to maximize the revenues generated from its extracting activities and the government agency in charge of regulating and taxing natural resources. We prove the existence of a Nash equilibrium and the convergence of an approximating scheme for the value functions. Furthermore, optimal extraction and fiscal policies that should be applied when the equilibrium is reached are derived.A numerical example is presented to illustrate these results.

Equations

$$\displaystyle\left\{\begin{array}[]{ll}\mathrm{d}X(t)&\displaystyle=\kappa\big{(}\mu(\alpha(t))-X(t)\big{)}\mathrm{d}t\\ &+\sigma(\alpha(t))\mathrm{d}W(t)\\ &+\int_{{\mathbb{R}}}\gamma(\alpha(t))z\bar{N}(dt,dz),\\ dY(t)&=-u_{1}(t)dt,\\ X(s)&=x,\,\,\,\,Y(s)=y\geq 0,\,\,0\leq s\leq t\leq T,\end{array}\right.$$

$$\begin{array}[]{rl}P(t,X(t),Y(t),u_{1}(t))=e^{X(t)}u_{1}(t)-\\ C(t,Y(t),u_{1}(t)).\end{array}$$

$$\begin{array}[]{rl}&\!\!\!\displaystyle L_{1}(t,X(t),Y(t),u_{1}(t),u_{2}(t))\\ &\!\!\!\displaystyle=\theta P(t,X(t),Y(t),u_{1}(t))(1-u_{2}(t)),\end{array}$$

$$\begin{array}[]{rl}&\!\!\!\displaystyle L_{2}(t,X(t),Y(t),u_{1}(t),u_{2}(t))\\ &\!\!\!\displaystyle=(1-\theta)P(t,X(t),Y(t),u_{1}(t))+\\ &\!\!\!\displaystyle u_{2}(t)\theta P(t,X(t),Y(t),u_{1}(t)).\end{array}$$

$$\begin{array}[]{rl}&\!\!\!\displaystyle J_{i}(s,x,y,\iota;u_{1},u_{2})\\ &\!\!\!\displaystyle=E\Bigg{[}\int_{s}^{T}e^{-r(t-s)}L_{i}(t,X(t),Y(t),u_{1}(t),u_{2}(t))dt\\ &\!\!\!\displaystyle+e^{-r(T-s)}\Phi_{i}(X(T),Y(T))\bigg{|}X(s)=x,Y(s)=y,\\ &\!\!\!\displaystyle\alpha(s)=\iota\Bigg{]}.\end{array}$$

$$\begin{array}[]{rl}&\!\!\!\displaystyle J_{1}(s,x,y,\iota;u_{1}^{*},u_{2}^{*})\geq J_{1}(s,x,y,\iota;u_{1},u_{2}^{*}),\\ &\!\!\!\displaystyle\hbox{ for all }u_{1}(\cdot)\in{\cal U}_{1}(s,x,y,\iota),\\ &\!\!\!\displaystyle J_{2}(s,x,y,\iota;u_{1}^{*},u_{2}^{*})\geq J_{2}(s,x,y,\iota;u_{1}^{*},u_{2}),\\ &\!\!\!\displaystyle\hbox{ for all }u_{2}(\cdot)\in{\cal U}_{2}(s,x,y,\iota).\end{array}$$

$$V_1(s,x,y,\iota )=\underset{u_1\in {\mathcal{U}}_1}{\sup }{J}_1(s,x,y,\iota ;u_1,u_2^{\ast })$$

$$V_2(s,x,y,\iota )=\underset{u_2\in {\mathcal{U}}_2}{\sup }{J}_2(s,x,y,\iota ;u_1^{\ast },u_2)$$

$$\displaystyle H_{1}\bigg{(}s,x,y,\iota,V,\frac{\partial V}{\partial x},\frac{\partial V}{\partial y},\frac{\partial^{2}V}{\partial x^{2}}\bigg{)}$$

$$\displaystyle=rV-\sup_{u_{1}\in U_{1}}\Bigg{(}\frac{1}{2}\sigma^{2}(\iota)\frac{\partial^{2}V}{\partial x^{2}}+\kappa\big{(}\mu(\iota)-x\big{)}\frac{\partial V}{\partial x}$$

$$\displaystyle-u_{1}\frac{\partial V}{\partial y}+\int_{{\mathbb{R}}}\bigg{(}V(s,x+\gamma(\iota)z,y,\iota)-V(s,x,y,\iota)$$

$$\displaystyle-{\bf 1}_{\{|z|<1\}}(z)\frac{\partial V}{\partial x}\gamma(\iota)z\bigg{)}\nu(dz)$$

$$\displaystyle+L_{1}(s,x,y,u_{1},u_{2}^{*})+QV(s,x,y,\cdot)(\iota)\Bigg{)},$$

$$\displaystyle H_{2}\bigg{(}s,x,y,\iota,V,\frac{\partial V}{\partial x},\frac{\partial V}{\partial y},\frac{\partial^{2}V}{\partial x^{2}}\bigg{)}$$

$$\displaystyle=rV-\sup_{u_{2}\in U_{2}}\Bigg{(}\frac{1}{2}\sigma^{2}(\iota)\frac{\partial^{2}V}{\partial x^{2}}+\kappa\big{(}\mu(\iota)-x\big{)}\frac{\partial V}{\partial x}$$

$$\displaystyle-u_{1}^{*}\frac{\partial V}{\partial y}+\int_{{\mathbb{R}}}\bigg{(}V(s,x+\gamma(\iota)z,y,\iota)-V(s,x,y,\iota)$$

$$\displaystyle-{\bf 1}_{\{|z|<1\}}(z)\frac{\partial V}{\partial x}\cdot\gamma(\iota)z\bigg{)}\nu(dz)$$

$$\displaystyle+L_{2}(s,x,y,u_{1}^{*},u_{2})+QV(s,x,y,\cdot)(\iota)\Bigg{)},$$

$$\displaystyle\left\{\begin{array}[]{ll}\frac{\partial V_{1}}{\partial s}=H_{1}\bigg{(}s,x,y,\iota,V_{1},\frac{\partial V_{1}}{\partial x},\frac{\partial V_{1}}{\partial y},\frac{\partial^{2}V_{1}}{\partial x^{2}}\bigg{)},\\ \frac{\partial V_{2}}{\partial s}=H_{2}\bigg{(}s,x,y,\iota,V_{2},\frac{\partial V_{2}}{\partial x},\frac{\partial V_{2}}{\partial y},\frac{\partial^{2}V_{2}}{\partial x^{2}}\bigg{)},\\ V_{1}(T,x,y,\iota)=\Phi_{1}(x,y)\\ V_{2}(T,x,y,\iota)=\Phi_{2}(x,y).\end{array}\right.$$

$$\displaystyle\displaystyle u_{1}^{*}(s)=\arg\max\bigg{(}-u_{1}\frac{\partial V_{1}(s,x,y,\iota)}{\partial y}$$

$$\displaystyle+L_{1}(s,x,y,u_{1},u^{*}_{2})\bigg{)},\quad s\in[0,T]$$

$$\displaystyle\displaystyle u_{2}^{*}(s)=\arg\max\bigg{(}\displaystyle L_{2}(s,x,y,u_{1}^{*},u_{2})\bigg{)},\,s\in[0,T].$$

$$\displaystyle\left\{\begin{array}[]{ll}\frac{\partial W_{1}}{\partial s}=H_{1}\bigg{(}s,x,y,\iota,W_{1},\frac{\partial W_{1}}{\partial x},\frac{\partial W_{1}}{\partial y},\frac{\partial^{2}W_{1}}{\partial x^{2}}\bigg{)}&\\ W_{1}(T,x,y,\alpha(T))=\Phi_{1}(x,y).\end{array}\right.$$

$$\begin{array}[]{rl}V(s,x,y,i)=\left(\begin{array}[]{ll}V_{1}(s,x,y,i)\\ V_{2}(s,x,y,i)\end{array}\right),\end{array}$$

$$\begin{array}[]{rl}&\!\!\!\displaystyle\displaystyle H(s,x,y,i,V,\frac{\partial V}{\partial x},\frac{\partial V}{\partial y},\frac{\partial^{2}V}{\partial x^{2}})\\ &\!\!\!\displaystyle=\left(\begin{array}[]{ll}\displaystyle H_{1}(s,x,y,i,V_{1},\frac{\partial V_{1}}{\partial x},\frac{\partial V_{1}}{\partial y},\frac{\partial^{2}V_{1}}{\partial x^{2}})\\ \displaystyle H_{2}(s,x,y,i,V_{2},\frac{\partial V_{2}}{\partial x},\frac{\partial V_{2}}{\partial y},\frac{\partial^{2}V_{2}}{\partial x^{2}})\end{array}\right),\end{array}$$

$$\begin{array}[]{rl}\Phi(x,y)=\left(\begin{array}[]{ll}\Phi_{1}(x,y)\\ \Phi_{2}(x,y)\end{array}\right).\end{array}$$

$$\displaystyle\left\{\begin{array}[]{ll}\displaystyle\frac{\partial V}{\partial s}=H\bigg{(}s,x,y,i,V,\frac{\partial V}{\partial x},\frac{\partial V}{\partial y},\frac{\partial^{2}V}{\partial x^{2}}\bigg{)},\\ V(T,x,y,i)=\Phi(x,y).\end{array}\right.$$

$$\begin{array}[]{rl}&\!\!\!\displaystyle{\rm I}f(s,x,y,i)=\int_{{\mathbb{R}}}f(s,x+\gamma(i)z,y,i)\nu(dz)\\ &\!\!\!\displaystyle-\frac{\partial f(s,x,y,i)}{\partial x}\int_{-1}^{1}\gamma(i)z\nu(dz)-f(s,x,y,i)\Gamma.\end{array}$$

$$\begin{array}[]{rl}&\!\!\!\displaystyle{\rm I}_{\xi}f(s,x,y,i)=\sum_{j=0}^{N_{\xi}}c_{j}f(s,x+\gamma(i)z_{j},y,i)-\\ &\!\!\!\displaystyle\frac{\partial f(s,x,y,i)}{\partial x}\sum_{j=0}^{M_{\xi}}d_{j}\gamma(i)z_{j}-f(s,x,y,i)\Gamma,\end{array}$$

$$H_{1,u_2^{\ast }}^{h,k,\xi }{V}_1(s,x,y,i)$$

$$H_{2,u_1^{\ast }}^{h,k,\xi }{V}_2(s,x,y,i)$$

$$\displaystyle\left\{\begin{array}[]{ll}\displaystyle 0=\Delta_{s}V-H_{u_{1}^{*},u_{2}^{*}}^{h,k,\xi}V(s,x,y,i),&\\ V(T,x,y,i)=\Phi(x,y),\end{array}\right.$$

$$\begin{array}[]{rl}H_{u_{1}^{*},u_{2}^{*}}^{j,k,\xi}V(s,x,y,i)=\left(\begin{array}[]{ll}H_{1,u_{2}^{*}}^{h,k,\xi}V_{1}(s,x,y,i)\\ H_{2,u_{1}^{*}}^{h,k,\xi}V_{2}(s,x,y,i)\end{array}\right).\end{array}$$

$${\mathcal{F}}_{\xi }(V)(s,x,y,i;h,k)$$

$$\begin{array}[]{rl}e(i;v)&\!\!\!\displaystyle=1-\epsilon\bigg{[}r+\frac{1}{k}+\frac{\sigma^{2}(i)}{h^{2}}+\frac{1}{h}\bigg{(}\kappa(\mu(i)-x)\\ &\!\!\!\displaystyle-\sum_{n=0}^{M_{\xi}}d_{n}\gamma(i)z_{n}-v\bigg{)}+\Gamma+\sum_{n\neq i}q_{in}\bigg{]},\end{array}$$

Full text

1 Introduction

Oil and natural gas have always been the main sources of revenues for a large number of developing countries as well as some industrialized countries around the world. Oil extraction policies vary from a country to another. In some countries, the extraction is done by a state-owned corporation in others it is done by foreign multinationals. The earliest work on the extraction of natural resources was done by Hotelling (1931) who derived an optimal extraction policy under the assumption that the commodity price is constant. Many economists have proposed various extensions of the Hotelling model by taking into account the uncertainty and randomness of commodity prices.

The majority of oil extraction contracts signed between multinational oil companies and governments of oil-rich nations are in the form of profit sharing agreements where each party will take a fraction of the profits. In addition, the host nation is also entitled to collect tax from all companies. This creates a very interesting dynamic for these two parties with converging as well as conflicting interests during the lifetime of the mining contract. We formulate this problem as a differential game where the two players are the multinational oil company and the government. To the best of our knowledge, this is the first time this approach is used to characterize the interplay between extraction and taxation of oil or natural gas. It is also self-evident that the price of oil in commodity exchange markets fluctuates following divers macroeconomic and global geopolitical forces. In this paper, we use the mean reverting regime switching Lévy processes to model the oil price. Oil prices also display a great deal of seasonality, jumps, and spikes due to various supply disruptions and political turmoils in oil-rich countries. We use regime-switching jump diffusions to capture all those effects. Thus our pricing model closely captures the instability of oil markets. The paper is organized as follows. In the next section, we formulate the problem under consideration. In Section 3, we prove the existence of a Nash equilibrium. And in section 4, we construct a finite difference approximation scheme and prove its convergence to the value functions. Finally, in section 5, we give a numerical example.

2 Problem formulation

Consider a multinational oil company who enters into a Production Sharing Agreement with the government of an oil-rich country with expiration $0<T<\infty$. Both parties will share the profits from the sales of the extracted oil on world markets following a simple rule where the company takes $100\theta$ percent and the government takes $100(1-\theta)$ percent of the profits, for some $\theta\in(0,1)$ . We assume that the market value of a barrel of oil at time $t$ is $S_{t}=e^{X_{t}}$. In fact, we assume that the oil price follows an exponential Lévy model, these models are natural extensions of the celebrated Black-Scholes model. Given that oil prices are very sensitive global macroeconomic and geopolitical shocks, we model $X_{t}$ as a mean reverting regime switching Lévy process with two states. Let $\alpha(t)\in{\cal M}=\{1,2\}$ be a finite state Markov chain that captures the state of the oil market: $\alpha(t)=1$ indicates the bull market at time $t$ and $\alpha(t)=2$ represents a bear market at time $t$. Let the matrix $Q$ be the generator of $\alpha(t)$. Let $(\eta_{t})_{t}$ be a Lévy process and let $N$ be the Poisson random measure of $(\eta_{t})_{t}$. Let $\nu$ be the Lévy measure of $(\eta_{t})_{t}$. The differential form of $N$ is denoted by $N(dt,dz)$, we define the differential $\bar{N}(dt,dz)$ as follows: $\bar{N}(dt,dz)=N(dt,dz)-\nu(dz)dt$ if $|z|<1$ and $\bar{N}(dt,dz)=N(dt,dz)$ if $|z|\geq 1.$ We assume that the Lévy measure $\nu$ has finite intensity, $\displaystyle\Gamma:=\int_{{\mathbb{R}}}\nu(dz)<\infty.$ Let $K<\infty$ be the total size of the oil field at the beginning of the lease, and let $Y(t)$ be the size of the remaining reserve of the oil field by time $t$, obviously $Y(t)\in[0,K]$. We model the evolution of the profit sharing agreement as a differential game where the two players are the oil company and the government. The state variables of our differential game are $X(t)\in{\mathbb{R}}$ and $Y(t)\in[0,K]$, and the state space is ${\mathbb{R}}\times[0,K]$, because the oil price $S_{t}$ is fully determined by its logarithm $X_{t}=\log(S_{t})\in{\mathbb{R}}$. We will refer the oil company as Player 1 and the government as Player 2. We assume that the processes $X(t)$ and $Y(t)$ follow the dynamics

[TABLE]

where $u_{1}(t)\in U_{1}=[0,\bar{u}_{1}]$ is the extraction rate chosen by the company and $u_{2}(t)\in U_{2}=[0,\bar{u}_{2}]$ is the tax rate chosen by the government. The constant $\bar{u}_{1}$ represents the maximum extraction rate and $\bar{u}_{2}$ is the maximum tax rate. The processes $u_{1}(t)$ and $u_{2}(t)$ are control variables, and $W(t)$ is the Wiener process defined on a probability space $(\Omega,{\cal F},P)$. Moreover, we assume that $W(t)$, $\alpha(t)$ and $\eta_{t}$ are independent. The parameter $\exp(\mu(\cdot))$ represents the equilibrium price of oil and $\kappa$ represents the coefficient of mean reversion. For each state $i\in\{1,2\}$ of the oil market, we assume that the corresponding equilibrium price $\exp(\mu(i))$ is known. Similarly $\sigma(\cdot)$ represents the volatility and $\gamma(\cdot)$ represents and the intensity of the jump diffusion. For each state $i\in\{1,2\}$ of the oil market, we assume that $\sigma(i),\,\gamma(i)$ are known nonzero constants. As a matter of fact, $\gamma(i)$ captures the frequencies and jump sizes of the oil price.

Definition 2.1**.**

The extraction and taxation rates $u_{1}(\cdot)$ and $u_{2}(\cdot)$, taking values on intervals $[0,\bar{u}_{1}]$ and $[0,\bar{u}_{2}]$ respectively, are called admissible controls with respect to the initial data $(s,x,y,i)\in{\cal D}:=[0,T]\times{\mathbb{R}}\times[0,K]\times{\cal M}$ if:

• •

Equation (6) has a unique solution with $X(s)=x$, $Y(s)=y$, $\alpha(s)=i$, and $X(t)\in{\mathbb{R}},\,\,Y(t)\in[0,K]$ for all $t\in[0,T]$.

• •

The processes $u_{1}(\cdot)$ and $u_{2}(\cdot)$ are $\{{\cal F}_{t}\}_{t\geq 0}$-adapted where ${\cal F}_{t}=\sigma\{\alpha(s),W(s),\eta_{s};s\leq t\}$.

We use ${\cal U}_{j}={\cal U}_{j}(s,x,y,i)$ to denote the set of admissible controls taking values in $U_{j}=[0,\bar{u}_{j}]$ such that $X(s)=x$, $Y(s)=y$, $\alpha(s)=i)$, for each $j\in\{1,2\}$.

Let $C(t,Y(t),u_{1}(t))$ be the extraction cost function per unit of time $t$. A typical example of extraction cost function is $C(t,y,u):=a+mu(c-by),$ where $a>0$ can be seen as the initial cost of setting up the oil field and $m$, $b$, and $c$ are constants such that $m>0$ and $b,c\geq 0$. The total profit rate for operating the mine is

[TABLE]

The total income tax the government levies on the oil company is $u_{2}(t)\theta P(t,X(t),Y(t),u_{1}(t))$. The post-tax profit rate of the company is

[TABLE]

and the government profit rate function is

[TABLE]

We assume that at the end of the lease there are no extraction revenues, therefore the profit rate for the oil company could be zero or equal to the cost of closing the mine. We will denote by $\Phi_{1}(X(T),Y(T))$ the terminal profit rate for the oil company. In most cases, $\Phi_{1}(x,y):=0$. However, the terminal profit rate of the government is the market value of the remaining reserve, $\Phi_{2}(X(T),Y(T))=Y(T)(\exp(X(T))-m)$, where $m$ is the cost of extracting one barrel. In sum, we will generally assume that the running profit rate and terminal profit rate functions $L_{i},\Phi_{i}$, $i=1,2$ are Lipschitz continuous on bounded sets. Given a discount rate $r>0$, the payoff functional of Player $i=1,2$ is

[TABLE]

Each player wants to maximize its own payoff. The company will try to maximize its payoff by adjusting the extraction rate $u_{1}(\cdot)$, while the government will maximize its payoff by changing the tax rate $u_{2}(\cdot)$. We, therefore, model this interaction as a noncooperative differential game. Our goal is to find a noncooperative Nash equilibrium $(u_{1}^{*},u_{2}^{*})$ such that

[TABLE]

In the next section we will prove the existence of a Nash equilibrium.

3 Nash Equilibrium

Definition 3.1**.**

Let $(u_{1}^{*},u_{2}^{*})$ be a Nash equilibrium of our differential game, the functions

[TABLE]

are called value functions of Player 1 and Player 2 respectively.

In order to find the optimal strategies $u_{1}^{*}$ and $u_{2}^{*}$ of a Nash equilibrium we first have to derive the value functions $V_{1}$ and $V_{2}$ of the differential game. Formally the value functions $V_{1}$ and $V_{2}$ should satisfy the following Hamilton-Jacobi-Isaacs equations. Assuming that we have a Nash equilibrium $(u_{1}^{*},u_{2}^{*})$ let us define the corresponding Hamiltonians:

[TABLE]

and

[TABLE]

with $\displaystyle QV(s,x,y,\cdot)(\iota)=\sum_{j\neq\iota}q_{\iota j}(V(s,x,y,j)-V(s,x,y,\iota))$. The corresponding Hamilton Jacobi Isaacs equations of this noncooperative game are

[TABLE]

The strategy for solving this differential game is to first find the solutions of the Isaacs equation (13) and then derive the optimal extraction and taxation policies from the Nash equilibrium. The next result gives the road map we will use to find a Nash equilibrium if we already have the value functions.

Theorem 3.2**.**

Assume that there exists $(u_{1}^{*},u_{2}^{*})\in{\cal U}_{1}\times{\cal U}_{2}$ such that the nonlinear Hamilton-Jacobi-Isaacs equations (13) have classical solutions $V_{i}(s,x,y,\iota)$, $i=1,2$ with

[TABLE]

and

[TABLE]

Then the pair $(u^{*}_{1},u^{*}_{2})$ is a Nash equilibrium and $J_{i}(s,x,y,\iota;u_{1}^{*},u_{2}^{*})=V_{i}(s,x,y,\iota)$, $i=1,2$.

Proof.

The proof relies on the fact that this problem can be uncoupled and solved as an optimal control problem. In fact, if we replace the control process $u_{2}(\cdot)$ by $u^{*}_{2}(\cdot)$ in (6) then, the differential game problem becomes an optimal control problem with the only control variable $u_{1}(\cdot)$. The HJB equation of this new control problem is

[TABLE]

Following the assumptions of this theorem, it is clear that the HJB equation (16) has a solution $V_{1}$ and the optimal policy of this new control problem is $u_{1}^{*}$. Therefore $u_{1}^{*}$ is in equilibrium with $u_{2}^{*}$ and $V_{1}$ is the value function of Player 1. A similar argument can be used to show that $u_{2}^{*}$ is in equilibrium with $u_{1}^{*}$ and that $V_{2}$ is the value function of Player 2.   \sqcup\hbox to0.0pt{\hss\sqcap}

Given that $L_{i},\Phi_{i},i=1,2$ are Lipschitz continuous, using standard methods from control theory it can be shown that the values functions $V_{1}$ and $V_{2}$ are the unique viscosity solutions of the Isaacs equations (13). The uniqueness of the viscosity solutions are obtained as in Barles and Imbert (2008) by applying nonlocal extensions of the Jensen-Ishii Lemma. For more on the derivation of the maximum principle for nonlocal operators, one can also refer to Biswas et al. (2010).

4 Numerical Approximation

In this section, we construct a finite difference scheme and show that it converges to the unique viscosity solutions of the Isaacs equation (13). We will use the following notations; we set $u=(u_{1},u_{2})$, $u^{*}=(u_{1}^{*},u_{2}^{*})$, $U=U_{1}\times U_{2}$, and

[TABLE]

[TABLE]

[TABLE]

The Isaacs equation (13) can be rewritten as follows

[TABLE]

Let $k\in(0,1)$ be the step size with respect to $s$, and $h\in(0,1)$ be the step size with respect to $x$ and $y$, we will use the standard finite difference operators $\Delta_{s}$, $\Delta_{x}$, $\Delta_{xx}$ and $\Delta_{y}$. Let ${\rm I}f$ denote the integral part of the Hamiltonians $H_{1}$ and $H_{2}$. We will approximate ${\rm I}f$ using the Simpson’s quadrature. Using the fact the Lévy measure is finite $\displaystyle\Gamma=\int_{{\mathbb{R}}}\nu(dz)<\infty$, we have

[TABLE]

We use the Simpson’s quadrature to approximate the integral part of the Hamiltonians. Let $\xi\in(0,1)$ be the step size of the Simpson’s quadrature, the corresponding approximation of the the integral part is

[TABLE]

where the $(c_{j})_{0\leq j\leq N_{\xi}}$ and $(d_{j})_{0\leq j\leq M_{\xi}}$ are the corresponding sequences of the coefficients of the Simpson’s quadrature and $z_{j}\in[-1,1]$, $j\in\{1,...,M_{\xi}\}$ are grid points of the interval $[-1,1]$. In fact, $\displaystyle\lim_{N_{\xi}\rightarrow\infty}\sum_{j=0}^{N_{\xi}}c_{j}=\Gamma$ and $\displaystyle\lim_{M_{\xi}\rightarrow\infty}\sum_{j=0}^{M_{\xi}}d_{j}=\int_{-1}^{1}\nu(dz)$. The corresponding discrete versions of the Hamiltonians $H_{1},H_{2}$ are defined as follows

[TABLE]

and

[TABLE]

Therefore the discrete version of (19) is

[TABLE]

with

[TABLE]

We have the following crucial Lemma.

Lemma 4.1**.**

Let $h\in(0,1)$ be small enough, for each $k,\xi\in(0,1)$, there exists a unique bounded function $V_{h,k,\xi}$ defined on ${\cal D}$ that solves equation (22).

Proof

Any solution $V$ of (22) should also satisfy $V=V+\epsilon(\Delta_{s}-H_{u_{1}^{*},u_{2}^{*}}^{h,k,\xi}V(s,x,y,i))$ for $\epsilon>0$ together with terminal condition $V(T,x,y,i)=\Phi(x,y)$. We define the operator ${\cal F}_{\xi}$ on bounded functions on ${\cal D}$ as follows

[TABLE]

where the coefficients $q_{in}$ are coefficients of the generator $Q$ and the quantities $a_{u,u^{*}}(i)$, $b(i)$ , $c_{u,u^{*}}(i)$ and $PV(s,x,y,\cdot)(i)$ are defined as follows

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

Note that equation (22) is equivalent to $V(s,x,y,i)={\cal F}_{\xi}(V)(s,x,y,i;h,k)$, it suffices to show the operator ${\cal F}_{\xi}$ has a fixed point. As a matter of fact, for $h$ small enough, it is clear that $g(i,v)>0$ for all $i\in{\cal M}$ and $v\in U_{1}$, thus both diagonal coefficients of $a_{u,u^{*}}(i)$ are positive. Moreover, we can choose $\epsilon$ such that $e(v,i)>0,$ for all $v\in U_{1},i\in{\cal M}$ and

[TABLE]

where $\delta$ is a constant. Using the fact that the difference of two suprema is less than the supremum of the difference. If we have two bounded functions $V,W$ defined on ${\cal D}$, it is clear that

[TABLE]

Therefore, the map ${\cal F}_{\xi}$ is a contraction on the space of bounded functions on ${\cal D}$, using the Banach’s Fixed Point Theorem we conclude the proof of the lemma.   \sqcup\hbox to0.0pt{\hss\sqcap}

Remark 4.2**.**

It is clear from Lemma 4.1 that the numerical scheme obtained from (22) is stable since the solution of the scheme is bounded independently of the step sizes $h,k,\xi\in(0,1)$ and obviously consistent because as the step sizes $h,k,\xi$ go to zero the finite difference operators converge to the actual partial differential operators. We have the following convergence theorem.

Theorem 4.3**.**

Let $V_{h,k,\xi}$ be the solution of the discrete scheme obtained in Lemma 4.1. Then as $(h,k,\xi)\rightarrow 0$ the sequence $V_{h,k,\xi}$ converges locally uniformly on ${\cal D}$ to the unique viscosity solution $V$ of (19).

This result is the standard method for approximating viscosity solutions, for more one can refer to Barles and Souganidis (1991).

5 Applications

Consider an oil company with a 10 years lease to extract oil from an oil field with an known capacity of K=10 billion barrels. We assume that the profit sharing agreement between the oil company and the government is such that the oil company takes 40% of profits and the government takes 60%, so $\theta=0.4$. The oil equilibrium price when the market is up is $\mu(1)=50$ and when the market is down is $\mu(2)=35$. The mean reversion coefficient is $\kappa=0.01$, the volatility when the market is up is $\sigma(1)=0.1$ and when the market is down $\sigma(2)=0.3$. And the jump intensity is $\lambda(1)=0.01$ when the market is up and $\lambda(2)=0.15$ when the market is down. The generator of the Markov chain is \displaystyle Q=\left(\begin{array}[]{ll}-0.005&0.005\\ 0.002&-0.002\end{array}\right). We assume that $L_{1}(t,x,y,u_{1},u_{2})=0.4(e^{x}u_{1}-(10+15u_{1}))(1-u_{2})$ , $\Phi_{1}(x,y)=0$, $L_{2}(t,x,y,u_{1},u_{2})=0.6(e^{x}u_{1}-(10+15u_{1}))+0.4u_{2}(e^{x}u_{1}-(10+15u_{1})),$ and $\Phi_{2}(x,y)=y(e^{x}-15)$. Moreover, we assume that the extraction $u_{1}(\cdot)\in[0,50000]$ and the top tax rate is 30% so $u_{2}(\cdot)\in[0,0.3]$. Keep in mind that, because the payoff rates are linear functions of each control variable $u_{1}(\cdot)$ and $u_{2}(\cdot)$, therefore using Theorem 3.2, the optimal strategies $u_{1}^{*}$ and $u_{2}^{*}$ are obtained by looking at the signs of the following functionals $\ F(t,x,y,u_{2},i)=-\frac{\partial V_{1}(t,x,y,i)}{\partial y}+0.4(e^{x}-15))(1-u_{2}(t))$, and $G(t,x,y,u_{1},i)=0.4(e^{x}u_{1}(t)-(10+15u_{1}(t))).$ The optimal strategies will only be attained at the endpoints of the intervals $U_{1}=[0,50000]$ and $U_{2}=[0,0.3]$, we have

[TABLE]

and

[TABLE]

In the next two figures we have the plots of the functions $F$ and $G$ when the market is up and when the market is down. For instance, in Figure 1 we have plots of the $F$ and $G$ when the market is bullish, the regions above the curves represent the domains where signs of $F$ and $G$ are positive, so in those regions it is always optimal to extract at full capacity or to tax at the maximal rate. And the regions below the curves represent the domains where the signs of $F$ and $G$ are negative, thus in those regions it is optimal not to extract at all and not to tax the oil company. In Figure 2, similar plots are given when the market is bearish.