A Joint Electricity Market-Clearing Mechanism for Flexible Ramping Products with a Convex Spot Market Model

By inergency On Mar 13, 2024

In the present study, we advance a mathematical depiction of a CESM, which will be based on the linear quotation approach in [29]. The time axis is discretized in this model, with a scheduling horizon of 24 h. We assign the index

h

to represent the hour of the day and construct the set

H : = \{h : h = 1, 2, \dots, 24\}

. We partition each hour h into K equal-length intra-hourly subperiods and give an index to each subperiod to signify its temporal granularity. At the same time, we define the set

T : = \{t : t = 1, 2 \dots K, \dots, 2 K, 2 K + 1, \dots, 24 \cdot K\}

. Taking the DAM clearing time scale (1 h) as an example, intra-hourly subperiods

K = 1

, and temporal granularity set

T : = \{t : t = 1, 2, \dots, 24\}

, we describe the power system bus index n and define the set of buses

N : = \{n : n = 1, 2 \dots N\}

, where the constant

N

denotes the aggregate number of buses within the system. At the same time, we describe the generator index g and define the set

G : = \{g : g = 1, 2 \dots G\}

. The constant

G

represents the number of generator sets in the system. To illustrate the generator set on the bus

i

, we describe the set

Ω (i)

. (Take the example of generator one and generator two on bus one,

Ω (1) : = \{g : g = 1, 2\}

2.2. CESM Vectorized Structure and Solution Methods

The CESM problem is a mixed-integer second-order cone problem. To evaluate the clearing price in the markets, we use the standard practice implemented by most independent operators [30]. First, we run the CESM model to obtain the optimal solution in binary variables, indicating the DPG’s startup and shutdown status [31]. We fix the DPG’s startup and shutdown status when we receive a continuous second-order cone programming problem (security-constrained economic dispatch). We utilize the Karush—Kuhn—Tucker (KKT) condition to solve the new convex problem.

To visualize the mathematical principles of the proposed model more intuitively, the method for solving for dual variables is described as follows.

We set the variables vector as follows:

$x = ({\{P_{γ, τ}\}}_{\forall γ \in G, τ \in T}, {\{Q_{γ, τ}\}}_{\forall γ \in G, τ \in T}, {\{S P C_{γ, τ}\}}_{\forall γ \in G, τ \in T}, {\{C_{α, β, τ}, S_{α, β, τ}\}}_{\forall α, β \in N, τ \in T})$

The vectorized CESM convex model can be expressed as

$M i n (b^{⊤} x + e^{⊤} c)$

(22)

$s . t . {\underline{f}}_{θ} \leq {‖A_{θ} x + d_{θ}‖}_{2} - a_{θ}^{⊤} x \leq {\bar{f}}_{θ}$

(23)

$k_{ξ}^{T} x = e^{⊤} h_{ξ}$

(24)

The problem’s parameters are $A \in ℝ^{(2 n^{2} + G) \times (2 n^{2} + G)}, d \in ℝ^{(2 n^{2} + G)}, a \in ℝ^{(2 n^{2} + G)}$ , $\bar{f}, \underline{f} \in ℝ, k \in ℝ^{(2 n^{2} + 3 G)}, h \in ℝ^{(2 n^{2} + 3 G)}, e = {[1, \dots, 1]}^{T}$ .

We derive the solution procedure based on the Karush—Kuhn—Tucker condition. Let

f (x) = b^{⊤} x + e^{⊤} c, g (x) = {‖A_{θ} x + d_{θ}‖}_{2} - a_{θ}^{⊤} x, h (x) = k_{ξ}^{T} x - e^{⊤} h_{ξ}

. The following function can describe the CESM vectorized convex model:

$\begin{array}{l} M i n & f (x) \\ s . t . & g (x) \leq 0 \\ h (x) = 0 \end{array}$

(25)

By slacking the problem in Equation (27), the slacked optimization problem can be obtained:

$\begin{array}{l} \min & f (x) \\ s . t . \end{array} h (x) = 0 g (x) - l - {\underline{f}}_{θ} = 0 g (x) + u - {\bar{f}}_{θ} = 0 l, u \geq 0$

(26)

In the optimization problem in Equation (28), the slack variable is $l = {[l_{1}, \dots, l_{c 1}]}^{T}, u = {[u_{1}, \dots, u_{c 2}]}^{T}$ .

By transforming the objective function, the barrier function is obtained:

$\begin{matrix} \min & f (x) - μ \sum_{θ = 1}^{c} \ln (l_{θ}) - μ \sum_{θ = 1}^{c} \ln (u_{θ}) \\ s . t . \end{matrix} h (x) = 0 g (x) - l - {\underline{f}}_{θ} = 0 g (x) + u - {\bar{f}}_{θ} = 0 μ \geq 0$

(27)

The Lagrange function of the optimization problem in Equation (29) is given by

$L = f (x) - h (x) λ - [g (x) - l - {\underline{f}}_{θ}] z - [g (x) + u - {\bar{f}}_{θ}] w - μ \sum_{θ = 1}^{c} \ln (l_{θ}) - μ \sum_{θ = 1}^{c} \ln (u_{θ})$

(28)

In Equation (30),

λ = {[λ_{1}, \dots, λ_{ξ}]}^{T}, z = {[z_{1}, \dots, z_{θ}]}^{T} \geq 0, w = {[w_{1}, \dots, w_{η}]}^{T} \geq 0

is the Lagrange multiplier. The necessary conditions for the existence of the extremum of the problem is that the partial derivatives of the Lagrange function for all variables and multipliers are zero. Therefore, we obtain the first-order necessary conditions as follows:

$\{\begin{array}{l} L_{x} = \frac{\partial L}{\partial x} = \nabla f (x) - \nabla h (x) λ - \nabla g (x) (z + w) = 0 \\ L_{y} = \frac{\partial L}{\partial y} = h (x) = 0 \\ L_{z} = \frac{\partial L}{\partial z} = g (x) - l - {\underline{f}}_{θ} = 0 \\ L_{w} = \frac{\partial L}{\partial w} = φ (x) + u - {\bar{f}}_{θ} = 0 \\ L_{l} = \frac{\partial L}{\partial l} = z - μ_{1} L^{- 1} e = 0 \Rightarrow L_{l}^{μ} = Z l - μ e = 0 \\ L_{u} = \frac{\partial L}{\partial u} = - w - μ_{2} U^{- 1} e = 0 \Rightarrow L_{u}^{μ} = W u + μ e = 0 \end{array}$

(29)

In Equation (31),

Z = d i a g (z), W = d i a g (w)

. With

L_{l}^{μ} = 0

and

L_{u}^{μ} = 0

, and we obtain the equation

$μ = l^{T} z - u^{T} w$

(30)

The duality gap is defined as

$G a p = l^{T} z - u^{T} w$

(31)

The authors of [32] proved that under certain conditions, if

x *

is a solution to the optimization problem in Equation (28), then

x (μ)

is a solution to the optimization problem in Equation (29) when

μ

is fixed. Then, when

G a p \to 0, μ \to 0

x (μ)

converges to

x *

In Equation (34),

σ \in (0, 1)

is called the center parameter. From [33], when the parameter μ in the objective function is taken as in Equation (32), the convergence of the algorithm is poor, and it is found that when the parameter

μ

in the objective function is taken according to Equation (33),

σ

is generally taken to be 0.1, and good convergence results can be obtained most often.

The necessary condition in Equation (29) for the extremum is a system of nonlinear equations which the Newton—Raphson method can solve. After linearizing Equation (29), we obtain the modified equation:

$\{\begin{cases} - L_{x} = & (\nabla^{2} f (x) - \nabla^{2} h (x) λ - \nabla^{2} g (x) z - \nabla^{2} g (x) w)) Δ x \\ - \nabla h (x) Δ λ - \nabla g (x) Δ z - \nabla g (x) Δ w) \\ - L_{y} = & \nabla h {(x)}^{T} Δ x \\ - L_{z} = & \nabla g {(x)}^{T} Δ x - Δ l \\ - L_{w} = & \nabla g {(x)}^{T} Δ x + Δ u \\ - L_{l}^{μ} = & Z Δ l + L Δ z \\ - L_{u}^{μ} = & W Δ u + U Δ w \end{cases}$

(34)

$\{\begin{cases} [\begin{matrix} H (\cdot) & \nabla h (x) & \nabla g (x) & \nabla g (x) & 0 & 0 \\ \nabla h {(x)}^{T} & 0 & 0 & 0 & 0 & 0 \\ \nabla g {(x)}^{T} & 0 & 0 & 0 & - I & 0 \\ \nabla g {(x)}^{T} & 0 & 0 & 0 & 0 & I \\ 0 & 0 & L & 0 & Z & 0 \\ 0 & 0 & 0 & U & 0 & W \end{matrix}] [\begin{matrix} Δ x \\ Δ λ \\ Δ z \\ Δ w \\ Δ l \\ Δ u \end{matrix}] = [\begin{matrix} L_{x} \\ - L_{λ} \\ - L_{z} \\ - L_{w} \\ - L_{l}^{μ} \\ - L_{u}^{μ} \end{matrix}] \\ H (\cdot) = - \nabla^{2} f (x) + \nabla^{2} h (x) λ + \nabla^{2} g (x) z + \nabla^{2} g (x) w) \end{cases}$

(35)

In Equation (35), the unit diagonal matrix and the column vectors that make it up are represented by $I = d i a g ([1, 1, \dots, 1]) = [I_{1}, I_{2}, \dots, I_{Γ}]$ .

Taking $g (x) = [\begin{matrix} g_{1} (x) & g_{2} (x) & \dots & g_{θ} (x) \end{matrix}]$ as an example, the vector function is calculated as follows.

The Jacobian matrix is obtained by

$\nabla g (x) = {[\begin{matrix} g_{1} (x) \\ g_{2} (x) \\ ⋮ \\ g_{θ} (x) \end{matrix}]}^{'} = [\begin{matrix} I_{1}^{T} A_{1}^{2} x + I_{1}^{T} (2 d_{1}^{T} A_{1} - a_{1}^{T}) & \dots & I_{θ}^{T} A_{1}^{2} x + I_{θ}^{T} (2 d_{1}^{T} A_{1} - a_{1}^{T}) \\ \begin{matrix} I_{1}^{T} A_{2}^{2} x + I_{1}^{T} (2 d_{2}^{T} A_{2} - a_{2}^{T}) \\ ⋮ \end{matrix} & \begin{array}{l} \dots \\ ⋱ \end{array} & \begin{matrix} I_{θ}^{T} A_{2}^{2} x + I_{θ}^{T} (2 d_{2}^{T} A_{2} - a_{2}^{T}) \\ ⋮ \end{matrix} \\ I_{1}^{T} A_{θ}^{2} x + I_{1}^{T} (2 d_{θ}^{T} A_{θ} - a_{θ}^{T}) & \dots & I_{θ}^{T} A_{θ}^{2} x + I_{θ}^{T} (2 d_{θ}^{T} A_{θ} - a_{θ}^{T}) \end{matrix}]$

(36)

\nabla^{2} g (x)

is a Hessian matrix. The derivation of the first line of

\nabla^{2} g (x)

is as follows:

$\nabla^{2} g (x) = \nabla g {(x)}^{'} \Rightarrow \nabla^{2} g_{1} (x) ≜ [\begin{matrix} I_{1}^{T} A_{1}^{2} I_{1} & I_{1}^{T} A_{1}^{2} I_{2} & \dots & I_{1}^{T} A_{1}^{2} I_{θ} \\ I_{2}^{T} A_{1}^{2} I_{1} & I_{2}^{T} A_{1}^{2} I_{2} & \dots & I_{2}^{T} A_{1}^{2} I_{θ} \\ ⋮ & ⋮ & ⋱ & ⋮ \\ I_{θ}^{T} A_{1}^{2} I_{1} & I_{θ}^{T} A_{1}^{2} I_{2} & \dots & I_{n}^{T} A_{1}^{2} I_{θ} \end{matrix}]$

(37)

The dual variable ( $λ$ ) of $h (x)$ is obtained by the above solution method. The dual variable in $λ$ corresponding to the bus power balance (BPB) constraint is the clearing price of electricity.