René Aïd

Pierre Gruet

Huyên Pham

An optimal trading problem in intraday electricity markets

Keywords: JEL Classification: G11, Q02, Q40 MSC Classification: 35Q93, 49J20, 60H30, 91G80 Optimal trading, power plant, intraday electricity markets, renewable energy, market impact, linear-quadratic control problem, jumps, delay

come

Introduction

The development of renewable energy sources in Europe as a response to global climate change has led to an increase of exchange in the intraday electricity markets. For instance, the exchanged volume on the European Energy Exchange (EEX) for Germany has grown from 2 TWh in 2008 to 25 TWh in 2013. This increase is mainly due to the level of forecasting error of wind production, which leads power producers owning a large share of wind production to turn more than ever to intraday markets in order to adjust their position and avoid penalties for their imbalances. The accuracy of forecasts for renewable power production from wind and solar may vary considerably depending on the agreggation level (local vs regional forecast) and the time horizon. For a complete survey on this problem, the reader can consult Giebel et al. [START_REF] Giebel | The state-of-the-art in short-term prediction of wind power. A literature overview[END_REF], and may have in mind that the root mean square error (RMSE) of the error forecast for the production of a wind farm in six hours can reach 20% of its installed capacity. Many different intraday markets have been designed and are subject to different sets of regulation. But, in all cases, intraday markets offer power producer the possibility to buy or sell power for the next (say) 9 hours to 32 hours (case of the French electricity market of EpexSpot). These trades can occur after the closing of the day-ahead market or during the clearing phase of the day-ahead market.

The problem of trading management in the intraday electricity market for a balancing purpose has already drawn the attention in the literature. Henriot [START_REF] Henriot | Market design with centralised wind power management: handling low-predictability in intraday markets[END_REF] studied the problem of how the intraday market can help a power producer to deal with the wind production error forecast in a stylized discrete time model. In his model, the power producer is a wind producer who is trying to minimize her sourcing cost on the intraday market while maintaining a balance position between her forecast production and her sales. Henriot's model takes into account the impact of the wind power producer on the intraday price with a deterministic inverse demand function, and the intraday price is not a risk factor. The only risk factor comes from the error forecast of the wind production and its auto-correlation. Garnier and Madlener [START_REF] Garnier | Balancing forecast errors in continuous-trade intraday markets[END_REF] studies the trade-off between entering into a deal in the intraday market right now and postponing it in a discrete time decision model where intraday prices follow a geometric Brownian model and wind production error forecast follows an arithmetic Brownian motion. In their framework, the power producer is supposed to have no impact on intraday prices. Liquidity risk is taken into account as a probability of not finding a counter-party at the next trading window.

In this paper, we consider a power producer having at disposal some renewable energy sources (e.g. wind and solar), and thermal plants (e.g. coal, gas, oil, and nuclear sources), and who can buy/sell energy in the intraday markets. Her purpose is to minimize the imbalance cost, i.e. the cost induced by the difference between the demand of her clients minus the electricity produced and traded, plus the production and trading costs. In contrast with thermal power plants whose generation can be controlled, the power generated from renewable sources is subject to non controllable fluctuations or risks (wind speed, weather forecast) and is then considered here as a random factor just like the demand. We then call the residual demand the demand minus the energy generated by renewable energy. Thus, the problem of the power producer is to minimize the imbalance costs arising from her residual demand by relying both on her own controllable thermal assets and on the intraday market. As in [START_REF] Garnier | Balancing forecast errors in continuous-trade intraday markets[END_REF], we assume that the power producer has access to a continuously updated forecast of the residual demand to be satisfied at terminal date T and that this forecast evolves randomly. Moreover, the intraday price for delivery at time T evolves also randomly and is correlated with the residual demand forecast. However, compared with [START_REF] Garnier | Balancing forecast errors in continuous-trade intraday markets[END_REF], the intraday market can be used for optimization purposes. We develop a model that allows us to study how power producers can take advantage of the interaction between the dynamics of the residual demand forecast and the dynamics of the intraday prices.

Our model shares some links with optimal order execution problems, as introduced in the seminal paper by Almgren and Chriss [START_REF] Almgren | Optimal execution of portfolio transactions[END_REF], and then largely studied in the recent literature, see e.g. the survey paper [START_REF] Schied | Some mathematical aspects of market impact modeling[END_REF]. In our context, the original feature with respect to this literature is the consideration of a random demand target and the possibility for the agent to use her thermal power production. This connection with optimal execution is fruitful in the sense that it allows us to take into account several features of intraday markets while maintaining the tractability of the model sufficiently high to allow analytical solutions. Hence, we take into account liquidity risk through a market impact, both permanent and temporary, on the electricity price generated by a power producer when trading in the intraday market. As in optimal execution problems, this impact is always in the adverse direction: when the producer sells, the price decreases and when she buys, the price increases. Our setting is a continuous-time decision problem representing the possibility for the producer to make a deal at each time she wants and not only at pre-specified windows. Moreover, it is general enough as it permits us to study the limiting cases of a pure retailer (no production function), a pure trader (no demand commitment) and an integrated player (player owning both clients and generation), small or large.

The main goal of this paper is to derive analytical results, which provide explicit solutions for the (approximate) optimal control, hence giving enlightening economic interpretations of the optimal trading strategies. In order to achieve such analytical tractability, we have to make some simplifying assumptions on the dynamics of the price process and of the residual demand forecast, as well as on the cost function, assumed to be of quadratic form meaning a simple linear growth of the marginal cost of production with respect to the production level. We first consider a simple model for a continuous price process with linear impact, and demand forecast driven by an arithmetic Brownian motion, and neglect in a first step the delay of production when using thermal power plants. We then study an auxiliary control problem by relaxing the nonnegativity constraint on the generation level, for which we are able to derive explicit solutions. The approximation error induced by this relaxation constraint is analyzed. In next steps, we consider more realistic situations and investigate two extensions: (i) On one hand, we incorporate the case where the residual demand forecast is subject to sudden changes, related to prediction error for wind or solar power production, which may be quite important due to the difficulties for estimating wind speed and forecasting weather, see [START_REF] Bensoussan | A generalized model approach to seasonal aspects of wind speed modeling[END_REF]. This is formalized by jumps in the dynamics of the demand process, and consequently also on the price process. Again, we are able to obtain explicit solutions. Actually, the key tool in the derivation of all these analytical results is a suitable treatment of the linear-quadratic structure of our stochastic control problem. (ii) On the other hand, we introduce natural delay constraints in the production, and show how the optimal decision problem can be explicitly solved by a suitable reduction to a problem without delay.

Our (approximate) optimal trading strategies present some remarkable properties. When the intraday price process is a martingale, the optimal trading rate inherits the martingale property, which implies in particular that the net position of electricity shares has a constant growth rate on average. Moreover, the optimal strategy consists in making at each time the forecast marginal cost equal to the forecast intraday price. This property follows the common sense of intraday traders. Consequently, if the producer has made sales or purchases on the day-ahead such that her forecast marginal cost equals the day-ahead price and if the initial condition of the intraday price is the day-ahead price, thus, on average, the producer optimal trading rate is zero. This fact is no longer true when the demand forecast and the price follow processes with jumps. In this case, the optimal trading rate is a supermartingale or a submartingale depending on the relative probability and size of positive and negative jumps on the price process. For this reason, contrary to the case without jumps, the power producer may need to have a non-zero initial trading rate even if she has made sales or purchases on the day-ahead such that her forecast marginal cost equals the day-ahead price and if the initial condition of the intraday price is the day-ahead price. We also quantify explicitly the impact of delay in production on the trading strategies. When the price process is a martingale, the net inventory in electricity shares grows linearly on average, with a change of slope (which is smaller) at the time decision for the production.

The outline of the paper is organized as follows. We formulate the optimal trading problem in Section 2. In Section 3, we study the optimal trading problem without delay. We first solve explicitly the auxiliary optimal execution problem, and then study the approximation on the solution to the original problem, by focusing in particular on the error asymptotics. We illustrate our results with some numerical tests and simulations. We extend in Section 4 our results to the case where jumps in demand forecast may arise. In Section 5, we show how the optimal trading problem with delay in production can be reduced to a problem without delay, and then leads to explicit solutions. Finally, the appendix collects the explicit derivations of our solutions, which are justified by verification theorems.

Problem formulation

We consider an agent on an intraday energy market, who is required to guarantee her equilibrium supply/demand for a given fixed time T : she has to satisfy the demand of her customers by purchase/sale of energy on the intraday market at time T and also by means of her thermal power generation. We denote by X t the net position of sales/purchases of electricity at time t ≤ T for a delivery at terminal time T , assumed to be described by an absolutely continuous trajectory up to time T , and by q t = Ẋt the trading rate: q t > 0 means an instantaneous purchase of electricity, while q t < 0 represents an instantaneous sale at time t:

X t = X 0 + t 0 q s ds, 0 ≤ t ≤ T. (2.1)
Given the trading rate, the transactions occur with a market price impact:

P t (q) = Pt + t 0 g(q s )ds + f (q t ).
Here, ( Pt ) t is the unaffected intraday electricity price process on a filtered probability space (Ω, F, F = (F t ) t∈[0,T ] , P), carrying some part of randomness of the market, and following the terminology in the seminal paper by Almgren and Chriss [START_REF] Almgren | Optimal execution of portfolio transactions[END_REF], the term f (q t ) refers to the temporary price impact, while t 0 g(q s )ds describes the permanent price impact. The price ( Pt ) t may be seen as a forward price, evolving in real time, for delivery at time T . Let us then denote by Y the intraday electricity price impacted by the past trading rate q of the agent, defined by:

Y t := Pt + t 0 g(q s )ds.
We assume that Y t is observable and quoted, which means actually that the agent is a large trader and electricity producer, whose actions directly impact the intraday electricity price. The case where the agent is a small producer can be also dealt with by simply considering a zero permanent impact function g ≡ 0. Notice that the transacted price is equal to the sum of the quoted price Y and the temporary price impact:

P t (q) = Y t + f (q t ).
(2.

2)

The residual demand D T is the consumption of clients of the agent minus the production from renewable energy at terminal date T , and we assume that the agent has access to a continuously updated forecast (D t ) t of the residual demand. The agent can use her thermal power production with a quantity ξ at cost c(ξ) in order to match as close as possible the target demand D T . In practice, generation of electricity can not be obtained instantaneously and needs a delay to reach a required level of production. Hence, the decision to produce a quantity ξ should be taken at time T -h, where h ∈ [0, T ] is the delay. Thus, for a controlled trading rate q = (q t ) t ∈ A, the set of real-valued F-adapted processes satisfying some integrability conditions to be precised later, a production quantity ξ ∈ L 0 + (F T -h ), the set of nonnegative F T -h -measurable random variables, the total cost is:

T 0 q t P t (q)dt + C(D T -X T , ξ) := T 0 q t P t (q)dt + c(ξ) + η 2 (D T -X T -ξ) 2 . (2.
3)

The first term in (2.3) represents the total running cost arising from the trading in the intraday electricity market, and the last term, where η > 0, represents the quadratic penalization when the net position in sales/purchases of electricity X T + ξ (including the production quantity ξ at cost c(ξ)) at terminal date T does not fit the effective demand D T . The objective of the agent is then to minimize over q and ξ the expected total cost:

minimize over q ∈ A, ξ ∈ L 0 + (F T -h ) E T 0 q t P t (q)dt + C(D T -X T , ξ) . (2.4) Remark 2.1 1)
The penalization term in the objective function above is a simplification of the effective penalization process that can be found in real electricity markets. For instance, the penalization of imbalances in the French electricity market depends both on the sign of the imbalance of the electricity system and on the price of imbalances (see [1, chap 2., sec. 2.2.1]). Nevertheless, the positive coefficient η captures the main objective of the penalization process. The agent has no incentive of being either too long or too short.

2) On real markets, trading ends some time before the date of delivery, at which the agent has to ensure equilibrium (e.g. on the French electricity market, there is a delay of 45 minutes). We do not include that practical fact in our framework, by considering that the delay is null for the sake of clarity. There is no mathematical consequence: it is enough to have in mind that the delivery and production do not really take place at T , but at T plus some delay.

3) The larger is η, the stronger is the incentive for the agent to be as close as possible to the equilibrium supply-demand. At the limit, when η goes to infinity, the agent is formally constrained to fit supply and demand. However, the limiting problem when η = ∞ is not mathematically well-posed since such perfect equilibrium constraint is in general not achievable. Indeed, the demand at terminal date T is random, typically modelled via a Gaussian noise, and the inventory X which is of finite variation, may exceed or underperform with positive probability the demand D T at terminal date T . Hence, in the scenario where X T > D T , and since by nature the production quantity ξ is nonnegative, it is not possible to realize the equilibrium X T + ξ = D T , even if there is no delay. In the sequel, we fix η > 0 (which may be large, but finite), and study the stochastic control problem (2.4).

4)

The optimization problem (2.4) shares somes similarities with the optimal execution problem in limit order book studied in the seminal paper by Almgren and Chriss [START_REF] Almgren | Optimal execution of portfolio transactions[END_REF], and then extended by many authors in the literature, see e.g. the survey paper [START_REF] Schied | Some mathematical aspects of market impact modeling[END_REF]. The main difference is that in the execution problem of equities, the target is to buy or sell a certain number of shares, i.e. lead X T to a fixed constant (meaning formally that η goes to infinity) while in our intraday electricity markets context, the target is to realize the equilibrium with the random demand D T , eventually with the help of production leverage ξ. However, in contrast with the case of constant target, it is not possible in presence of random target D T to achieve perfectly the equilibrium, which justifies the introduction of the penalty factor η as pointed out above.

The main aim of this paper is to provide explicit (or at least approximate explicit) solutions to the optimization problem (2.4), which are easily interpreted from an economic point of view, and also allow to measure the impact of the various parameters of the model. In order to achieve this goal, we shall adapt our modeling as close as possible to the linearquadratic framework of stochastic control, and make the following assumptions: The energy production cost function is in the quadratic form:

c(x) = β 2 x 2 ,
for some β > 0. Although simple, a quadratic cost function represents the increase of the marginal cost of production with the level of production.

Remark 2.2 (Pure retailer) In the limiting case when β goes to infinity, meaning an infinite cost of production, this corresponds to the framework where the agent never uses the production leverage and only trades in the intraday-market by solving the optimal execution problem:

minimize over q ∈ A E T 0 q t P t (q)dt + C(D T -X T , 0) . (2.5)
As in Almgren and Chriss, we assume that the price impact (both permanent and temporary) is of linear form, i.e. g(q) = νq, f (q) = γq, for some constants ν ≥ 0 and γ > 0. The unaffected intraday electricity price is taken as a Bachelier model:

Pt = P0 + σ 0 W t , (2.6)
where W is a standard Brownian motion, and σ 0 > 0 is a positive constant. Such assumption might seem a shortcoming at first sight since it allows for negative values of the unaffected price. However, in practice, for our intraday execution problem within few hours, negative prices occur only with negligible probability. The martingale assumption is also standard in the market impact literature since drift effects can often be ignored due to short trading horizon. The quoted price Y , impacted by the past trading rate q ∈ A, is then governed by the dynamics:

dY t = νq t dt + σ 0 dW t . (2.7)
The dynamics of the residual demand forecast is given by

dD t = µdt + σ d dB t , (2.8)
where µ, σ d are constants, with σ d > 0, and B is a Brownian motion correlated with W :

d < W, B > t = ρdt, ρ ∈ [-1, 1]. From (2.
2), one can then define the value function associated to the dynamic version of the optimal execution problem (2.4) by:

v(t, x, y, d) := inf q∈At,ξ∈L 0 + (F T -h ) J(t, x, y, d; q, ξ) (2.9) with J(t, x, y, d; q, ξ) := E T t q s (Y t,y s + γq s )ds + C(D t,d T -X t,x T , ξ) 2 , (2.10)
for (t, x, y, d) ∈ [0, T ] × R × R × R
, where A t denotes the set of real-valued processes q = (q s ) t≤s≤T s.t. q s is F s -adapted and E T t q 2 s ds] < ∞, D t,d is the solution to (2.8) starting from d at t, and given a control q ∈ A t , Y t,y denotes the solution to (2.7) starting from y at time t, and X t,x is the solution to (2.1) starting from x at t.

In a first step, we shall consider the case when there is no delay in the production, and then show in the last section of this paper how to reduce the problem with delay to a no delay problem. We shall also study the case when there are jumps in the residual demand forecast.

Optimal execution without delay in production

In this section, we consider the case when there is no delay in production, i.e. h = 0. In this case, we notice that the optimization over q and ξ in (2.4) is done separately. Indeed, the production quantity ξ ∈ L 0 + (F T ) is chosen at the final date T , after the decision over the trading rate process (q t ) t∈[0,T ] is achieved (leading to an inventory X T ). It is determined optimally through the optimization a.s. at T of the terminal cost C(D T -X T , ξ), hence in feedback form by ξ * T = ξ+ (D T -X T ) where

ξ+ (d) := arg min ξ≥0 C(d, ξ) = arg min ξ≥0 β 2 ξ 2 + η 2 (d -ξ) 2 = η η + β d1 d≥0 . (3.1)
The value function of problem (2.9) may then be rewritten as

v(t, x, y, d) = inf q∈At E T t q s (Y t,y s + γq s )ds + C + (D t,d T -X t,x T ) 2 , (3.2)
where

C + (d) := C(d, ξ+ (d)) = 1 2 ηβ η + β d 2 1 d≥0 + η 2 d 2 1 d<0 . (3.3)
and the optimal trading rate q * is derived by solving (3.2). Due to the indicator function in C + , caused by the non negativity constraint on the production quantity, there is no hope to get explicit solutions for the problem (3.2), i.e. solve explicitly the associated dynamic programming Hamilton-Jacobi-Bellman (HJB) equation. We shall then consider an auxiliary execution problem by relaxing the sign constraint on the production quantity, for which we are able to provide explicit solution. Next, we shall see how one can derive an approximate solution to the original problem in terms of this auxiliary explicit solution, and we evaluate the error and illustrate the quality of this approximation by numerical tests.

Auxiliary optimal execution problem

We consider the optimal execution problem with relaxation on the non negativity constraint of the production leverage, and thus introduce the auxiliary value function

ṽ(t, x, y, d) := inf q∈A,ξ∈L 0 (F T ) J(t, x, y, d; q, ξ), for (t, x, y, d) ∈ [0, T ] × R × R × R.
By same arguments as for the derivation of (3.2), we have

ṽ(t, x, y, d) = inf q∈A E T t q s (Y t,y s + γq s )ds + C(D t,d T -X t,x T ) , (3.4)
where

ξ(d) := arg min ξ∈R C(d, ξ) = η η + β d, C(d) := C(d, ξ(d)) = 1 2 ηβ η + β d 2 =: 1 2 r(η, β)d 2 . (3.5)
The function in (3.5) can be interpreted as a reduced cost function. Because the production cost function and the penalization are both quadratic, they can be reduced to a single production function where the imbalances are internalized by the producer. By exploiting the linear-quadratic structure of the stochastic control problem (3.4), we can obtain explicit solutions for this auxiliary problem.

Theorem 3.1 The value function to (3.4) is explicitly equal to:

ṽ(t, x, y, d) = r(η, β)( ν 2 (T -t) + γ) (r(η, β) + ν)(T -t) + 2γ (d -x) 2 + 2µ(T -t)(d -x) + T -t (r(η, β) + ν)(T -t) + 2γ - y 2 2 + r(η, β)µ(T -t)y + r(η, β)(T -t) (r(η, β) + ν)(T -t) + 2γ (d -x)y + γ σ 2 0 + σ 2 d r 2 (η, β) -2ρσ 0 σ d r(η, β) r(η, β) + ν 2 ln 1 + (r(η, β) + ν)(T -t) 2γ + σ 2 d r(η, β)ν + 2ρσ 0 σ d r(η, β) -σ 2 0 2 r(η, β) + ν (T -t) + r(η, β)µ 2 (T -t) 2 ( ν 2 (T -t) + γ) (r(η, β) + ν)(T -t) + 2γ , for (t, x, y, d) ∈ [0, T ] × R × R × R
, with an optimal trading rate given in feedback form by:

qs = q T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s , t ≤ s ≤ T q(t, d, y) := r(η, β)(µt + d) -y (r(η, β) + ν)t + 2γ . (3.6)
Here ( Xt,x,y,d , Ŷ t,x,y,d , D t,d ) denotes the solution to (2.1)-(2.7)-(2.8) when using the feedback control q, and starting from (x, y, d) at time t. Finally, the optimal production leverage is given by:

ξT = ξ(D t,d T -Xt,x,y,d T ) = η η + β D t,d T -Xt,x,y,d T . (3.7)
Skech of proof. We look for a candidate solution to (3.4) in the quadratic form:

w(t, x, y, d) = A(T -t)(d -x) 2 + B(T -t)y 2 + F (T -t)(d -x)y + G(T -t)(d -x) + H(T -t)y + K(T -t),
for some deterministic functions A, B, F , G, H and K. Plugging this ansatz into the Hamilton-Jacobi-Bellman (HJB) equation associated to the stochastic control problem (3.4), we find that these deterministic functions should satisfy a system of Riccati equations, which can be explicitly solved. Then, by a classical verification argument, we check that this ansatz w is indeed equal to ṽ, with an optimal feedback control derived from the argmax in the HJB equation. The details of the proof are reported in Appendix.

Remark 3.1 (Pure trader) By sending β to infinity in the expression of the value function ṽ and of the optimal feedback control q, and observing that r(η, β) goes to η, we obtain the solution to the optimal execution problem (2.5) without leverage production:

v N P (t, x, y, d) := inf q∈A E T t q s (Y t,y s + γq s )ds + C(D t,d T -X t,x T , 0) (3.8) = η( ν 2 (T -t) + γ) (η + ν)(T -t) + 2γ (d -x) 2 + 2µ(T -t)(d -x) + T -t (η + ν)(T -t) + 2γ - y 2 2 + ηµ(T -t)y + η(T -t) (η + ν)(T -t) + 2γ (d -x)y + γ σ 2 0 + σ 2 d η 2 -2ρσ 0 σ d η η + ν 2 ln 1 + (η + ν)(T -t) 2γ + σ 2 d ην + 2ρσ 0 σ d η -σ 2 0 2 η + ν (T -t) + ηµ 2 (T -t) 2 ( ν 2 (T -t) + γ) (η + ν)(T -t) + 2γ , for (t, x, y, d) ∈ [0, T ] × R × R × R
, with an optimal trading rate given in feedback form by:

qNP s = qNP T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s , t ≤ s ≤ T qNP (t, d, y) := η(µt + d) -y (η + ν)t + 2γ .
Interpretation:

1. The optimal trading rate qs at time s ∈ [t, T ], given in feedback form by (3.6), is decomposed in two terms: the first one

r(η, β) (r(η, β) + ν)(T -t) + 2γ µ(T -s) + D t,d s -Xt,x,y,d
s is related to the trading rate in order to follow the trend of the demand, and to the incentive to invest when the forecast of the residual demand is larger than the current inventory. The second term

- 1 (r(η, β) + ν)(T -t) + 2γ Ŷ t,x,y,d s
represents the negative impact of the quoted price on the investment strategy: the higher the price is, the more the agent decreases her trading rate until she reaches negative value meaning a resale of electricity shares. These effects are weighted by the constant denominator term depending on the penalty factor η, the marginal cost production factor β, the temporary and permanent price impact parameters γ, ν, and the time to maturity T -t.

2. By introducing the marginal cost function: c ′ (x) = βx, and the process

ξs := η η + β D t,d s + µ(T -s) -Xt,x,y,d s -qs (T -s) , t ≤ s ≤ T,
which is interpreted as the forecast production for the final time T (recall expression (3.7) of the final production), we notice from the expression of the optimal trading rate that the following relation holds:

Ŷ t,x,y,d s + ν qs (T -s) + 2γ qs = c ′ ( ξs ), t ≤ s ≤ T. (3.9)
This relation means that at each time, the optimal trading rate is to make the forecast intraday price plus marginal temporary impact (left hand side), which can be seen as the marginal cost of electricity on the intraday market at time T , equal to the forecast marginal cost of production. Here, the instantaneous impact γ appears as a marginal cost of buying or selling, and the forecast at time s supposes that the optimal trading rate qs is held constant between s and T .

We complete the description of the optimal trading rate by pointing out a remarkable martingale property.

Proposition 3.1 The optimal trading rate process (q s ) t≤s≤T in (3.6) is a martingale.

Proof. By applying Itô's formula to qs

= q(T -s, D t,d s -Xt,x, y,d s , Ŷ t,x,y,d s
), t ≤ s ≤ T , and since q is linear in d and y, we have:

dq s = - ∂ q ∂t + (µ -q) ∂ q ∂d + ν q ∂ q ∂y (T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s )ds + ∂ q ∂d (T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s )σ d dB s + ∂ q ∂y (T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s )σ 0 dW s ,
from the dynamics (2.1), (2.8), and (2.7) of Xt,x,y,d , D t,d and Ŷ t,x,y,d . Now, from the explicit expression of the function q(t, y, d), we see that

- ∂ q ∂t + (µ -q) ∂ q ∂d + ν q ∂ q ∂y = 0,
and so:

dq s = r(η, β)σ d (r(η, β) + ν)(T -s) + 2γ dB s - σ 0 (r(η, β) + ν)(T -s) + 2γ dW s , (3.10)
which shows the required martingale property.

Remark 3.2 Recall that in the classical optimal execution problem as studied in [START_REF] Almgren | Optimal execution of portfolio transactions[END_REF], the optimal trading rate is constant. We retrieve this result in their framework which corresponds to the case where σ d = 0 (constant demand target), β = ∞ (there is no production), and η = ∞ (constraint to lead X T to the fixed target), see Remark 2.1 4). Indeed, in these limiting regimes, we see from (3.10) that dq s = 0, meaning that {q s , t ≤ s ≤ T } is constant. In our framework, this is generalized to the martingale property of the optimal trading rate process, which implies that the optimal inventory { Xt,x,y,d s , t ≤ s ≤ T } has a constant growth rate in mean, i.e. dE[ Xt,x,y,d s ] ds is constant equal to the initial trading rate at time t given by q(T -t, d -x, y).

As a consequence of this martingale property, if the producer already satisfies the relation (3.9) in the day-ahead market, and if the initial intraday price is the day-ahead price, her initial trading rate on the intraday market will be zero. And thus, on average, her trading rate will be zero.

The martingale property of the trading rate process is actually closely related to the martingale dynamics of the unaffected price P in (2.6). As we shall see in Section 4 where we consider jumps on price, making P a sub or super martingale, the optimal trading rate will inherit the converse sub or super martingale property.

Approximate solution

We go back to the original execution problem with the non negativity constraint on the production quantity. As pointed out above, there is no explicit solution in this case, due to the form of the terminal cost function C + . The strategy is then to use the explicit control consisting in the trading rate q derived in (3.6), and of the truncated nonnegative production quantity:

ξ * T := ξT 1 ξT ≥0 = ξ+ (D t,d T -Xt,x,y,d T ), (3.11)
with ξT defined in (3.7) from the auxiliary problem. In other words, we follow the trading rate strategy q determined from the problem without constraint on the final production quantity, and at the terminal date use the production leverage if the final inventory Xt,x,y,d . The aim of this section is to measure the relevance of this approximate strategy (q, ξ * T ) ∈ A × L 0 + (F T ) with respect to the optimal execution problem (2.9) by estimating the induced error:

E 1 (t, x, y, d) := J(t, x, y, d; q, ξ * T ) -v(t, x, y, d), for (t, x, y, d) ∈ [0, T ] × R × R × R.
We also measure the approximation error on the value functions:

E 2 (t, x, y, d) := v(t, x, y, d) -ṽ(t, x, y, d).
Notice that if ξT ≥ 0 a.s., i.e. D t,d T ≥ Xt,x,y,d T a.s. (which is not true), and so ξ * T = ξT , then clearly (q, ξT ) would be the solution to (2.9), and so

E 1 (t, x, y, d) = E 2 (t, x, y, d) = 0.
Actually, these errors depend on the probability of the event: { Xt,x,y,d

T > D t,d
T }, and we have the following estimate:

Proposition 3.2 For all (t, x, y, d) ∈ [0, T ] × R × R × R, we have 0 ≤ E i (t, x, y, d) ≤ ηr(η, β) 2β V (T -t)ψ m(T -t, d -x, y) V (T -t) , i = 1, 2, (3.12)
where

ψ(z) := (z 2 + 1)Φ(-z) -zφ(z), z ∈ R,
with φ = Φ ′ the density of the standard normal distribution, and

m(t, d, y) := (νt + 2γ)(µt + d) + yt (r(η, β) + ν)t + 2γ , (3.13) V (t) := t 0 σ 2 0 s 2 + σ 2 d (νs + 2γ) 2 + 2ρσ 0 σ d s(νs + 2γ) (r(η, β) + ν)s + 2γ 2 ds ≥ 0. (3.14)
Proof. By definition of the value functions v and ṽ, recalling that (q, ξT ) is an optimal control for ṽ, and since (q, ξ * T ) ∈ A × L 0 + (F T ), we have:

J(t, x, y, d; q, ξT ) = ṽ(t, x, y, d) ≤ v(t, x, y, d) ≤ J(t, x, y, d; q, ξ * T ), for all (t, x, y, d) ∈ [0, T ] × R × R × R.
E(t, x, y, d) = E C(D t,d T -Xt,x,y,d T , ξ * T ) -C(D t,d T -Xt,x,y,d T , ξT ) = E C(D t,d T -Xt,x,y,d T , ξ+ (D t,d T -Xt,x,y,d T )) -C(D t,d T -Xt,x,y,d T , ξ(D t,d T -Xt,x,y,d T )) = E C + (D t,d T -Xt,x,y,d T ) -C(D t,d T -Xt,x,y,d T ) = ηr(η, β) 2β E D t,d T -Xt,x,y,d T 2 1 D t,d T -Xt,x,y,d T <0 , (3.15)
from the definitions and expressions of C + and C in (2.3), (3.3) and (3.5). Now, from (3.10) and by integration, we obtain the explicit (path-dependent) form of the optimal trading rate control:

qs = qt + s t r(η, β)σ d (r(η, β) + ν)(T -u) + 2γ dB u - s t σ 0 (r(η, β) + ν)(T -u) + 2γ dW u , t ≤ s ≤ T, with qt = q(T -t, d -x, y).
We then obtain the expression of the final spread between demand and inventory:

D t,d T -Xt,x,y,d T = d -x + µ(T -t) + T t σ d dB s - T t qs ds = m(T -t, d -x, y) + T t σ d (ν(T -s) + 2γ) (r(η, β) + ν)(T -s) + 2γ dB s + T t σ 0 (T -s) (r(η, β) + ν)(T -s) + 2γ
E(t, x, y, d) = ηr(η, β) 2β V (T -t)ψ m(T -t, d -x, y) V (T -t) ,
while the probability that the final inventory is larger than the terminal demand is:

P D t,d T -Xt,x,y,d T < 0 = Φ - m(T -t, d -x, y) V (T -t) . (3.16)
Error asymptotics. We now investigate the accuracy of the upper bound in (3.12)

Ē(T -t, d -x, y) := ηr(η, β) 2β V (T -t)ψ m(T -t, d -x, y) V (T -t) .
It is well-known (see e.g. Section 14.8 in [START_REF] Williams | Probability with martingales[END_REF]) that

zΦ(-z) ≤ φ(z), ∀z ∈ R, (3.17)
from which we easily see that ψ is non increasing, convex, and ψ(∞) = 0. Thus, Ē(Tt, d -x, y) decreases to zero for large m(T -t, d -x; y) or small V (T -t). We shall study its asymptotics in three limiting cases (i) the time to maturity T -t is small, (ii) the initial demand spread d -x is large, (iii) the initial quoted price y is large. We prove that the error bound Ē(T -t, d -x, y), and thus E 1 (t, x, y, d), E 2 (t, x, y, d), converge to zero at least with an exponential rate of convergence in these limiting regimes:

Proposition 3.3 (i) For all (x, y, d) ∈ R × R × R with d > x, we have lim sup T -t↓0 (T -t) ln Ē(T -t, d -x, y) ≤ - 1 2 d -x σ d 2 . (3.18) (ii) For all (t, y) ∈ [0, T ) × R, we have lim sup d-x→∞ 1 (d -x) 2 ln Ē(T -t, d -x, y) ≤ - 1 2 m 2 ∞ (T -t) V (T -t) , (3.19)
where

m ∞ (t) = νt + 2γ r(η, β) + ν t + 2γ (iii) For all (t, x, d) ∈ [0, T ) × R × R, we have lim sup y→∞ 1 y 2 ln Ē(T -t, d -x, y) ≤ - 1 2 n 2 ∞ (T -t) V (T -t) , (3.20)
where

n ∞ (t) = t (r(η, β) + ν)t + 2γ .
Proof. From (3.17), we have:

0 ≤ ψ(z) ≤ z -1 φ(z), ∀z > 0.
Notice that in the three asymptotic regimes (i) (with d -x > 0), (ii), and (iii), the quantity m(T -t, d -x, y) is positive, and we thus have:

Ē(T -t, d -x, y) ≤ ηr(η, β) 2β V (T -t) 3 2 m(T -t, d -x, y) φ m(T -t, d -x, y) V (T -t) . (3.21) (i) For small time to maturity T -t, we see that m(T -t, d -x, y) converges to d -x > 0, while V (T -t) ∼ σ 2 d (T -t), i.e. V (T -t)/σ 2 d (T -t) converges to 1.
This shows from (3.21) that, when T -t goes to zero, the error bound Ē(T -t, d -x, y), converges to zero at least with an exponential rate of convergence, namely the one given by (3.18).

(ii) For large demand spread d -x, we see that

m(T -t, d -x, y) ∼ m ∞ (T -t)(d -x), i.e. the ratio m(T -t, d -x, y)/m ∞ (T -t)(d -x)
converges to 1 when d -x goes to infinity. This shows from (3.21) that, when d -x goes to infinity, the error bound Ē(T -t, d -x, y), converges to zero at least with an exponential rate of convergence, namely the one given by (3.19).

(iii) For large y, we see that m(T -t, d -x, y) ∼ n ∞ (T -t)y, i.e. the ratio m(T -t, dx, y)/n ∞ (T -t)y converges to 1 when y goes to infinity. This shows from (3.21) that, when d -x goes to infinity, the error bound Ē(T -t, d -x, y) converges to zero at least with an exponential rate of convergence, namely the one given by (3.20).

Interpretation. Recall from (3.16) that

P D t,d T < Xt,x,y,d T = Φ - m(T -t, d -x, y) V (T -t) ,
and thus following the same arguments as in the above proof, we have:

(i) lim sup T -t↓0 (T -t) ln P D t,d T < Xt,x,y,d T = - 1 2 d -x σ d 2 , (3.22) for all (x, y, d) ∈ R × R × R with d > x.
We observe that the rate in the rhs of (3.18) (or (3.22)) depends only on the demand volatility σ d and the initial demand spread d -x. Moreover, it is all the larger, the smaller σ d is, and the larger d -x is. This means that the terminal demand will stay with very high probability above the final inventory once we are near from the maturity with a low volatile demand, initially larger than the inventory, in which case, the explicit strategy (q, ξ * T ) approximates very accurately the optimal strategy (q * , ξ * T ).

(ii)

lim sup d-x→∞ 1 (d -x) 2 ln P D t,d T < Xt,x,y,d T = - 1 2 m 2 ∞ (T -t) V (T -t) , (3.23)
for all (t, y) ∈ [0, T ) × R, The rate in the rhs of (3.19) (or (3.23)) is all the larger, the smaller the volatilities σ 0 and σ d of the electricity price and demand are. Again, we have the same interpretation than in the asymptotic regime (i), and this means that the explicit strategy (q, ξ * T ) approximates very accurately the optimal strategy (q * , ξ * T ) in the limiting regime when the initial demand spread is large, and the volatilities are small.

(iii)

lim sup y→∞ 1 y 2 ln P D t,d T < Xt,x,y,d T = - 1 2 n 2 ∞ (T -t) V (T -t) , (3.24)
for all (t, x, d) ∈ [0, T ) × R × R. In the limiting regime where the initial quoted price y is large, the agent has a strong incentive to sell energy on the intraday market, which leads to a final inventory staying under the final demand with high probability, and thus to a very accurate approximate strategy (q, ξ * T ). As in case (ii), this accuracy is strengthened for small volatilities σ 0 and σ d of the electricity price and demand.

Numerical results

Numerical tests

We measure quantitatively the accuracy of the error bound derived in the previous paragraph with some numerical tests. Let us fix the following parameter values:

σ 0 = 1/60 e•(MW) -1 • s -1/2 , σ d = 1000/60 MW•s -1/2 , β = 0.002 e•(MW) -2 , η = 200 e•(MW) -2 , µ = 0 MW•s -1 , ν = 10 -10 e•(MW) -2 , γ = 10 -10 e•s•(MW) -2 and ρ = 0.8.
We start from the initial time t = 0, with a zero inventory X 0 = 0, and vary respectively the maturity T , the initial demand D 0 and the initial price Y 0 . We compute the probability for the final inventory to exceed the final demand P[ XT > D T ], the approximate value function ṽ(0, X 0 , Y 0 , D 0 ), and the error bound Ē(T, D 0 -X 0 , Y 0 ). The results are reported in Table 1 when varying T , in Table 2 when varying D 0 and in Table 3 when varying Y 0 . 1 shows that for time to maturity less than T = 24h, the probability for the final inventory to exceed the final demand is very small, and consequently the error bound is rather negligible. When the time horizon increases, the agent has the possibility to spread over time her trading strategies for reducing the price impact, and purchase more energy, in which case the probability for the final inventory to exceed the demand increases.

T (h) P[ XT > D T ] ṽ(0, X 0 , Y 0 , D 0 ) (e) Ē(T, D 0 -X 0 , Y 0 ) (e) 1 < 10
Table 2 shows that the probability for the final inventory to exceed the final demand, and the error bound are not much sensitive to the variations of the initial positive demand D 0 . Actually, the main impact is caused by the initial stock price, as observed in Table 3.

For small initial electricity price Y 0 , the agent will buy more energy in the intraday market and produce less. Therefore, the inventory will overtake with higher probability the demand, in which case the approximate value function can be significantly different from the original one, as observed from the error bound in Table 3 for Y 0 = 20.

D 0 (MW) P[ XT > D T ] ṽ(0, X 0 , Y 0 , D 0 ) (e) Ē(T, D 0 -X 0 , Y 0 ) (e) 500 <
Y 0 (e•(MW) -1 ) P[ XT > D T ] ṽ(0, X 0 , Y 0 , D 0 )(e) Ē(T, D 0 -X 0 , Y 0 ) (e) 500

Simulations

We plot trajectories of some relevant quantities that we simulate with the following set of parameters:

σ 0 = 1/60 e•(MW) -1 •s -1/2 , σ d = 1000/60 MW•s -1/2 , β = 0.002 e•(MW) -2 , η = 100 e•(MW) -2 , µ = 0 MW•s -1 , ρ = 0.8, ν = 4.00 × 10 -5 e•(MW) -2 , γ = 2.22e•s•(MW) -2 , T = 24h, X 0 = 0, D 0 = 50000 MW and Y 0 = 50 e•(MW) -1 .
For such parameter values, the probability P[ XT > D T ] is bounded above by 10 -16 , the error Ē(0, D 0 -X 0 , Y 0 ) is bounded by 2.82 × 10 -10 e, and ṽ(0, X 0 , Y 0 , D 0 ) = 1916700e.

The executed strategy (q, ξ * T ) can then be considered as very close to the optimal strategy. Figure 1 represents the evolution of the trading rate control (q t ) t∈[0,T ] derived in (3.6) for a given trajectory of price and demand, and this is consistent with the martingale property as shown in Proposition 3.1. Figure 2 represents a simulation of the quoted price Ŷt with impact and of the unaffected price Pt . Due to the buying strategy, i.e. positive q, we observe that the quoted price Ŷ is larger than P . In Figure 3, we plot the evolution of the optimal inventory ( Xt ) t∈[0,T ) , and of the forecast residual demand (D t ) t∈[0,T ] . We see that Xt is increasing, with a growth rate which looks constant as pointed out in Remark 3.2. At final time, if XT < D T (which is the case in our simulation), the agent uses her production leverage ξT , and achieves a final inventory: XT + ξT , which is represented by the peak at time T . From the expression (3.7) of ξT , the final imbalance cost is equal to

D T -XT -ξT = β η + β (D T -XT ),
and is then positive, as shown in Figure 3.

Jumps in the residual demand forecast

In this section, we incorporate the case where the residual demand forecast is subject to sudden changes induced by prediction errors on renewable production, which may be quite large. Our aim is to study the impact on the strategies obtained in the previous section, and we shall also neglect the delay in thermal plants production. The sudden changes in the demand forecast are modeled via a compound Poisson process N t = (N + t , N - t ) t≥0 with intensity λ > 0, where N + t is the counting process associated to positive jumps of the demand forecast with size δ + > 0, occurring with probability p + ∈ [0, 1], while N - t is the counting process associated to negative jumps of the demand forecast with size δ -< 0, occurring with probability p -= 1 -p + . We denote by δ := δ + p + + δ -p - the mean of the jump size of the demand forecast. The dynamics of the residual demand forecast D is then given by:

dD t = µdt + σ d dB t + δ + dN + t + δ -dN - t , (4.1)
where we add a jump component with respect to the model in (2.8). Moreover, as soon as a jump in the residual demand forecast occurs, this is impacted into the intraday electricity price since the main producers are assumed to have access to the whole updated forecast. We thus model the unaffected electricity price by:

Pt = P0 + σ 0 W t + π + N + t + π -N - t , (4.2)
where we add with respect to the Bachelier model in (2.6) a jump component of size π + > 0 (resp. π -< 0) when the jump on residual demand is positive (resp. negative), which means that a higher (resp. lower) demand induces an increase (resp. drop) of price. We denote by π := π + p + + π -p -the mean of the jump size of the intraday price. Given a trading rate q ∈ A, the dynamics of the quoted price Y is then governed by

dY t = νq t dt + σ 0 dW t + π + dN + t + π -dN - t . (4.3)
By considering this simplified modeling of demand forecast subject to sudden shift in terms of a Poisson process, we do not have additional state variables with respect to the no jump case of the previous section. Let us then denote by v = v (λ) (t, x, y, d) the value function to the optimal execution problem (2.4) with cost functional J = J (λ) (t, x, y, d, q, ξ), where we stress the dependence in λ for taking into account jumps in demand forecast. The value function in the no jump case derived in the previous section is denoted by v = v (0) .

As in the case with no jumps, there is no explicit solution to v (λ) due to the non negativity constraint on the final production: we shall first study the auxiliary execution problem without sign constraint on the final production, then provide an approximate solution to the original one with an estimation of the induced error approximation, and with some numerical illustrations. We compare the results with the no jump case by focusing on the impact of the jump components.

Auxiliary optimal execution problem

Similarly as in Subsection 3.1, we consider the optimal execution problem without non negativity constraint on the final production, denoted by ṽ = ṽ(λ) (t, x, y, d).

As in Theorem 3.1 for the case of the value function ṽ(0) without jumps, we have an explicit solution to this auxiliary problem. Theorem 4.1 The value function to the auxiliary optimization problem is explicitly given by: ṽ(λ) (t, x, y, d) = ṽ(0) (t, x, y, d)

+ λ 2 r(η, β)(T -t) π(T -t) + 2δ(ν(T -t) + 2γ) r(η, β) + ν (T -t) + 2γ (d -x) - λ 2 (T -t) 2 π -2r(η, β)δ r(η, β) + ν (T -t) + 2γ y +λγ p + (π + -r(η, β)δ + ) 2 + p -(π --r(η, β)δ -) 2 r(η, β) + ν 2 ln 1 + (r(η, β) + ν)(T -t) 2γ - λ 2 p + ((π + ) 2 -r(η, β)δ + (2π + + νδ + )) + p -((π -) 2 -r(η, β)δ -(2π -+ νδ -)) r(η, β) + ν (T -t) + λr(η, β) 2 2νµδ + λ((p + ) 2 δ + (π + + νδ + ) + (p -) 2 δ -(π -+ νδ -)) r(η, β) + ν (T -t) 2 +λ 2 γr(η, β) r(η, β)δ 2 + 2νp + p -δ + δ --((p + ) 2 δ + π + + (p -) 2 δ -π -) (r(η, β) + ν) (r(η, β) + ν)(T -t) + 2γ (T -t) 2 + 2λγr 2 (η, β)µδ (r(η, β) + ν) (r(η, β) + ν)(T -t) + 2γ (T -t) 2 - λ 2 π 2 48γ (T -t) 3 + λ 2 p + p -r(η, β) 2 2νδ + δ -+ δ -π + + δ + π - (r(η, β) + ν)(T -t) + 2γ (T -t) 3 + 1 8 4r(η, β)µλπ -λ 2 π 2 (r(η, β) + ν)(T -t) + 2γ (T -t) 3 , for (t, x, y, d) ∈ [0, T ] × R × R × R
, with an optimal trading rate given in feedback form by:

q(λ) s = q(λ) (T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s ), t ≤ s ≤ T q(λ) (t, d, y) := q(0) (t, d, y) + λ r(η, β)δt + π 4γ (r(η, β) + ν)t 2 (r(η, β) + ν)t + 2γ = q(0) (t, d + λδt, y + λ 2 πt) + λπ 4γ t, (4.4)
where q(0) is the optimal trading rate given in (3.6) in the case with no jump in the demand forecast. Here ( Xt,x,y,d , Ŷ t,x,y,d , D t,d ) denotes the solution to (2.1)-(4.3)-(4.1) when using the feedback control q(λ) , and starting from (x, y, d) at time t. Finally, the optimal production quantity is given by: ξ(λ)

T = η η + β D t,d T -Xt,x,y,d T . (4.5)
Proof. See Appendix.

Interpretation. The expression of the optimal trading rate q(λ) s , s ∈ [t, T ], as q(λ)

s = q(0) s + λ r(η, β)δ(T -s) + π 4γ (r(η, β) + ν)(T -s) 2 (r(η, β) + ν)(T -s) + 2γ
, where q(0) s = q(0) (T -s, D ) represents the optimal trading rate that the agent would use if she believes that the demand forecast will not jump, shows that under the information knowledge about jumps, the agent will purchase more (resp. less) electricity shares and this impact is all the larger, the larger the intensity λ of jumps, and the positive (resp. negative) mean δ and π of jump size in demand forecast and price are. On the other hand, the expression of q(λ) s as the sum of two terms:

q(λ) s = q(0) T -s, D t,d s + λδ(T -s), Ŷ t,x,y,d s + λ 2 π(T -s) + λπ 4γ (T -s), (4.6)
can be interpreted as follows. The first term is analog to the optimal trading rate in the no jump case, with an adjustment λδ(T -s) in the demand, which represents the expectation of the demand jump size up to the final horizon, and an adjustment λ 2 π(T -s) on the price, which represents half of the expectation of the price jump size up to the final horizon. The second term, λπ 4γ (T -s), is deterministic, and linear in time, and we shall see on the simulations for some parameter values that it can be dominant with respect to the first stochastic term. Moreover, notice that the equilibrium relation (3.9) in the no jump case between forecast intraday price and forecast marginal cost of production does not hold anymore in the presence of jumps, except at terminal date T :

Ŷ t,x,y,d T + 2γ q(λ) T = c ′ ( ξ(λ) T ). (4.7)
The unaffected price P in (4.2) is no more a martingale in presence of jumps, except when π = 0. It is actually a supermartingale when π < 0 (predominant negative jumps), and submartingale when π > 0 (predominant positive jumps). The next result shows that the optimal trading rate inherits the converse submartingale or supermartingale property of the price process.

Proposition 4.1 The optimal trading rate process (q (λ) s ) t≤s≤T in (4.4) is a supermartingale if π > 0, and a submartingale if π < 0. More precisely, the process {q

(λ) s + λπ 2γ (s -t), t ≤ s ≤ T } is a martingale. Proof. Notice that N ±
t is a Poisson process with intensity λp ± , and let us introduce the compensated martingale Poisson process Ñ ± t = N t -λp ± t. By applying Itô's formula to the trading rate process q(λ)

s = q(λ) (T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s
), t ≤ s ≤ T , and from the dynamics (2.1), (4.1) and (4.3), we have:

dq (λ) s = - ∂ q(λ) ∂t + (µ -q(λ) ) ∂ q(λ) ∂d + ν q(λ) ∂ q(λ) ∂y + λp + q(λ) (., . + δ + , . + π + ) -q(λ) + λp -q(λ) (., . + δ -, . + π -) -q(λ) (T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s )ds + ∂ q(λ) ∂d (T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s )σ d dB s + ∂ q(λ) ∂y (T -s, D t,d s -Xt,x,y,d s , Ŷ t,x,y,d s )σ 0 dW s + q(λ) (T -s, D t,d s -+ δ + -Xt,x,y,d s , Ŷ t,x,y,d s - + π + ) -q(λ) (T -s, D t,d s --Xt,x,y,d s , Ŷ t,x,y,d s - ) d Ñ + s + q(λ) (T -s, D t,d s -+ δ --Xt,x,y,d s , Ŷ t,x,y,d s - + π -) -q(λ) (T -s, D t,d s --Xt,x,y,d s , Ŷ t,x,y,d s - ) d Ñ - s .
Now, from the expression (4.4) of q(λ) (t, d, y), we see that:

- ∂ q(λ) ∂t + (µ -q(λ) ) ∂ q(λ) ∂d + ν q(λ) ∂ q(λ) ∂y +λ p + q(λ) (., . + δ + , . + π + ) + p -q(λ) (., . + δ -, . + π -) -q(λ) = - λπ 2γ ,
and then:

dq (λ) s = - λπ 2γ ds + r(η, β)σ d (r(η, β) + ν)(T -s) + 2γ dB s - σ 0 (r(η, β) + ν)(T -s) + 2γ dW s + r(η, β)δ + -π + (r(η, β) + ν)(T -s) + 2γ d Ñ + s + r(η, β)δ --π - (r(η, β) + ν)(T -s) + 2γ d Ñ - s . (4.8)
This proves the required assertions of the proposition.

Remark 4.1 The above supermartingale (or submartingale) property implies in particular that the mean of the optimal trading rate process (q

(λ)
s ) 0≤s≤T is decreasing (or increasing) in time, and so that the trajectory of the optimal inventory mean E[ X0,x,y,d s ], 0 ≤ s ≤ T , is concave (or convex). Moreover, from the martingale property of q(λ) s + λπ 2γ s, 0 ≤ s ≤ T , we have:

E[q (λ) s ] = q(λ) (T, d -x, y) -λπ 2γ s for 0 ≤ s ≤ T . Fix d, x
, y, and let us then denote by s(λ) := 2γ λπ q(λ) (T, d -x, y), which is explicitly written as:

s(λ) = T 2 + 1 λπ r(η, β)µ + λ(r(η, β)δ -π 2 ) T + r(η, β)(d -x) -y 1 + (r(η,β)+ν)T 2γ
We have the following cases:

• s(λ) ≤ 0 and π > 0: this may arise for large y, or d << x, or r(η, β)δ << π/2. In this

extreme case, dE[ X0,x,y,d s ] ds = E[q (λ)
s ] ≤ 0 for 0 ≤ s ≤ T , i.e. the trajectory of E[ X0,x,y,d s ], 0 ≤ s ≤ T , is decreasing, which means that the agent will "always" sell electricity shares since she takes advantage of high price, in order to decrease her inventory for approaching the demand, and because in average, the jump size of the demand is much lower than the positive jump size of the price.

• s(λ) ≤ 0 and π < 0: this may arise for small y, or d >> x, or r(η, β)δ >> π/2. In this

extreme case, dE[ X0,x,y,d s ] ds = E[q (λ)
s ] ≥ 0 for 0 ≤ s ≤ T , i.e. the trajectory of E[ X0,x,y,d s ], 0 ≤ s ≤ T , is increasing, which means that the agent will "always" buy electricity shares since she takes advantage of low price, in order to increase her inventory for approaching the demand, and because in average, the jump size of the price is much lower than the jump size of the demand.

• s(λ) ≥ T and π > 0: this may arise for r(η, β)δ >> π/2, d >> x or small y. In this other extreme case, the trajectory of E[ X0,x,y,d s ], 0 ≤ s ≤ T , is increasing, which means that the agent will "always" buy electricity shares at low price in order to approach the residual demand at final time.

• s(λ) ≥ T and π < 0: this may arise for r(η, β)δ << π/2, d << x or large y. The trajectory of E[ X0,x,y,d s ], 0 ≤ s ≤ T , is decreasing, which means that the agent will "always" sell electricity shares at high price in order to approach the residual demand at final time.

• 0 < s(λ) < T : in this regular case, it is interesting to comment on the two subcases:

if π > 0, the trajectory of s → E[ X0,x,y,d s ] is increasing for s ≤ s(λ) and then decreasing for s(λ) < s ≤ T . This means that the agent starts by purchasing electricity shares for taking profit of the positive price jumps (which have more impact than the negative price jumps as p + π + + p -π -> 0), and then resells shares in order to achieve the equilibrium relation (4.7).

if π < 0, i.e. the negative jumps have more impact than the positive ones: the agent starts by selling electricity shares and then purchases shares.

Approximate solution

We turn back to the original optimal execution problem with the non negativity constraint on the final production, and as in Section 3.2, we use the approximate strategy consisting in the trading rate q(λ) derived in (4.4), and of the truncated nonnegative final production:

ξ(λ), * T := ξ(λ) T 1 ξ(λ) T ≥0 = ξ+ (D t,d T -Xt,x,y,d T ),
with ξ(λ) T given in (4.5). We measure the relevance of this strategy (q (λ) , ξ(λ), * T ) ∈ A×L 0 + (F T ) by estimating the induced error:

E (λ) 1 (t, x, y, d) := J (λ) (t, x, y, d; q(λ) , ξ(λ), * T ) -v (λ) (t, x, y, d), for (t, x, y, d) ∈ [0, T ] × R × R × R,
and also measure the approximation error on the value functions:

E (λ) 2 (t, x, y, d) := v (λ) (t, x, y, d) -ṽ(λ) (t, x, y, d). Proposition 4.2 For all (t, x, y, d) ∈ [0, T ] × R × R × R, we have 0 ≤ E (λ) i (t, x, y, d) ≤ ηr(η, β) 2β V (T -t)E ψ m (λ) (T -t, d -x, y) + Σ -,t T V (T -t) (4.9)
for i = 1, 2, where ψ, m, V are defined in Proposition 3.2,

m (λ) (t, d, y) = m t, d, y + λ π 2 -r(η, β)δ t (4.10) + λ r(η, β)δ -π r(η, β) + ν t - 2γ r(η, β) + ν ln 1 + r(η, β) + ν 2γ t , and
Σ -,t T = T t δ -(ν(T -s) + 2γ) + π -(T -s) (r(η, β) + ν)(T -s) + 2γ dN - s ≤ 0, a.s.
Proof. By the same arguments as in Proposition 3.2, we have

0 ≤ E (λ) i (t, x, y, d) ≤ E (λ) (t, x, y, d) := J (λ) (t, x, y, d; q(λ) , ξ(λ), * T ) -J (λ) (t, x, y, d; q(λ) , ξ(λ)
T ), for i = 1, 2, and

E (λ) (t, x, y, d) = ηr(η, β) 2β E D t,d T -Xt,x,y,d T 2 1 D t,d T -Xt,x,y,d T <0 , (4.11)
for (t, x, y, d) ∈ [0, T ] × R × R × R. Now, recall from (4.8) that:

dq (λ) s = -λ π 2γ + r(η, β)δ -π (r(η, β) + ν)(T -s) + 2γ ds + r(η, β)σ d (r(η, β) + ν)(T -s) + 2γ dB s - σ 0 (r(η, β) + ν)(T -s) + 2γ dW s + r(η, β)δ + -π + (r(η, β) + ν)(T -s) + 2γ dN + s + r(η, β)δ --π - (r(η, β) + ν)(T -s) + 2γ dN - s ,
where we write the dynamics directly in terms of the Poisson processes N ± . By integration, we deduce the (path-dependent) expression of q(λ) s , t ≤ s ≤ T :

q(λ) s = q(λ) t - λπ 2γ (s -t) + λ r(η, β)δ -π r(η, β) + ν ln (r(η, β) + ν)(T -s) + 2γ (r(η, β) + ν)(T -t) + 2γ + s t r(η, β)σ d (r(η, β) + ν)(T -u) + 2γ dB u - s t σ 0 (r(η, β) + ν)(T -u) + 2γ dW u + s t r(η, β)δ + -π + (r(η, β) + ν)(T -u) + 2γ dN + u + s t r(η, β)δ --π - (r(η, β) + ν)(T -u) + 2γ dN - u , with q(λ) t = q(λ) (T -t, d -x, y).
We thus obtain the expression of the final spread between demand and inventory:

D t,d T -Xt,x,y,d T = d -x + µ(T -t) + T t σ d dB s + T t δ + dN + s + T t δ -dN - s - T t q(λ) s ds = m (λ) (T -t, d -x, y) + T t σ d (ν(T -s) + 2γ) (r(η, β) + ν)(T -s) + 2γ dB s + T t σ 0 (T -s) (r(η, β) + ν)(T -s) + 2γ dW s + T t δ + (ν(T -s) + 2γ) + π + (T -s) (r(η, β) + ν)(T -s) + 2γ dN + s + T t δ -(ν(T -s) + 2γ) + π -(T -s) (r(η, β) + ν)(T -s) + 2γ dN - s , (4.12)
by Fubini's theorem, and where

m (λ) (t, d, y) := d + µt -tq (λ) (t, d, y) + λπ 2γ t 0 sds - λ r(η, β)δ -π r(η, β) + ν t 0 ln (r(η, β) + ν)s + 2γ (r(η, β) + ν)t + 2γ ds,
is explicitly written as in (4.10) after some straightforward calculation. Denoting by

∆ t,x,y,d T the continuous part of D t,d T -Xt,x, y,d
T consisting in the three first terms in the rhs of (4.12), and by Σ +,t T , Σ -,t T the jump parts consisting in the two last terms of (4.12), so that

D t,d T -Xt,x,y,d T = ∆ t,x,y,d T + Σ +,t T + Σ -,t T ,
we notice that ∆ t,x,y,d T follows a normal distribution law with mean m (λ) (T -t, d -x, y) and variance V (T -t), independent of Σ ±,t T . Then, conditionally on Σ ±,t T , D t,d T -Xt,x,y,d

T follows a normal distribution law with mean m (λ) (T -t, d -x, y) + Σ +,t T + Σ -,t
T , and variance V (T -t), and this implies from (4.11) that:

E (λ) (t, x, y, d) = ηr(η, β) 2β V (T -t)E ψ m (λ) (T -t, d -x, y) + Σ +,t T + Σ -,t T V (T -t) ≤ ηr(η, β) 2β V (T -t)E ψ m (λ) (T -t, d -x, y) + Σ -,t T V (T -t) ,
since Σ +,t T ≥ 0 a.s. and ψ is non-increasing.

Comments on the approximation error. Let us discuss about the accuracy of the upper bound in (4.9):

Ē(λ) (T -t, d -x, y) := ηr(η, β) 2β V (T -t)E ψ m (λ) (T -t, d -x, y) + Σ -,t T V (T -t) , First, notice that m (λ) (T -t, d -x, y) + Σ -,t T ∼ m(T -t, d -x, y) a.
s. in the limiting regimes where T -t goes to zero, d -x or y goes to infinity. Therefore, by dominated convergence theorem, Ē(λ) (T -t, d -x, y) converges to zero in these limiting regimes as in the no jump case. However, we are not able to derive an asymptotic limit as in the no jump case of Proposition 3.3, except when Σ -,t T = 0, i.e. δ -= π -= 0, fo which we get the same asymptotic limit. Actually, in the presence of negative jumps on the demand, it is intuitively clear that our approximation should be less accurate than in the no jump case since the probability for the residual demand to stay above the final inventory is decreasing. Anyway, the explicit strategies (q (λ) , ξ(λ), *

T

) still provide a very accurate approximation of the optimal strategies at least in these limiting regimes, as illustrated in the next paragraph.

Numerical results

We plot trajectories of some relevant quantities that we simulate with the same set of parameters as in Paragraph 3.3.2:

σ 0 = 1/60 e•(MW) -1 • s -1/2 , σ d = 1000/60 MW•s -1/2 , β = 0.002 e•(MW) -2 , η = 200 e•(MW) -2 , µ = 0 MW•s -1 , ρ = 0.8, ν = 4.00 • 10 -5 e•(MW) -2 , γ = 2.22e•s•(MW) -2 , T = 24h, X 0 = 0, D 0 = 50000 MW and Y 0 = 50 e•(MW) -1 .
Moreover, we fix the probability of positive jumps, p + = 1 (then all jumps are positive: p -= 0), and the following values for the jump components: λ = 1.5/(3600 • 24) s -1 , π + = 10 e•(MW) -1 , δ + = 1500 MW.

For such parameter values, we observe two occurrences of jumps on the trajectories of the demand of price. Moreover, the probability P[ XT > D T ] is bounded above by 2.92 × 10 -16 , the error Ē(λ) (0, D 0 -X 0 , Y 0 ) is bounded by 2.66 × 10 -5 e, and ṽ(λ) (0, X 0 , Y 0 , D 0 ) = 2020950e.

The executed strategy (q (λ) ,

) can then be considered as very close to the optimal strategy. This has to be compared with the numerical result obtained in the previous section in the no jump case where we obtained a lower expected total cost: ṽ(0, X 0 , Y 0 , D 0 ) = 1916700e.

Figure 4 represents the evolution of the trading rate (q

(λ) t ) t∈[0,T ] ,
and we see that it is decreasing consistently with the supermartingale property in Proposition 4.1. Actually, we observe that the deterministic part in (4.6), which is linear in time, dominates the stochastic part. The interpretation of the strategy is the following: since positive price jumps are expected, the agent purchases a large number of shares in electricity with the hope to sell it later at a higher price thanks to the possible occurrence of a positive jump. At the price jump times, which can be visualized in Figure 5, we notice that the control q(λ) reacts by a decrease in the trading rate. The reaction to the second jump is more sensible than to the first jump since it occurs a short time before the final horizon T , where the objective is also to achieve the equilibrium relation (4.7) between price and marginal cost. Finally, we observe clearly in Figure 6 the concavity of the trajectory of the optimal inventory process ( Xt ) t∈[0,T ) , as expected from Remark 4.1. This emphasizes the double objective of the agent: on one hand, the purchase of electricity shares for taking profit of the positive price jumps, and on the other hand the resale of electricity shares for attaining the equilibrium relation between price and marginal cost at terminal date. We also plot the production ξT at the final time T on Figure 6, and observe as in the no jump case that the imbalance cost D T -XT -ξT is positive. Next, we plot trajectories with the same set of parameters, but with p + = 0.3 (i.e. p -= 0.7), π -= -10 e•(MW) -1 , δ -= -1500 MW. There are, in average, more negative than positive jumps. Now ṽ(λ) (0, X 0 , Y 0 , D 0 ) = 1756330e.

Figure 7 shows that the trading rate (q (λ) t ) t∈[0,T ] is increasing, which is consistent with the submartingale property in Proposition 4.1: the deterministic part in (4.6) dominates the stochastic part. Since negative jumps are more expected than negative jumps are, the agent first sells a large number of shares in electricity with the hope to buy it later at a lower price thanks to the possible occurrence of jumps, that should be mainly negative. Here, the control reacts to the negative price jumps by an increase in the trading rate. Finally, in Figure 9 we observe the convexity of the trajectory of the optimal inventory ( Xt ) t∈[0,T ) process, as expected from Remark 4.1. We also plot the production ξT at the final time T on that figure.

Delay in production

In this section, we consider the more realistic situation when there is delay in the production, assumed to be fixed equal to h ∈ [0, T ], and we denote by v = v h the value function to the associated optimal execution problem, as defined in (2.9), where we stress the dependence in the delay h. Our aim is to show how one can reduce the problem with delay to a suitable problem without delay, and then solve it explicitly. We shall consider the problem without jumps on demand forecast and price, but the same argument also works for the case with jumps.

Explicit solution with delay

For simplicity of presentation, and without loss of generality, we shall focus on the derivation of the value function v h (t, x, y, d) for an initial time t = 0, and fixed

(x, y, d) ∈ R × R × R.
Given a control trading rate q ∈ A, and from pathwise uniqueness for the solution to the dynamics (2.1), (2.7), (2.8), we observe that for any ξ ∈ L 0 (F T -h ):

   X 0,x T + ξ = X T -h,X 0,x T -h +ξ T a.s. Y 0,y T = Y T -h,Y 0,y T -h T , D 0,d T = D T -h,D 0,d T -h T a.s. (5.1)
To alleviate notations, we shall omit the dependence in the fixed initial conditions (x, y, d), and simply write X s = X 0,x s , Y s = Y 0,y s , D s = D 0,d s , for s ≥ 0, v h = v h (0, x, y, d), and J(0; q, ξ) = J(0, x, y, d; q, ξ) for the the cost functional in (2.10). By the tower property of conditional expectations and from (5.1), the cost functional can be written, for all q ∈ A, ξ ∈ L 0 (F T -h ), as:

J(0; q, ξ) = E T -h 0 q s Y s + γq s )ds + c(ξ) + J(T -h, X T -h + ξ, Y T -h , D T -h ; q, 0) (5.2) ≥ E T -h 0 q s Y s + γq s )ds + c(ξ) + v N P (T -h, X T -h + ξ, Y T -h , D T -h ) ,
by definition (3.8) of the value function v N P for the optimal execution problem without production, i.e. the pure retailer problem. Since q is arbitrary in A, this shows that:

inf q∈A J(0; q, ξ) (5.3) ≥ inf q∈A E T -h 0 q s Y s + γq s )ds + c(ξ) + v N P (T -h, X T -h + ξ, Y T -h , D T -h ) ,
for all ξ ∈ L 0 (F T -h ). Now, given q ∈ A, and ξ ∈ L 0 (F T -h ), let us consider the trading rate qNP,ξ in A T -h solution to the pure retailer problem:

v N P (T -h, X T -h + ξ, Y T -h , D T -h ),
hence starting at time T -h from an inventory X T -h + ξ. By considering the process q ∈ A defined by: qs = q s for 0 ≤ s < T -h, and qs = qNP,ξ s , for T -h ≤ s ≤ T , we then obtain from (5.2): (5.4) which proves together with (5.3) the equality: inf q∈A J(0; q, ξ)

J(0; q, ξ) = E T -h 0 q s Y s + γq s )ds + c(ξ) + v N P (T -h, X T -h + ξ, Y T -h , D T -h ) ,
(5.5)

= inf q∈A E T -h 0 q s Y s + γq s )ds + c(ξ) + v N P (T -h, X T -h + ξ, Y T -h , D T -h ) ,
for all ξ ∈ L 0 (F T -h ). Therefore, v h = inf q∈A,ξ∈L 0 + (F T -h ) J(0; q, ξ) can be written as:

v h = inf q∈A,ξ∈L 0 + (F T -h ) E T -h 0 q s Y s + γq s )ds + c(ξ) + v N P (T -h, X T -h + ξ, Y T -h , D T -h ) . (5.6)
In other words, the original problem with delay in production is formulated as an optimal execution problem without delay, namely with final horizon T -h, and terminal cost function:

C h (x, y, d, ξ) := c(ξ) + v N P (T -h, x + ξ, y, d).
Notice from the explicit expression of v N P in Remark 3.1 that this cost function C h does not depend on T , and is in the form:

C h (x, y, d, ξ) = C h (0, y, d -x -ξ, 0) = c(ξ) + v N P (T -h, 0, y, d -x -ξ).
The optimization over q and ξ in (5.6) is done separately: the production ξ ∈ L 0 + (F T -h ) is decided at time T -h, after the choice of the trading rate (q s ) for 0 ≤ s ≤ T -h (leading to an inventory X T -h ), and is determined optimally from the optimization a.

= η η + β (νh + 2γ)(µh + d) + hy (r(η, β) + ν)h + 2γ . ( 5.7)
The problem (5.6) is then rewritten as

v h = inf q∈A E T -h 0 q s Y s + γq s )ds + C + h (D T -h -X T -h , Y T -h ) , (5.8)
where

C + h (d, y) := C h (0, y, d -ξh,+ (d), 0).
Notice that when h = 0, we retrieve the expressions in the no delay case: ξ0,+ = ξ+ in (3.1), C + 0 = C + in (3.3) and v 0 = v in (3.2). As in the no delay case, there is no explicit solution to the HJB equation associated to the stochastic control problem (5.8). We then consider the approximate control problem where we relax the non negativity constraint on the production, i.e. ṽh = inf q∈A,ξ∈L 0 (F T -h ) J(0; q, ξ). Therefore by following the same arguments as above, the corresponding value function is written as:

ṽh = inf q∈A E T -h 0 q s Y s + γq s )ds + Ch (D T -h -X T -h , Y T -h ) , (5.9)
where Ch (d, y) := C h (0, y, d -ξh (d), 0).

From the explicit expressions of ξh in (5.7) and v N P in Remark 3.1, it appears after some tedious but straightforward calculations that the auxiliary terminal cost function Ch simplifies remarkably into:

Ch (d, y) = ṽ0 (T -h, 0, y, d) + K h ,
where ṽ0 is the auxiliary value function without delay explicitly obtained in Theorem 3.1, and K h is a constant depending only on the delay h and the parameters of the model, given explicitly by

K h = η 2 2 σ 2 0 + σ 2 d ν 2 + 2ρσ 0 σ d ν (η + β)(η + ν)(r(η, β) + ν) h + γ σ 2 0 + σ 2 d η 2 -2ρσ 0 σ d η (η + ν) 2 ln 1 + (η + ν)h 2γ -γ σ 2 0 + σ 2 d r 2 (η, β) -2ρσ 0 σ d r(η, β) (r(η, β) + ν) 2 ln 1 + (r(η, β) + ν)h 2γ .
One easily checks that K h = 0 for h = 0, and K h is increasing with h (actually the derivative of K h w.r.t. h is positive), hence in particular K h is nonnegative. Plugging into (5.9), we then get

ṽh = inf q∈A E T -h 0 q s Y s + γq s )ds + ṽ0 (T -h, X T -h , Y T -h , D T -h ) + K h . (5.10)
Therefore, by using the dynamic programming principle for the control problem ṽ0 = ṽ0 (0, x, y, d) in (3.4), we obtain this remarkable relation

ṽh = ṽ0 + K h , (5.11)
which explicitly relates the (approximate) value function with and without delay. As expected from the very definition of ṽh , this relation implies that ṽh -ṽ0 is nonnegative, and is increasing in h. This is consistent with the intuition that when making the production choice in advance, we do not take into account the future movements of the price and of the residual demand, which should therefore lead to an average positive correction of the cost. More precisely, the relation (5.11) gives an explicit quantification of the delay impact via the term K h (which does not depend on the state variables x, y, d) in function of the various model parameters. Moreover, the optimal control of the stochastic control problem (5.10) over [0, T -h) is explicitly given by the optimal control (q s ) 0≤s≤T -h of problem ṽ0 without delay in Theorem 3.1.

Let us now consider the following strategy (q h,+ , ξh, * T -h ) ∈ A × L 0 + (F T -h ) for the original problem v h with delay:

• Before T -h, follow the trading strategy qh,+ s = qs , s < T -h, corresponding to the solution of the auxiliary problem without delay as if production choice is made at time T , and leading to an inventory XT -h , and an impacted price ŶT -h .

• At time T -h, choose the production quantity: In order to estimate the quality of this approximate strategy with respect to the optimal trading problem v h , measured by

ξh, * T -h := ξh,+ (D T -h -XT -h , ŶT -h ).
E h 1 := J(0; qh,+ , ξh, * T -h ) -v h ,
we shall compare it with the following strategy (q h , ξh

T -h ) ∈ A × L 0 (F T -h ):
• Before T -h, follow the trading strategy qh s = qs , s < T -h, corresponding to the solution of the auxiliary problem without delay as if production choice is made at time T , and leading to an inventory XT -h , and an impacted price ŶT -h .

• At time T -h, choose the "production" quantity (which can be negative): Then, by construction and following the arguments (see in particular (5.4), (5.9), (5.10)) leading to the expression (5.11) of ṽh , we see that (q h , ξh T -h ) is the optimal solution for ṽh , i.e. ṽh = J(0; qh , ξh T -h ). On the other hand, since ṽh ≤ v h ≤ J(0; qh,+ , ξh, * T -h ), we deduce that

ξh T -h = ξh (D T -h -XT -h , ŶT -h ).
max(v h -ṽh , E h 1 ) ≤ Ēh := J(0; qh,+ , ξh, * T -h ) -J(0; qh , ξh T -h ).
Now, from the expression (5.4) of J, and by same arguments as in the proof of Proposition 3.2 (see the derivation of relation (3.15)), we have

Ēh = E C + h (D T -h -XT -h , ŶT -h ) -Ch (D T -h -XT -h , ŶT -h ) = E v N P (T -h, 0, ŶT -h , D T -h -XT -h -ξh, * T -h ) + c( ξh, * T -h ) -v N P (T -h, 0, ŶT -h , D T -h -XT -h -ξh T -h ) -c( ξh T -h ) = ηr(η, β) 2β (r(η, β) + ν)h + 2γ (η + ν)h + 2γ V h (T )ψ m(T, d -x, y) V h (T )
where m and ψ are defined as in (3.12), and

V h (T ) = T h σ 2 0 s 2 + σ 2 d (νs + 2γ) 2 + 2ρσ 0 σ d s(νs + 2γ) (r(η, β) + ν)s + 2γ 2 ds.
We recover when h = 0 the expression in Proposition 3.2 of the error in the no delay case, and notice that Ēh decreases when the delay increases: indeed, the error comes from the trading procedure before deciding how much to produce, which is dictated by the auxiliary problem, in which the final "production" can be negative. After T -h, the followed control is optimal, as there remains no production decision at some further date. The shorter the period before making the production decision is, the weaker the error is.

Let us finally discuss some properties of the (approximate) optimal trading strategy qh,+ . Recalling from Proposition 3.1 that the optimal trading rate is a martingale in the no delay case, we see by construction of (q h,+ s ) 0≤s≤T that it is a martingale on [0, T -h) and a martingale on [T -h, T ]. Moreover, for any s ∈ [T -h, T ], and t ∈ [0, T -h), we have

E qh,+ s |F t = E E qNP, ξh, * T -h s |F T -h |F t = E qNP, ξh, * T -h T -h |F t = E η(µh + D T -h -XT -h -ξh, * T -h ) -ŶT -h (η + ν)h + 2γ F t = E η(µh + D T -h -XT -h -ξh T -h ) -ŶT -h (η + ν)h + 2γ F t + η (η + ν)h + 2γ E ξh T -h -ξh, * T -h F t = E r(η, β)(µh + D T -h -XT -h ) -ŶT -h (r(η, β) + ν)h + 2γ F t + η (η + ν)h + 2γ E ξh T -h -ξh, * T -h F t = E qT -h |F t + η (η + ν)h + 2γ E ξh T -h -ξh, * T -h F t = qt + η (η + ν)h + 2γ E ξh T -h 1 ξh T -h <0 F t ≤ qt = qh,+ t .
(5.12)

where we used the tower rule for conditional expectations, the martingale property and the explicit expression of q N P, ξh, * T -h in Remark 3.1, the definition of ξh T -h , the martingale property and explicit expression of q in Theorem 3.1, and finally the fact that ξh, * T -h = ξh T -h 1 ξh T -h ≥0 . This shows in particular the supermartingale property of qh,+ over the whole period [0, T ]. Notice that the same arguments as for the derivation of (5.12) shows the martingale property over the whole period [0, T ] of the optimal trading strategy qh associated to the auxiliary problem ṽh . Moreover, by the martingale property of qh,+ on [0, T -h), and relation (5.12), we see that the (approximate) optimal inventory process Xh,+ with trading rate qh,+ has on average, a growth rate dE[ Xh,+ s ] ds , which is piecewise constant, equal to:

E[q h,+ s ] = q0 , for 0 ≤ s < T -h q(h) 0 := q0 + η (η+ν)h+2γ E ξh T -h 1 ξh T -h <0 < q0 , for T -h ≤ s ≤ T,
with q0 = r(η,β))(µT +d-x)-y (r(η,β)+ν)T +2γ , and q(h) 0 = q0 -ηr(η, β) β((η + ν)h + 2γ)

V h (T ) ψ m(T, d -x, y) V h (T ) , where ψ(z) := φ(z) -zΦ(-z), z ∈ R is a nonnegative function, as pointed out in (3.17).

Numerical results

We plot figures showing relevant trajectories with the same parameters as in Section 3.3.2. We add a delay h = 4 hours: the production choice has to be made four hours before the end of the trading period. We have ṽh (0, X 0 , Y 0 , D 0 ) = 1925460e, which is slightly higher than the value ṽ0 (0, X 0 , Y 0 , D 0 ) = 1916700e without delay. On Figure 10, we see that at time T -h, the positive production choice ξh, * T -h is made, and then we go on buying shares on the intraday market in order to go nearer to the demand forecast, with a smaller slope of trading rate. On Figure 11, which represents the control process without the last hour of trading (because oscillations then become overwhelming), we see that after date T -h, as we do not plan to use final production leverage any more, the approximate optimal control process qh,+ oscillates a lot as we are approaching the end of trading time. We can compare with Figure 1 to assert that qualitatively, the control in the problem with no production oscillates more than the one in the problem with final production, as in the former problem, the intraday market is the only way to seek to reach the equilibrium. with the initial conditions A(0) = 1 2 r(η, β), B(0) = 0, F (0) = 0, G(0) = 0, H(0) = 0, K(0) = 0. We first solve the Riccati system relative to the triple (A, B, F ), and obtain: Notice that q is linear, and Lipschitz in x, y, d, uniformly in time t, and so given an initial state (x, y, d) at time t, there exists a unique solution ( Xt,x,y,d , Ŷ t,x,y,d , D ds] < ∞, hence q ∈ A t . We now call on a classical verification theorem (see e.g. Theorem 3.5.2 in [START_REF] Pham | Continuous-time stochastic control and optimization with financial applications[END_REF]), which shows that w is indeed equal to the value function ṽ, and q is an optimal control. Finally, once the optimal trading rate q is determined, the optimal production is obtained from the optimization over ξ ∈ R of the A.2 Proof of Theorem 4.1

A(t) =
The Hamilton-Jacobi-Bellman (HJB) integro-differential equation arising from the dynamic programming associated to the stochastic control problem ṽ = ṽ(λ) with jumps in the dynamics of Y and D is:

        
∂ṽ (λ) ∂t + inf q∈R q ∂ṽ (λ) ∂x + νq ∂ṽ (λ) ∂y + µ ∂ṽ (λ) ∂d + 1 2 σ 2 0 ∂ 2 ṽ(λ)

∂y 2 + 1 2 σ 2 d
∂ 2 ṽ(λ) ∂d 2 + ρσ 0 σ d ∂ 2 ṽ(λ) ∂y∂d + q(y + γq)

+ λ p + ṽ(λ) (t, x, y + π + , d + δ + ) + p -ṽ(λ) (t, x, y + π -, d + δ -) -ṽ(λ) (t, x, y, d) = 0 ṽ(λ) (T, x, y, d) = C(d -x) = 1 2 r(η, β)(d -x) 2 .

where G and H are given from the no jump case (A.4). Finally, after some tedious but straightforward calculations, we explicitly obtain K λ from the last equation: Again, notice that q(λ) is linear, and Lipschitz in x, y, d, uniformly in time t, and so given an initial state (x, y, d) at time t, there exists a unique solution ( Xt,x,y,d , Ŷ t,x,y,d , D t,d ) t≤s≤T to (2.1)-(4.3)-(4.1) with the feedback control q(λ)

K λ (t

T

  is below the terminal demand D t,d T , by choosing a quantity proportional to this spread D t,d T -Xt,x,y,d T

Figure 2 :Figure 3 :

 23 Figure 2: Simulation of the quoted impacted price Ŷ and of the unaffected price P .

Figure 4 :Figure 5 :Figure 6 :

 456 Figure 4: Evolution of the trading rate control q(λ)

Figure 7 :Figure 8 :-Figure 9 :

 789 Figure 7: Evolution of the trading rate control q(λ)

  s. at T -h of the terminal cost C h (X T -h , Y T -h , D T -h , ξ). It is then given in feedback form by ξ * T -h = ξh,+ (D T -h -X T -h , Y T -h )where ξh,+ (d, y) := arg min ξ≥0 C h (0, y, d -ξ, 0) = arg min ξ≥0 c(ξ) + v N P (T -h, 0, y, d -ξ) , hence explicitly given from the expression of v N P in Remark 3.1 by: ξh,+ (d, y) = ξh (d, y)1 ξh (d,y)≥0 , ξh (d, y) :

•

  Between time T -h and T , follow the trading strategy qh,+ s = qNP, ξh, * T -h s , T -h ≤ s ≤ T , corresponding to the solution of the problem without production, and starting at T -h from an inventory XT -h + ξh, * T -h .

•

  Between time T -h and T , follow the trading strategy qh s = qNP, ξh T -h s , T -h ≤ s ≤ T , corresponding to the solution of the problem without production, and starting at T -h from an inventory XT -h + ξh T -h .

Figure 10 :

 10 Figure 10: Evolution of the inventory X (with production choice at time T -h) and of the forecast residual demand D

2 t

 2 r(η,β)( ν 2 t+γ) (r(η,β)+ν)t+2γ , B(t) = -1 (r(η,β)+ν)t+2γ , F (t) = r(η,β)t (r(η,β)+ν)t+2γ .(A.3)Then we solve the first-order linear system of ODE relative to the pair (G, H), which leads to the explicit solution:G(t) = 2µtA(t),andH(t) = -2r(η, β)µtB(t). (A.4)Finally, we explicitly obtain K from the last equation:K(t) = γ σ 2 0 + σ 2 d r 2 (η, β) -2ρσ 0 σ d r(η, β) r(η, β) + ν 2 ln 1 + (r(η, β) + ν)t 2γ + σ 2 d r(η, β)ν + 2ρσ 0 σ d r(η, β) -σ 2 0 2 r(η, β) + ν t + r(η, β)µ 2 t 2 ( ν 2 t + γ) (r(η, β) + ν)t + 2γ . (A.5)By construction, w in (A.2) with A, B, F , G, H and K explicitly given by (A.3)-(A.4)-(A.5), is a smooth solution with quadratic growth condition to the HJB equation (A.1). Moreover, the argmin in HJB equation for w is attained for q(t, x, y, d) = -, β)(µ(T -t) + d -x) -y (r(η, β) + ν)(T -t) + 2γ =: q(T -t, d -x, y).

  terminal cost C(D t,d T -Xt,x,y,d T , ξ), hence given by: ξT = η η+β (D t,d T -Xt,x,y,dT).

  , β)(µ(T -t) + d -x) -y (r(η, β) + ν)(T -t) + 2γ + λ r(η, β)δ(T -t) + π 4γ (r(η, β) + ν)(T -t) 2 (r(η, β) + ν)(T -t) + 2γ =: q(λ) (T -t, d -x, y).

  This clearly implies that both errors E 1 and E 2 are nonnegative, and

max(E 1 (t, x, y, d), E 2 (t, x, y, d)) ≤ E(t, x, y, d) := J(t, x, y, d; q, ξ * T ) -J(t, x, y, d; q, ξT ).

We now focus on the upper bound for E. By definition of J in (2.10), ξT and ξ * T in (3.7) and (3.11), we have

Table 1 :

 1 Y 0 = 50 e•(MW) -1 and D 0 = 50000 MW

		-16	1.88 × 10 6	< 10 -16
	8	< 10 -16	1.88 × 10 6	< 10 -16
	24	< 10 -16	1.89 × 10 6	4.16 × 10 -12
	50	7.72 × 10 -13	1.90 × 10 6	2.48 × 10 -4
	Table

Table 2 :

 2 T = 24 h and Y 0 = 50 e•(MW) -1

		10 -16	-5.86 × 10 5	4.16 × 10 -12
	5000	< 10 -16	-3.62 × 10 5	4.16 × 10 -12
	50000	< 10 -16	1.89 × 10 6	4.16 × 10 -12
	500000	< 10 -16	2.44 × 10 7	4.16 × 10 -12

Table 3 :

 3 T = 24 h and D 0 = 50000 MW

		< 10 -16	2.51 × 10 6	< 10 -16
	50	< 10 -16	1.89 × 10 6	4.16 × 10 -12
	40	9.51 × 10 -15	1.61 × 10 6	3.80 × 10 -4
	30	4.57 × 10 -10	1.29 × 10 6	1.30 × 10 -2
	20	2.23 × 10 -5	9.13 × 10 5	1.26 × 10 3

  ) = K(t) + λγ p + (π + -r(η, β)δ + ) 2 + p -(π --r(η, β)δ -) 2 + ((π + ) 2 -r(η, β)δ + (2π + + νδ + )) + p -((π -) 2 -r(η, β)δ -(2π -+ νδ -)) + ) 2 δ + (π + + νδ + ) + (p -) 2 δ -(π -+ νδ -)) r(η, β) + ν t 2 + λ 2 γr(η, β) r(η, β)δ 2 + 2νp + p -δ + δ --((p + ) 2 δ + π + + (p -) 2 δ -π -) (r(η, β) + ν) (r(η, β) + ν)t + 2γ t 2 + 2λγr(η, β) 2 µδ (r(η, β) + ν) (r(η, β) + ν)t + 2γ t 2 -λ 2 π 2 48γ t 3 + λ 2 p + p -r(η, β) 2 2νδ + δ -+ δ -π + + δ + π - (r(η, β) + ν)t + 2γ t 3 , β)µλπ -λ 2 π 2 (r(η, β) + ν)t + 2γ t 3 ,with K in (A.5). The function w(λ) in (A.7) may thus be rewritten as the sum of w in (A.2) and another function of t, d -x and y, and is by construction a smooth solution with quadratic growth condition to the HJB equation (A.6). Moreover, the argmin in HJB equation for w(λ) is attained for

			r(η, β) + ν	2	ln 1 +	(r(η, β) + ν)t 2γ
	-	λ 2	p r(η, β) + ν	t
	+ 2νµδ + λ((p + λr(η, β) 2 8 1 4r(η

This study was supported by FiME (Finance for Energy Market Research Centre) and the "Finance et Développement Durable -Approches Quantitatives" EDF -CACIB Chair. The authors would like to thank Marc Ringeisen, Head of EDF R&D Osiris Department for insightful discussion on trading and intraday market.

† EDF R&D and Finance for Energy Market Research

The Hamilton-Jacobi-Bellman (HJB) equation arising from the dynamic programming associated to the stochastic control problem (3.4) is:

The argmin in HJB is attained for q(t, x, y, d) = -1 2γ ∂ṽ ∂x + ν ∂ṽ ∂y + y , and the HJB equation is rewritten as:

We look for a candidate solution to HJB in the form

for some deterministic functions A, B, F , G, H and K. Plugging the candidate function w into equation (A.1), we see that w is solution to the HJB equation iff the following system of ordinary differential equations (ODEs) is satisfied by A, B, F , G, H and K:

Notice that with respect to the no jump case, there is in addition a linear integro-differential term in the HJB equation (which does not depend on the control), and the argmin is attained as in the no jump case for

The HJB equation is then rewritten as

(A.6) We look again for a candidate solution to (A.6) in the form

for some deterministic functions A λ , B λ , F λ , G λ , H λ and K λ . Plugging the candidate function w(λ) into equation (A.6), we see that w(λ) is solution to the HJB equation iff the following system of ordinary differential equations (ODEs) is satisfied by

We first solve the Riccati system relative to the triple (A λ , B λ , F λ ), which is the same as in the no jump case, and therefore obtain: A λ = A, B λ = B, F λ = F as in (A.3). Then we solve the first-order linear system of ODE relative to the pair (G λ , H λ ), which involves the jump parameters λ, π and δ, and get: