Optimal Bidding in Wind Farm Management

Abstract

We study the problem of wind farm management, in which the manager commits himself to deliver energy in some future time. He reduces the consequences of uncertainty by using a storage facility (a battery, for instance). We consider a simplified model in discrete time, in which the commitment is for the next period. We solve an optimal control problem to define the optimal bidding decision. Application to a real dataset is done, and the optimal size of the battery (or the overnight costs) for the wind farm is determined. We then describe a continuous time version involving a delay between the time of decision and the implementation.

1. Introduction

A higher penetration level of the wind energy into electric power systems plays a part in the reduction of CO₂ emissions. In the meantime, traditional operational management of power systems is transformed by taking into consideration this fluctuating and intermittent resource. Smart grids and storage systems have been developed to overcome these challenges.

For wind power plants, storage is a straightforward solution to reduce renewable variability. It can be used to store electricity in hours of high production and inject electricity in the grid later on. The performance of the operational management can be therefore improved by considering simple charge-discharge plans based on short-term forecasts of the renewable production [1]. For instance, optimal management of wind farms associated with hydropower pumped storage showing economic benefit and increasing the controllability have been studied in [2, 3, 4]. Other examples are the sizing of a distributed battery in order to provide frequency support for a grid-connected wind farm [5] and the optimal operation of a wind farm equipped with a storage unit [6, 7].

For the specific case of isolated systems, which is the aim of our paper, it is necessary to think about distributed energy storage as battery [8], ultra-capacitors [9], or flywheels [10]. In this setting, the question of economic viability in isolated islands without additional reserves arises. Here, the storage unit allows wind farms to respect the scheduled production.

The storage costs will represent a large part of the overnight capital costs and motivate the different researches on storage. Generally the sizing of the storage device is reduced to a minimization problem of the fixed and variable costs of the storage and its application (see [11, 12], for a complete analysis of the cash flow of the storage unit).

In this paper, we present a simplified model in discrete time, in which the commitment is for the next period. We solve an optimal control problem to define the optimal bidding decision. The mathematical setting of the problem is described in Section 2. The main result is detailed in Section 3. Application on a real dataset is described in Section 4. The continuous version of the problem is also described in Section 5. A conclusion ends the paper.

2. Setting of the problem

2.1 General description

In our problem, the manager has to announce an energy production to be delivered to the next period. Considering the $k^{\mathit{th}}$ period, we may think that the announcement is made at the beginning of the period and the delivery at the end of the period. Of course the real delivery will be split along the ${\left(k+1\right)}^{\mathit{th}}$ period. This splitting will be omitted in this stylized model. It is convenient to consider the full delivery at the end of the $k^{\mathit{th}}$ period which is the beginning of the ${\left(k+1\right)}^{\mathit{th}}$ period. So, at the beginning of the $k^{\mathit{th}}$ period (day or hour), the manager commits himself to deliver $v_k$ units of energy (kWh or MWh). To simplify, we discard margins of tolerance. To decide, he knows the amount of energy stored in the battery, called $y_k$. The second element concerns the windfarm. The energy produced by the windfarm is a stochastic process $Z_k$. More precisely, $Z_k$ is the energy produced during the ${\left(k-1\right)}^{\mathit{th}}$ period, which we consider to be available at the beginning of the $k^{\mathit{th}}$ period. So $Z_k$ and all previous values are known. However to fulfill his commitment, the manager will rely on $Z_{k+1}$, the energy produced during the $k^{\mathit{th}}$ period, which we consider, with our convention, to occur at the end of the $k^{\mathit{th}}$ period, which is the beginning of the ${\left(k+1\right)}^{\mathit{th}}$ period. So it is not known by the manager, when he takes his decision. We model the process $Z_k$ as a Markov chain with transition probability density $f_k\left(\zeta |z\right).$ A key issue concerns the choice of this density which is discussed in the application in Section 4. Precisely, although formally

$$\text{Prob}\left(Z_{k+1}=\zeta |Z_k=z\right)=f_k\left(\zeta |z\right)$$E1

In fact, $Z_k$ is obtained through the power law, operating on another Markov chain, the wind speed (see [13, 14] for examples of such Markov chains). We denote by $F_k\left(\zeta |z\right)$ the CDF of the transition probability. We also set ${\overline{F}}_k\left(\zeta |z\right)=1-F_k\left(\zeta |z\right).$

In the language of stochastic control, the decision $v_k$ (control) is measurable with respect to ${\mathcal{F}}^k=\sigma \left(Z_1,\cdots ,Z_k\right).$ The storage $y_k$ is also adapted to ${\mathcal{F}}^k.$ The evolution of $y_k$ is defined by the equation

$$y_{k+1}={\left(min\left(M,y_k+Z_{k+1}-v_k\right)\right)}^{+}$$E2

Indeed, the available energy at the end of the $k^{\mathit{th}}$ period is $y_k+Z_{k+1}.$ If this quantity is smaller than $v_k$, then the manager cannot fulfill his commitment; he delivers what he has, namely, $y_k+Z_{k+1}$; and the storage becomes empty. If the available energy $y_k+Z_{k+1}$ is larger than $v_k,$ then the manager delivers his commitment $v_k$ and tries to store the remainder $y_k+Z_{k+1}-v_k.$ This is possible only when this quantity is smaller than $M,$ which represents the maximum storage of the battery. If $y_k+Z_{k+1}-v_k>M,$ then he charges up to $M,$ and the quantity $y_k+Z_{k+1}-v_k-M$ is lost (given free to the grid). This results in formula (2). In this equation, we do not consider the constraint to keep a minimum reserve. We also are considering the battery as a reservoir of kWh, which we can reduce or increase as done in inventory of goods. Finally, we neglect the losses in the battery.

2.2 The payoff

We want now to write the payoff to be optimized. During the period $k$, if the manager delivers his commitment $v_k$, he receives the normal income ${\mathit{pv}}_k.$ If he fails and delivers only $y_k+Z_{k+1}<v_k,$ there is a penalty. In the sequel, we have chosen the following penalty $\varpi {\left(y_k+Z_{k+1}-v_k\right)}^{-},$ which is common in inventory theory. The parameter $\varpi$ can be adjusted, for instance, to compare with the conditions on the spot market.

In our set up, the pair $y_k,Z_k$ is a Markov chain. So we have a controlled Markov chain and the state is two-dimensional. We introduce a discount factor denoted by $\alpha$, which is useful if we consider an infinite horizon. We can have $\alpha =1$, otherwise. Initial conditions are given at time $n$ and denoted by $x,z$. We call $V=\left(v_{n,}\cdots v_N\right)$ the control, where $N$ is the horizon. Finally, we want to maximize the functional

$$J_n\left(x,z;V\right)=\sum _{k=n}^N{\alpha }^{k-n}E\left[pmin\left(v_k,y_k+Z_{k+1}\right)-\varpi {\left(y_k+Z_{k+1}-v_k\right)}^{-}\right]$$E3

3. Dynamic programming

3.1 Bellman equation

The value function is defined by

$$U_n\left(x,z\right)=\underset{V}{\text{sup}}{J}_n\left(x,z;V\right)$$E4

Writing

$$min\left(v_k,y_k+Z_{k+1}\right)=v_k-{\left(y_k+Z_{k+1}-v_k\right)}^{-}$$

we get also

$$J_n\left(x,z;V\right)=\sum _{k=n}^N{\alpha }^{k-n}E\left[{\mathit{pv}}_k-\left(p+\varpi \right){\left(y_k+Z_{k+1}-v_k\right)}^{-}\right]$$E5

We can then write Bellman equation

$$\begin{array}{l}{U}_n\left(x,z\right)=\underset{v>0}{\text{sup}}\{\mathit{pv}+E[-\left(p+\varpi \right){\left(x+Z_{n+1}-v\right)}^{-}+\\ +\alpha U_{n+1}\left({\left(min\left(M,x+Z_{n+1}-v\right)\right)}^{+},Z_{n+1}\right)\mid Z_n=z\left]\right\}\end{array}$$E6

It is convenient to make the change of variables $v-x=y$ and obtain

$$\begin{array}{l}{U}_n\left(x,z\right)=\mathit{px}+\underset{x+y>0}{\text{sup}}\{\mathit{py}+E[-\left(p+\varpi \right){\left(Z_{n+1}-y\right)}^{-}+\\ +\alpha U_{n+1}\left({\left(min\left(M,Z_{n+1}-y\right)\right)}^{+},Z_{n+1}\right)\mid Z_n=z\left]\right\}\end{array}$$E7

with final equation

$$U_{N+1}\left(x\right)=0$$

We have $x\in \left[0,M\right]$ and $z>0.$

3.2 Main result

We state the following proposition:

Proposition 1. The solution of (7) is of the form

$$U_n\left(x,z\right)=\mathit{px}+K_n\left(z\right)$$E8

Proof. For $n=N,$ we have

$$U_N\left(x,z\right)=\mathit{px}+\underset{x+y>0}{\text{sup}}\left\{\mathit{py}-\left(p+\varpi \right)E\left[{\left(Z_{N+1}-y\right)}^{-}|Z_N=z\right]\right\}$$

Consider the function

$${\phi }_{N+1}\left(y\right)=\mathit{py}-\left(p+\varpi \right)E\left[{\left(Z_{N+1}-y\right)}^{-}|Z_N=z\right]$$

then, for $y<0,$ we have ${\phi }_{N+1}\left(y\right)=\mathit{py}.$ It is monotone increasing, so the maximum cannot be reached at a point $y<0.$ It follows that (8) is satisfied at $n=N,$ with

$$K_N\left(z\right)=\underset{y>0}{\text{sup}}\left\{\mathit{py}-\left(p+\varpi \right)E\left[{\left(Z_{N+1}-y\right)}^{-}|Z_N=z\right]\right\}$$E9

We have, for $y>0$

$${\phi }_{N+1}\left(y\right)=\mathit{py}-\left(p+\varpi \right)\underset{0}{\overset{y}{\int }}\left(y-\zeta \right)f\left(\zeta |z\right)\mathit{d\zeta }$$

and ${\phi }_{N+1}\left(y\right)$ is concave, since

$$\begin{array}{c}{\phi }_{N+1}^{\prime }\left(y\right)=p-\left(p+\varpi \right)F\left(y|z\right)\\ \phi {"}_{N+1}\left(y\right)=-\left(p+\varpi \right)f\left(y|z\right)<0\end{array}$$

and ${\phi }_{N+1}^{\prime }\left(0\right)=p,{\phi }_{N+1}^{\prime }\left(+\infty \right)=-\varpi .$ So the maximum is uniquely defined. Assume now (8) for $n+1.$ We insert it in (7) to obtain

$$\begin{array}{l}{U}_n\left(x,z\right)=\mathit{px}+\underset{x+y>0}{\text{sup}}\{\mathit{py}+E[-\left(p+\varpi \right){\left(Z_{n+1}-y\right)}^{-}+\\ +\mathit{\alpha p}{\left(min\left(M,Z_{n+1}-y\right)\right)}^{+}\mid Z_n=z\left]\right\}+\mathit{\alpha E}\left[K_{n+1}\left(Z_{n+1}\right)|Z_n=z\right]\end{array}$$

Consider now the function

$$\begin{array}{l}{\phi }_{n+1}\left(y|z\right)=\mathit{py}+E[-\left(p+\varpi \right){\left(Z_{n+1}-y\right)}^{-}+\\ +\mathit{\alpha p}{\left(min\left(M,Z_{n+1}-y\right)\right)}^{+}\mid Z_n=z]\end{array}$$

For $y<0,$ it reduces to

$${\phi }_{n+1}\left(y|z\right)=\mathit{py}+\mathit{\alpha p}\left[M{\overline{F}}_{n+1}\left({\left(y+M\right)}^{+}\right)+\underset{0}{\overset{{\left(y+M\right)}^{+}}{\int }}\left(\zeta -y\right)f_{n+1}\left(\zeta |z\right)\mathit{d\zeta }\right]$$

and

$$\begin{array}{c}{\phi }_{n+1}^{\prime }\left(y|z\right)=p-\alpha {\mathit{pF}}_{n+1}\left({\left(y+M\right)}^{+}\right)\\ \ge p-\alpha {\mathit{pF}}_{n+1}\left(M\right)>0\end{array}$$

and thus cannot reach a maximum for $y<0.$ Considering $y>0,$ we have

$$\begin{array}{l}{\phi }_{n+1}\left(y|z\right)=\mathit{py}-\left(p+\varpi \right)\underset{0}{\overset{y}{\int }}\left(y-\zeta \right)f_{n+1}\left(\zeta |z\right)\mathit{d\zeta }+\\ +\mathit{\alpha p}\left[M{\overline{F}}_{n+1}\left(y+M\right)+\underset{y}{\overset{y+M}{\int }}\left(\zeta -y\right)f_{n+1}\left(\zeta |z\right)\mathit{d\zeta }\right]\end{array}$$

Again, this function is concave and

$${\phi }_{n+1}^{\prime }\left(y|z\right)=p-\left(p\left(1-\alpha \right)+\varpi \right)F_{n+1}\left(y|z\right)-\alpha {\mathit{pF}}_{n+1}\left(y+M|z\right)$$

$${\phi }_{n+1}^{″}\left(y|z\right)=-\left(p\left(1-\alpha \right)+\varpi \right)f_{n+1}\left(y|z\right)-\mathit{\alpha p}{f}_{n+1}\left(y+M|z\right)<0$$

$${\phi }_{n+1}^{\prime }\left(0|z\right)=p-\alpha {\mathit{pF}}_{n+1}\left(M|z\right)>0$$

$${\phi }_{n+1}^{\prime }\left(+\infty |z\right)=-\varpi$$

and the property (8) is proven with

$$\begin{array}{l}{K}_n\left(z\right)=\alpha \underset{0}{\overset{+\infty }{\int }}{K}_{n+1}\left(\zeta \right)f_{n+1}\left(\zeta |z\right)\mathit{d\zeta }+\underset{y>0}{max}\{\mathit{py}-\left(p+\varpi \right)\underset{0}{\overset{y}{\int }}\left(y-\zeta \right)f_{n+1}\left(\zeta |z\right)\mathit{d\zeta }+\\ +\mathit{\alpha p}\left[M{\overline{F}}_{n+1}\left(y+M\right)+\underset{y}{\overset{y+M}{\int }}\left(\zeta -y\right)f_{n+1}\left(\zeta |z\right)\mathit{d\zeta }\right]\}\end{array}$$E10

The proof is completed. $■$

3.3 Optimal feedback

We define by $S_n\left(z\right)$ the point at which ${\phi }_{n+1}\left(y|z\right)$ attains its maximum. It is positive and uniquely defined. It follows that the optimal feedback in Bellman equation (6) is ${\widehat{v}}_n\left(x,z\right)=x+S_n\left(z\right)$ and $S_n\left(z\right)$ is the unique solution of

$$\begin{array}{c}p-\left(p\left(1-\alpha \right)+\varpi \right)F_{n+1}\left(S_n|z\right)-\alpha {\mathit{pF}}_{n+1}\left(S_n+M|z\right)=0\\ p-\left(p+\varpi \right)F_{N+1}\left(S_N|z\right)=0\end{array}$$E11

The recursion (10) writes

$$\begin{array}{c}{K}_n\left(z\right)=\alpha \underset{0}{\overset{+\infty }{\int }}{K}_{n+1}\left(\zeta \right)f_{n+1}\left(\zeta |z\right)\mathit{d\zeta }+p\underset{0}{\overset{S_n\left(z\right)}{\int }}{\overline{F}}_{n+1}\left(\zeta |z\right)\mathit{d\zeta }+\mathit{\alpha p}\underset{S_n\left(z\right)}{\overset{S_n\left(z\right)+M}{\int }}{\overline{F}}_{n+1}\left(\zeta |z\right)\mathit{d\zeta }\\ -\varpi \underset{0}{\overset{S_n\left(z\right)}{\int }}{F}_{n+1}\left(\zeta |z\right)\mathit{d\zeta }\\ K_N\left(z\right)=p\underset{0}{\overset{S_N\left(z\right)}{\int }}{\overline{F}}_{N+1}\left(\zeta |z\right)\mathit{d\zeta }-\varpi \underset{0}{\overset{S_N\left(z\right)}{\int }}{F}_{N+1}\left(\zeta |z\right)\mathit{d\zeta }\end{array}$$E12

It is worth emphasizing that the function $K_n\left(z\right)$ increases with $M$, as can be expected. The feedback has an easy interpretation. The bidding is the level of inventory plus a fixed amount depending on the last value of energy produced by the turbine. It is interesting to note that the quantity $S_n\left(z\right)$ decreases with $M.$ This is not so obvious. Clearly, the larger $M,$ the better it is, as captured by the increase of $K_n\left(z\right).$ We can understand why $S_n\left(z\right)$ decreases with $M,$ as follows: When $M$ is large, the risk of wasting energy by lack of storage is reduced, so it makes sense to focus on the other risk, to overbid and pay the penalty. Hence it makes sense to reduce the bid.

4. Application

We describe in this section an application on a wind farm project financed by EREN on a French island with national tender process. First we set the energy price $p=230$ EUR/MWh (Official Journal from March 8, 2013) and the discount factor $\alpha =1$. In this first application, $N=48$ which is the number of periods of 30 min in a day.

In the sequel, we have chosen the penalty $\varpi {\left(y_k+Z_{k+1}-v_k\right)}^{-}$ which is common in inventory theory. The parameter $\varpi$ can be adjusted.

Some analysts would prefer the penalty $\frac{p}{2}\left(y_k+Z_{k+1}\right){11}_{v_k>y_k+Z_{k+1}}$. This is rather strange, because it is fixed, whatever be the level of failure. If the failure is very small, the damage is not big, and still the penalty is high. Conversely, if the failure is big, the damage is big, and still the penalty does not change. Even more surprising, for a given level of commitment, the bigger the failure, the lower the penalty.

The production over a period of 30 min is presented on Figure 1. It is worth mentioning that we used directly the data proposed from July 26, 2005 to March 9, 2008 captured by a measurement mast.

Figure 1.

Histogram of the production over a period of 30 min.

For this first application, stationary law is considered as Gaussian. Mean and variance of the model are similar to those of the empirical distribution in Figure 1. This model allows to construct closed-form cumulative distribution function $F_k$. One-step forecasting error is 24% of the mean and 11% of the nominal power.

But this process does not take into account the stylized facts of the production on a period of 30 min (positive values below nominal power limit, atom for zero production, intraday seasonality, etc.). Consequently, in the optimal control problem, we use the corresponding truncated Gaussian distribution (between 0 and 7 MWh).

Finally, the penalty is fixed (geometrically) to $\varpi =\frac{3}{4}p$.

With these assumptions, the payoff with respect to the size of the storage is given in Figure 2 for an empty storage $x=0$, and $z$ is the average production as initial conditions. Some simple economic models penalizing the size of the battery with its costs would reveal an optimal size of the storage unit between 4 and 6 MWh.

Figure 2.

On the left, daily payoff ${\mathit{U}}_n\left(\mathit{x},\mathit{z}\right)$ in terms of the size of the storage $\mathit{M}$ for $\mathit{\alpha }=\mathit{0.9}$, $\mathit{p}=\mathit{230}$, and $\mathit{\varpi }=\frac{\mathit{3}}{\mathit{4}}\mathit{p}$. Here $\mathit{z}=\mathit{3}$ MWh and the initial storage is empty with $\mathit{x}=\mathit{0}$ MWh. On the right, part of the decision $S_n\left(z\right)$ (from the direct wind production) in terms of the size of the storage $M$.

5. Continuous time version

In the last section, we present a continuous version of the aforementioned problem. This new problem exhibits interesting questions in control theory when there is a delay between the decision and the application of the decision.

5.1 A continuous time model

We model the wind speed by a diffusion

$$\begin{array}{c}\mathit{dz}=g\left(z\right)\mathit{dt}+\sigma \left(z\right)\mathit{dw}\\ z\left(0\right)=z\end{array}$$E13

where $w\left(t\right)$ is a standard one-dimensional Wiener process, built on a probability space $\mathrm{\Omega },\mathcal{A},P$, and we denote by ${\mathcal{F}}^t$ the filtration generated by the Wiener process. This is the unique source of uncertainty in the model. We suppose that the model has a positive solution.

The energy produced per unit of time at time $t$ is $\phi \left(z\left(t\right)\right)$ where the function $\phi$ is the power law. So the energy produced on an interval of time $\left(0,t\right)$ is ${\int }_0^t\phi \left(z\left(s\right)\right)\mathit{ds}.$ We assume that the manager bids for a delivery program with a fixed delay $h.$ In other words, if he decides a level $v\left(t\right)$ per unit of time at time $t,$ the delivery will be at $t+h.$ On the interval of time $\left(0,t\right)$ the delivery is ${\int }_h^{t}v\left(s-h\right)\mathit{ds},$ provided $t>h$, otherwise it is $0.$

Define

$$\begin{array}{c}\eta \left(t\right)=x+\underset{0}{\overset{t}{\int }}\phi \left(z\left(s\right)\right)\mathit{ds}-\underset{h}{\overset{t}{\int }}v\left(s-h\right)\mathit{ds},t>h\\ \eta \left(t\right)=x+\underset{0}{\overset{t}{\int }}\phi \left(z\left(s\right)\right)\mathit{ds},t\le h\end{array}$$E14

which represents the excess of production of energy over the delivery on the interval $\left(0,t\right).$ The initial value $x$ represents the initial amount of energy in the storage unit. We have $0\le x\le M.$ We should have similarly $0\le \eta \left(t\right)\le M$, $\forall t.$ Indeed one cannot deliver more than one produces, and the storage of the excess production is limited by $M.$ This constraint is more complex to handle than in the discrete time case. To simplify we shall treat the constraints with penalties and not impose them. In particular, for coherence, we remove the condition $0\le x\le M,$ which is a purely mathematical extension. The control is the process $v\left(.\right)$, which is adapted to the filtration ${\mathcal{F}}^t$ and just positive. We then define the payoff. The payoff will include the penalty terms and the profit from selling the energy. We assume that the manager can sell his production up to his commitment at a fixed price per unit of time and unit of energy $p.$ We note that $\eta \left(t\right)<0$ captures the situation in which the manager delivers less than his commitment, and there is a penalty for it. The payoff is now written as

$$\begin{array}{c}{J}_{\mathit{xz}}\left(v\left(.\right)\right)=\mathit{pE}\underset{h}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}min\left(\phi \left(z\left(s\right)\right),v\left(s-h\right)\right)\mathit{ds}\\ -\mathit{\varpi E}\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{\eta }^{-}\left(s\right)\mathit{ds}-\mathit{\pi E}\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{\left(\eta \left(s\right)-M\right)}^{+}\mathit{ds}\end{array}$$E15

5.2 Rewriting the payoff functional

Because of the delay, we cannot consider the pair $z\left(t\right),\eta \left(t\right)$ as the state of a dynamic system and apply dynamic programming. In fact, we shall see that it is possible to rewrite the payoff (15) in terms of the pair $z\left(t\right),x\left(t\right)$ with

$$x\left(t\right)=x+\underset{0}{\overset{t}{\int }}\phi \left(z\left(s\right)\right)\mathit{ds}-\underset{0}{\overset{t}{\int }}v\left(s\right)\mathit{ds}$$E16

and the standard reasoning of dynamic programming will become applicable. The first transformation concerns the term

$$I=E\underset{h}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}min\left(\phi \left(z\left(s\right)\right),v\left(s-h\right)\right)\mathit{ds}$$

We have

$$\begin{array}{cc}I& =exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}min\left(\phi \left(z\left(s+h\right)\right),v\left(s\right)\right)\mathit{ds}\\ & =exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}E\left[min\left(\phi \left(z\left(s+h\right)\right),v\left(s\right)\right)|{\mathcal{F}}^s\right]\mathit{ds}\end{array}$$

and we need to compute $E\left[min\left(\phi \left(z\left(s+h\right)\right),v\left(s\right)\right)|{\mathcal{F}}^s\right]$. We remember that $v\left(s\right)$ is ${\mathcal{F}}^s$ measurable and that $z\left(s\right)$ is a stationary Markov process. Let us introduce the transition probability density $m\left(\eta ,s;z\right)$ representing

$$m\left(\eta ,s;z\right)=\text{Prob}\left[z\left(s\right)=\eta |z\left(0\right)=z\right]$$

The function $m\left(\eta ,s;z\right)$ is the solution of Fokker-Planck equation

$$\begin{array}{c}\frac{\partial m}{\partial s}-\frac{1}{2}\frac{{\partial }^2}{\partial {\eta }^2}\left({\sigma }^2\left(\eta \right)m\right)+\frac{\partial }{\partial \eta }\left(g\left(\eta \right)m\right)=0\\ m\left(\eta ,0;z\right)=\delta \left(\eta -z\right)\end{array}$$E17

Then by stationarity of the Markov process $z\left(s\right)$, we can write

$$E\left[min\left(\phi \left(z\left(s+h\right)\right),v\left(s\right)\right)|{\mathcal{F}}^s\right]={\int }_{R}min\left(\phi \left(\xi \right),v\left(s\right)\right)m\left(\xi ,s;z\left(s\right)\right)\mathit{d\xi }$$E18

Therefore

$$I=exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{\int }_{R}min\left(\phi \left(\xi \right),v\left(s\right)\right)m\left(\xi ,s;z\left(s\right)\right)\mathit{d\xi }$$E19

The next transformation concerns

$$\begin{array}{c}\mathit{II}=E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{\eta }^{-}\left(s\right)\mathit{ds}\\ \mathit{II}=E\underset{0}{\overset{h}{\int }}exp-\mathit{\alpha s}{\eta }^{-}\left(s\right)\mathit{ds}+E\underset{h}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{\eta }^{-}\left(s\right)\mathit{ds}\\ ={\mathit{II}}_1+{\mathit{II}}_2.\end{array}$$E20

The first integral does not depend on the control and is $0,$ when $x\ge 0,$ as it will be the case in practice. The second integral is written as

$${\mathit{II}}_{2.}=exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{\eta }^{-}\left(s+h\right)\mathit{ds}$$

We note that

$$\eta \left(s+h\right)=x\left(s\right)+\underset{s}{\overset{s+h}{\int }}\phi \left(z\left(\tau \right)\right)\mathit{d\tau }$$

Recalling the definition of $x\left(s\right),$ see (16). We then can write

$$\begin{array}{cc}E{\eta }^{-}\left(s+h\right)& =EE\left[{\eta }^{-}\left(s+h\right)|{\mathcal{F}}^s\right]=\\ & =\mathit{E\theta }\left(x\left(s\right),z\left(s\right);h\right)\end{array}$$

with

$$\theta \left(x,z;s\right)=E{\left(x+\underset{0}{\overset{s}{\int }}\phi \left(z\left(\tau \right)\right)\mathit{d\tau }\right)}^{-}$$E21

The argument $x$ is a real number and $z\left(0\right)=z$. So

$${\mathit{II}}_2=exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}\mathit{E\theta }\left(x\left(s\right),z\left(s\right);h\right)\mathit{ds}$$

We can also write

$${\mathit{II}}_1=\underset{0}{\overset{h}{\int }}exp-\mathit{\alpha s}\theta \left(x,z;s\right)\mathit{ds}$$

so we have

$$\mathit{II}=\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}\theta \left(x,z;s\right)\mathit{ds}+exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}\mathit{E\theta }\left(x\left(s\right),z\left(s\right);h\right)\mathit{ds}$$E22

We can give a similar formula for the second penalty term

$$\mathit{III}=E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{\left(\eta \left(s\right)-M\right)}^{+}\mathit{ds}$$E23

We introduce the function

$$\chi \left(x,z;s\right)=E{\left(x+\underset{0}{\overset{s}{\int }}\phi \left(z\left(\tau \right)\right)\mathit{d\tau }-M\right)}^{+}$$E24

and we can write

$$\mathit{III}=\underset{0}{\overset{h}{\int }}exp-\mathit{\alpha s}\chi \left(x,z;s\right)\mathit{ds}+exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}\mathit{E\chi }\left(x\left(s\right),z\left(s\right);h\right)\mathit{ds}$$E25

Combining results, we obtain the formula

$$\begin{array}{cc}{J}_{\mathit{xz}}\left(v\left(.\right)\right)& =-\varpi \underset{0}{\overset{h}{\int }}exp-\mathit{\alpha s}\theta \left(x,z;s\right)\mathit{ds}-\pi \underset{0}{\overset{h}{\int }}exp-\mathit{\alpha s}\theta \left(x,z;s\right)\mathit{ds}+\\ & +exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}\left(p{\int }_{R}min(\phi \left(\xi \right),v\left(s\right)\left)m\right(\xi ,s;z\left(s\right))\mathit{d\xi }-\right.\\ & -\left.\left(\mathit{\varpi \theta }\right(x\left(s\right),z\left(s\right);h)+\mathit{\pi \chi }(x\left(s\right),z\left(s\right);h\left)\right)\right)\end{array}$$E26

or

$$J_{\mathit{xz}}\left(v\left(.\right)\right)={\rho }_h\left(x,z\right)+exp-\mathit{\alpha h}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{l}_h\left(x\left(s\right),z\left(s\right),v\left(s\right)\right)\mathit{ds}$$E27

with

$$\begin{array}{c}{l}_h\left(x,z,v\right)=p{\int }_{R}min\left(\phi \left(\xi \right),v\right)m\left(\xi ,s;z\right)\mathit{d\xi }-\\ -\left(\mathit{\varpi \theta }\left(x,z;h\right)+\mathit{\pi \chi }\left(x,z;h\right)\right)\end{array}$$E28

The stochastic control problem becomes

$$\begin{array}{l}\frac{\mathit{dx}}{\mathit{dt}}=\phi \left(z\left(t\right)\right)-v\left(t\right)\\ \mathit{dz}=g\left(z\right)\mathit{dt}+\sigma \left(z\right)\mathit{dw},x\left(0\right)=x,z\left(0\right)=z\\ \underset{v\left(.\right)}{\text{sup}}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{l}_h\left(x\left(s\right),z\left(s\right),v\left(s\right)\right)\mathit{ds}\end{array}$$E29

which is a standard stochastic control problem. To avoid singularities, we impose a bound on the control $v\left(t\right).$ So we impose

$$0\le v\left(t\right)\le \phi \left(z\left(t\right)\right)+\overline{v}\mathrm{a}.\mathrm{s}.$$E30

in which $\overline{v}$ is a fixed constant.

5.3 Dynamic programming

Let us define the value function

$$\mathrm{\Phi }\left(x,z\right)=\underset{\left\{v\left(.\right)|0\le v\left(t\right)\le \phi \left(z\left(t\right)\right)+\overline{v}\right\}}{\text{sup}}E\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}{l}_h\left(x\left(s\right),z\left(s\right),v\left(s\right)\right)\mathit{ds}$$E31

Then it is easy to write the Bellman equation for the value function, namely,

$$\begin{array}{c}\alpha \mathrm{\Phi }=\phi \left(z\right)\frac{\mathrm{\partial \Phi }}{\partial x}+g\left(z\right)\frac{\mathrm{\partial \Phi }}{\partial z}+\frac{1}{2}{\sigma }^2\left(z\right)\frac{{\partial }^2\mathrm{\Phi }}{\partial z^2}+\\ +\underset{0\le v\le \phi \left(z\right)+\overline{v}}{\text{sup}}\left[l_h\left(x,z,v\right)-v\frac{\mathrm{\partial \Phi }}{\partial x}\right]\end{array}$$E32

A priori $x\in R$ and $z\ge 0$ (we may assume that $\sigma \left(0\right)=0).$ We can add growth conditions to get a problem which is well posed. The optimal feedback ${\widehat{v}}_h\left(x,z\right)$ is obtained by taking the sup in the bracket, with respect to the argument $v.$.

5.4 The case $h=0$

The case $h=0$ has a trivial solution. We note that

$$l_0\left(x,z,v\right)=pmin\left(\phi \left(z\right),v\right)-\varpi x^{-}-\pi {\left(x-M\right)}^{+}$$E33

The optimal feedback is then

$${\widehat{v}}_0\left(x,z\right)=\left|\begin{array}{cc}0& \text{if}x<0\\ \phi \left(z\right)& \text{if}0\le x\le M\\ \phi \left(z\right)+\overline{v}& \text{if}x>M\end{array}\right.$$E34

so (32) becomes

$$\alpha \mathrm{\Phi }=\mathit{\varpi x}+\phi \left(z\right)\frac{\mathrm{\partial \Phi }}{\partial x}+g\left(z\right)\frac{\mathrm{\partial \Phi }}{\partial z}+\frac{1}{2}{\sigma }^2\left(z\right)\frac{{\partial }^2\mathrm{\Phi }}{\partial z^2},\text{if}x<0$$E35

$$\alpha \mathrm{\Phi }=\mathit{p\phi }\left(z\right)+g\left(z\right)\frac{\mathrm{\partial \Phi }}{\partial z}+\frac{1}{2}{\sigma }^2\left(z\right)\frac{{\partial }^2\mathrm{\Phi }}{\partial z^2},\text{if}0<x<M$$E36

$$\alpha \mathrm{\Phi }=\mathit{p\phi }\left(z\right)-\pi \left(x-M\right)-\overline{v}\frac{\mathrm{\partial \Phi }}{\partial x}+g\left(z\right)\frac{\mathrm{\partial \Phi }}{\partial z}+\frac{1}{2}{\sigma }^2\left(z\right)\frac{{\partial }^2\mathrm{\Phi }}{\partial z^2},\text{if}x>M$$E37

The solution for $0<x<M$ does not depend on $x$ and has an easy probabilistic interpretation

$$\mathrm{\Phi }\left(x,z\right)=\mathrm{\Phi }\left(z\right)=\mathit{pE}\underset{0}{\overset{+\infty }{\int }}exp-\mathit{\alpha s}\phi \left(z\left(s\right)\right)\mathit{ds},$$E38

For $x<0$ or $x>M$, we have to solve parabolic problems, considering $x$ as a time, backward and forward. We define the values $\mathrm{\Phi }\left(0,z\right)$ and $\mathrm{\Phi }\left(M,z\right)$ by $\mathrm{\Phi }\left(z\right)$ defined by (38).

5.5 Analytical problems for $\theta$ and $\chi$

The functions $\theta \left(x,z,s\right)$ and $\chi \left(x,z,s\right)$ are solutions of the parabolic PDE

$$\begin{array}{c}-\frac{\partial \theta }{\partial s}+\phi \left(z\right)\frac{\partial \theta }{\partial x}+g\left(z\right)\frac{\partial \theta }{\partial z}+\frac{1}{2}{\sigma }^2\left(z\right)\frac{{\partial }^2\theta }{\partial z^2}=0\\ \theta \left(x,z,0\right)=x^{-}\end{array}$$E39

$$\begin{array}{c}-\frac{\partial \chi }{\partial s}+\phi \left(z\right)\frac{\partial \chi }{\partial x}+g\left(z\right)\frac{\partial \chi }{\partial z}+\frac{1}{2}{\sigma }^2\left(z\right)\frac{{\partial }^2\chi }{\partial z^2}=0\\ \chi \left(x,z,0\right)={\left(x-M\right)}^{+}\end{array}$$E40

This allows to compute $\theta \left(x,z,h\right)$ and $\chi \left(x,z,h\right)$.

6. Conclusions

The problem of the optimal delivery for wind energy in some future time with a storage facility (a battery for instance) is considered. We solve an optimal control problem to define the optimal bidding decision in a simple discrete stochastic problem and apply it to real data. Optimal size of the battery and the overnight costs are discussed.

Acknowledgments

The research is supported by EREN-GROUPE, the grant NSF-DMS 161 2880 and the grant ANR-15-CE05-0024.