MODELING AND COMPUTATION OF MEAN FIELD GAME WITH COMPOUND CARBON ABATEMENT MECHANISMS

. In this paper, we present a mean ﬁeld game to model the impact of the coexistence mechanism of carbon tax and carbon trading (we call it compound carbon abatement mechanism) on the production behaviors for a large number of producers. The game’s equilibrium can be presented by a system which is composed of a forward Kolmogorov equation and a backward Hamilton-Jacobi-Bellman (HJB) partial diﬀerential equation. An implicit and fractional step ﬁnite diﬀerence method is proposed to discretize the resulting partial diﬀerential equations, and a strictly positive solution is obtained for a non-negative initial data. The eﬃciency and the usefulness of this method are illustrated through the numerical experiments. The sensitivity analysis of the parameters is also carried out. The results show that an agent under concentrated carbon emissions tends to choose emission levels diﬀerent from other agents, and the choices of agents with uniformly distributed emission level will be similar to their initial distribution. Finally, we ﬁnd that for the compound carbon abatement mechanism carbon tax has a greater impact on the permitted emission rights than carbon trading price does, while carbon trading price has a greater impact on carbon emissions than carbon tax.


1.
Introduction. Global climate change is one of the most significant environmental issue of our time. Some studies have shown that if we fail to adequately reduce greenhouse gas emissions in a decade or two, we may face major disasters by the end of the century. Several carbon dioxide reduction policies have been proposed and implemented to address the challenges of climate change. Among these policies, cap-and-trade and carbon tax mechanisms are two of the most popular ones, which use the market-based policy instrument. Currently, carbon taxes are implemented in Denmark, Finland, Italy, the Netherlands, Norway, and Sweden in line with the requirements of the Kyoto Protocol [6,29,8,4]. Simultaneously, the cap-and-trade mechanisms for carbon dioxide emissions have been implemented by the European Union. Several emission permit trading markets have emerged and continued to develop successfully. Extensive studies on carbon dioxide reduction policies are available in the literature, see, for example, [36,9,12,20,32].
The cap-and-trade and carbon tax mechanisms have different characteristics both in practice and in theory. First, by imposing an overall cap on the allowable level of permitted carbon emissions, the cap-and-trace policy can guarantee environmental benefits, but it also shows uncertainty in economic costs. In contrast, the carbon tax policy provides accurate economic costs, but shows ambiguity about environmental impacts. Second, the focus of these two policies is different. The main focus of the cap and trade policy is the distribution and trading of carbon dioxide allowances, while the carbon tax policy focuses on generating revenue. Third, the total allowable emissions are subject to the upper bound constraint of the cap-andtrade mechanism. This feature will disappear in the carbon tax scheme, because enterprisers can increase their emissions to any level they want by paying enough taxes [19,37,39].
Some studies focus on the mixed scheme of the two systems. In [31], it is proposed to reduce social costs by combining these two policies (i.e., quantity-based-approach and price-based-approach ). In [34] based on the basis of command-control mechanisms, it is proposed how to combine tradable permits with other policy instruments, such as taxes and subsidy mechanisms. In [30], it is found that the mixed path mechanism is slightly better than the price mechanism. However, it is better than the quantity mechanism. In [35], it is suggested that the policy mix should be carefully designed and transparent. The coexisting instruments may increase the overall cost of abatement, but no further contribution to emission reductions. Therefore the objectives and trade-offs in the policy mix must be clear. In [28], the focus is on the situation where a part of the economy is regulated by cap-and-trade program, while the rest is subject to emissions taxes. Results show that under optimal conditions, the size of the sector being taxed, increases with the relative steepness of the total aggregate marginal abatement cost function. In [38], it is suggested that carbon trading and carbon tax mechanisms should be combined to effectively regulate China's environmental problems. Based on the analysis of China's energy, economy and environment by using dynamic CGE model, the finding in [33] suggests that the two emission reduction mechanisms should be combined. The results show that carbon taxes should be implemented in decentralized enterprises and carbon trading should be implement in large-scale highly centralized enterprises.
Inspired by the ideas reported in [33,10], we use the mean field game model to simulate a system which has a large number of anonymous and homogeneous enterprises. These enterprises adjust their behaviors according to the mixed carbon abatement polices. We choose the mean field game model based on the following reasons. First, the fact that economic behavior, includes competitions and cooperation, has made the game theory framework popular in recent years, especially in complex economic systems. Second, in all game models, the mean field game model studies the indistinguishable interactions between individuals and the surrounding population [24,25,26]. It is designed to address problems involving a large number of agents. This is because carbon abatement policies affect almost every country, every enterprise, and even every one in the world. Third, considering the randomness of emissions makes the model more realistic for the emission permit trading market [9,12,20,32,11].
Here, we assume that each agent in the mean field game uses a tax mechanism to pay for carbon dioxide lower than a fixed permission level. The excess amount of the permitted level comes from the second mechanism, namely the permission trade. Because the trading cost is higher than the tax rate, agents have the incentive to negotiate for higher permitted levels. Furthermore, we use the dynamic processes to describe actual and permitted emissions. By setting up the revenues of all agents, we study production, negotiation and trading behavior within the framework of the mean field game. In addition, our model extended the work reported in [10], which only considers the permission trading policy.
Theory of mean field game (MFG for short) was first reported in [24,25,26]. MFG models have been applied to a wide range of disciplines, such as econmics, physics, biology and network engineering see [1,15,14,16,21,27,40,17,18]). In this paper, we use the standard structure of the MFG, which includes the expected collective revenue function of all agents, and the emissions and permitted levels that are modeled as stochastic processes. By using control theory and maximum principle, the game can be converted into systems composting of two partial differential equations (PDEs) ( [7]). The first PDE is the Hamilton-Jacobi-Bellman equation, which returns the optimal individual response to the population mean-field emission behavior. The second PDE is the Fokker-Planck equation which describes the density evolution of the players as a results of the implementation of the individual optimal production strategies [24]. It does not appear possible to derive the explicit solution of this problem due to the complex nature of these two PDEs. Therefore, numerical methods will be developed for solving the problem. There exist some numerical methods in the literature for solving mean field game equations, such as those proposed in [2,3] based on the finite differential scheme and and those in [10] based on the finite volume scheme. In this paper, we adopt the method proposed by [23], in which the fractional step finite difference method in time approximation and the finite element method in space are employed. This paper has three main contributions. First, a mean field game is introduced to represent the compound emission reduction strategy. Second, the results show that, due to negative externalities, an agent under concentrated carbon emissions tends to choose different emission levels than other agents. On the other hand, the choices of agents with uniformly distributed emission level will be similar to their initial distribution. Third, we find that carbon tax has a greater impact on the permitted emission rights than carbon trading price, while carbon trading price has a greater impact on the mixed carbon abatement polices than carbon tax.
The remainder of this paper is organized as follows. We present the model in Section 2. The numerical methods and algorithm are described in Section 3. Section 4 discusses the efficiency and usefulness of our method in several numerical experiments and the effects of parameters on population density. Finally, some concluding remarks are presented in Section 5.
2. Model. Our model is defined on the time interval [0, T ], where T denotes the maturity of the game. The agents of the game are the producers of emissions. Each producer is associated with a production at time t denoted by q(t). On the other hand, q(t) produces an amount of emissions E(t). Generally, an increase in production results in additional emissions and vice versa. Thus, emission is treated as a state variable rather than production [10,11]. Its dynamics is described by the following stochastic differential equation.
where l(t) is the level of emission reduction and σ 1 dW 1 (t) is the stochastic disturbance of emission, which is due to technological innovations and market fluctuations (or other uncertain factors). Here, dW 1 (t) is a standard Brownian motion and σ 1 is a noise parameter. The symbol dN t (E(t)) is the reflection part which guarantees that the process remains in [E min , E max ]. It is determined by the production capacity of an agent [22], [13]. Another dynamic state of the system is the amount of the permitted emission x(t). The agents are assumed to exert their best efforts to negotiate with each other to seek a high level of permitted emission. The effort level is represented by µ(t). The permitted emission x(t) evolves according to the stochastic process given below: Here, σ 2 dW 2 (t) is the stochastic disturbance in which dW 1 (t) and dW 2 (t) are independent of each other, dN t (x(t)) is the reflection part which guarantees that the process remains in [x min , x max ]. Note that dN t (x(t)) is determined by the policymakers of the cap-and-trade mechanism.
In the proposed mean field game model, all agents are assumed to show the same ability in production and negotiation but make different decisions. The states E and x of the producers at time t are distributed in the interval [E min , E max ]×[x min , x max ], following the probability density function m(t, E, x). The initial density denoted by m 0 (E, x) is given. A given agent cannot influence the distribution of all player states, but can generate information, thereby affecting the decisions made by other agents.
The revenue function of the agent consists of five parts: the production revenue, the cost of emission reduction, the cost of purchasing permitted emissions, the cost of the carbon tax and carbon trading.
The production revenue of the agent is a quadratic function of the amount of the emissions. It is strictly concave. The form of this function can be expressed by where c 1 , and c 2 are positive constants. The marginal revenue of the agent is positive and decreasing. Due to the density m(t, E, p), the production revenue is decreasing. This property is related to the economic concept, known as "negative externality". That is, an agent should face more intense competition and receive fewer revenues if this agent selects the same state as that of other agents.
The costs of emission reduction and negotiation for the permitted emission are given, respectively, by Here, d 1 and d 2 are positive constants. The quadratic form of these cost functions guarantees that the marginal costs are increasing.
The cost of carbon tax is computed as follows: where P a (t) denotes the tax rate at time t, and x is the permission level of carbon emissions. If the emission amount is less than x, then the agent should pay the carbon tax of his emission. If the emission amount is more than x, then the agent should pay the carbon tax of the basic portion. The permission of the excess portion should be brought by the carbon trading mechanism, as expressed by where P b (t) denotes the price of the emission permits at time t and P b (t) > P a (t).
In summary, the net revenue of the agent at time t is The agents can adjust their strategies l(t) and µ(t) to maximize the expected discounted revenue over the time horizon [0, T ]. When the probability density function is independent of the control variable, the optimal control problem involving the agent is as follows: subject to where r is the risk-free discount rate. The evolution of the law of stochastic process indicates that m(t, E, x) is related to l(t) and µ(t). The density function m(t, E, x) satisfies the forward Kolmogorov equation as given below: with Neumann boundary condition and the initial condition m(0, E, x) = m 0 (E, x). See [15]. Clearly, this problem is described in the framework of mean field game theory.
Then the problem is expressed as follows: subject to the partial differential equation (4) with the Neumman boundary (5) and the initial condition m(0, E, x) = m 0 (E, x).
To derive the optimality conditions, set Using the first order condition, we can obtain with Neumann boundary condition (5).
c means a function set in which the derivatives of any order are smooth and the boundary values of the derivatives of any order are zero). Clearly, this system consists of the backward HJB and forward Kolmogorov equations. The solution to this system is the mean field equilibrium of our game. The existence and uniqueness of its solution is given in the following theorem. Proof. The proof is based on Theorem 2.5 in [26]. We only have to verify that for any m 1 ≥ 0 and m 2 ≥ 0, the condition Since E max E − 1 2 E 2 > 0, c i > 0 and m i ≥ 0, i = 1, 2, the conclusion follows readily.

Numerical experiments.
3.1. Approximation of the state equation. In this section, we present a numerical method to solve the above two coupled partial differential equations appeared in the previous section. LetΩ = ε∈T h ε be a decomposition of Ω = [E min , E max ] × [x min , x max ] into a family T h of nonoverlapping rectangular h E × h x finite elements. Moreover, V h ⊂ H 1 (Ω) is the finite element space of the continuous and piecewise bilinear functions. ByV 0h we denote the subspace of V h × V h of the vector-functionsᾱ h such thatᾱ wheren = (n 1 , n 2 ) is unit outword normal vector to ∂Ω. We use several quadrature formulas for approximating the integrals of a continuous function g(x). Namely, be the vertices of a finite element ε, respectively. We approximate the components of the form a α (m, ν) by the following bilinear forms of m h and v h ∈ V h depending onᾱ ∈V 0h : are, respectively, the trapezoidal quadrature formula and right and left rectangular quadrature formulas in E-direction and x-direction. Here α + Eh and α + xh are positive parts, α − Eh and α − xh are negative parts of corresponding mesh functions. Let ω t = {t j = jτ, j = 0, · · · , M, M τ = T } be the uniform mesh in time on the segment [0, T ] and ω 0t = ω t /{t = 0}. For a function u(t) : ω t → V h we denote by u k = u(t k ) ∈ V h its value on a time level k.
We consider the implicit finite difference scheme of the state equation. Let the initial value m 0 = m 0 ∈ V h and a mesh control functions l and µ ∈V 0h be given. We search m(t) : ω t → V h , which satisfies for all υ ∈ V h and k = 1, 2, · · · , M : We approximate the objective functional J(m, l, µ) by using the quadrature formulas in space and the time variables: For the mesh functions m k ∈ V h and l k , The gradient ∇ m k J h (m, l, µ) with respect to m k is defined by the equality In the case of implicit approximation the adjoint problem is given by the following equation: where g k = f k + m k ∂f k ∂m and the initial condition v M +1 = 0. Proof. See Theorem 2 in [23].
We use the iterative method given in [23] to solve m, v, l and µ. The algorithm is presented as an algorithm given below. We set n as an iteration number.

Algorithm
Step 1: An initial estimate is provided for l and µ, and which are set as l 0 and µ 0 . A tolerance threshold is set, where T ol > 0 and n = 0; Step 2: Equation (8) is solved to obtain m n .
Step 3: Equation (9) is solved to obtain v n ; Step 4: v n is used to compute l n and µ n by (6) in the discretized form. Compute n = µ n+1 − µ n + l n+1 − l n ; Step 5: If n ≤ T ol, then let l = l n , µ = µ n m = m n , v = v n , and stop. Otherwise, let n = n + 1, and return to Step 2.

Numerical results.
To solve our game model numerically, the following parameter values are used.
Parameters: T = 1, E min = 1, E max = 5, x min = 1, x max = 3, c 1 = 10, c 2 = 0.1, First, we consider the convergence rate of our discretization method to show its accuracy and efficiency. Owing to the limitation of space, we only test the adjoint state v. In addition, we consider the solution obtained by selecting M E = M x = 48 in space and M τ = 512 in time as the "exact" solution for v because the closedform solution of Equation (6) cannot be determined. We compute the errors in the discrete L ∞ -norm at the initial computational time step t = 0 on a sequence of meshes with M E = M x = 2 n and M τ = 2 n+1 for a positive integer n from n = 2 to a maximum n = 5. The discrete L ∞ -norm is defined as follows: where v h denotes the numerical solution. The log-log plots of the computed maximum errors and the linear fitting are depicted in Figure 1. This figure illustrates that the rate of convergence of v h in the discrete L ∞ norm is of the order O(h 0.4034 ). This figure also displaces the numerical solution of equation (6). Numerical studies show the evolution of different distributions as summarized in Figures 2, 3 and 4. The states e and x of the continuum of agents at time t are distributed in the interval [E min , E max ]×[x min , x max ]. In the first example the probability density function m(t, E, x) follows a normal distribution. The distribution considered in the second example is as in [5], where the probability density function is For the third example, we consider uniform distribution. The evolution of density functions at time t = 0, t = 1/3, t = 2/3 and t = 1 are given in Figures 2, 3 and 4, respectively.
In the first case, the density function of agents is concentrated initially at the medium levels, and dispersed evenly on the intervals as time progress. In the second case, the density function focuses initially on the left, and over time, it gradually becomes flattened within the permitted emission range. The third case involves a uniform distribution, the density function changes slightly over.  In term of emission concentration, we can see from Figure 2 and Figure 3 that each agent tends to choose different emission levels from other agents due to negative externality. In addition, we find that with the implementation of the compound strategy, the number of low-emission agents will increase. These changes suggest that the compound emission reduction strategy can effectively promote the emission reduction when agents are centralized. However, it has little impact on small and medium-sized agents when they are decentralized.
In this paper, we mainly focus on the impact of compound strategy on optimal production behaviors of a very large number of producers. We provide a numerical experiment to show the influence of carbon tax and carbon trading price on the optimal strategies of producers. We introduce the marginal density functions, as follows: To compensate the errors of the computation, we use the adjusted marginal density functionsm 1 (t, E) andm 2 (t, x), expressed bỹ   The evolutions of average emissionsĒ(E) and permitted emissionĒ(x) (in the sense of expectation) with carbon tax and carbon trading price are given in Figure  5. We can see thatĒ(E) decreases with the increases of carbon tax P a and carbon trading price P b , whileĒ(x) decrease with the increase of carbon tax P a and increases with the increase of carbon trading price P b . In fact, the increase of carbon tax will increase the producer's production cost, thus they will reduce emissions to save costs and increase profits. At the same time, the increase of carbon tax will entail additional cost, which will weaken the producer's desire to pursue high emission rights and henceĒ(x) will also decline. On the contrary, with the increase of carbon trading price P b , the producer's desire to pursue high emission rights becomes stronger under the condition of P b > P a , soĒ(x) will rises. In addition, Figure 5 also revealed that carbon tax has a greater impact on the permitted emission rights, and carbon trading price has a greater impact on carbon emissions in compound mechanism.
From Figures 4 and 5, it can be seen that for areas where carbon emissions are concentrated (take some industrial areas as an example), carbon emissions will be dispersed under the compound emission reduction mechanism, and the total carbon emissions will be reduced at the same time. For areas where carbon emission are dispersed, carbon emissions have a minor change under the compound emission reduction mechanism. On the other hand, we find that carbon tax and carbon trading strategy influence each other under the compound emission reduction mechanism, increasing carbon tax will reduce agents desire for carbon emission rights, thus reduce the total carbon emissions, and increasing carbon trading prices will prompt  agents to reduce carbon emissions. The compound strategy can make up for each other's shortcomings of a single emission reduction mechanism, and it is an effective emission reduction strategy compared with a single carbon tax or carbon trading mechanism.
4. Concluding remarks. We present a mean field game model to study the behavior of agents in a compound carbon abatement scheme. Our model consists of a large munber of agents for which each producer is homogeneous. The equilibrium of this game can be represented by a system composed of a forward Kolmogorov and backward HJB equations. Then an implicit and fractional step finite difference method is proposed to discretize the resulting partial differential equations, a strictly positive solution is obtained for a non-negative initial data. The efficiency and the usefulness of this method are illustrated in several numerical experiments. The sensitivity analysis of the parameter was also presented. It is found that, each agent with concentrated carbon emissions level tends to choose a different emission levels from other agents because of negative externality, agents with uniformly distributed emission level will chose similar to their initial distribution. Finally, we find that under the compound mechanism, carbon tax has a greater impact on the permitted emission rights, while carbon trading price has a greater impact on carbon emissions.