2. State of the Art

2.1. Comprehensive Structure of Energy Distribution Network

Currently, a comprehensive energy distribution network typically comprises various components, including the comprehensive energy system, lithium battery ESS, public power grid, central control unit, and electrical equipment, as depicted in Figure 1.

Figure 1 System structure of the integrated energy distribution network.

Within the structure of the comprehensive energy distribution network, we uphold the principle of restricting power supply to electrical equipment through comprehensive energy sources. This strategy significantly curtails the power requirements of electrical equipment from the grid.

2.2. Application and Optimization Principle of the ESS

The mathematical model of the integrated energy distribution network is illustrated in Figure 2. This system has the capability to supply multiple forms of energy in the form of heat and power to users. The distribution network encompasses a higher-level grid, an ES device, a photovoltaic power generation system, and a natural engine power generation system. The ESS performs peak shaving and valley filling based on time-sharing tariffs to optimize resource utilization and reduce costs. When the power generated by the new energy generation system exceeds the electrical equipment load, the ES device is charged. Conversely, when the power is insufficient, the ES device discharges and additional power is procured from the superior grid.

Figure 2 IES structure.

The optimization model of the ESS is depicted in Figure 3. Initially, the output power of the various power supply systems within the IES is regulated in response to the electrical load demand. Subsequently, the system operation model is optimized to identify the most efficient configuration among multiple options for the comprehensive energy system, improving the system’s economic performance in the shortest time. Simultaneously, an optimal ES configuration model is developed and the ES capacity is determined based on the output power of different power supply networks and time-of-use (TOU) pricing. As a result, it enhances the annual income of the comprehensive energy system while promoting a slower degradation of ES capacity.

Figure 3 Optimization model structure.

2.3. Load Analysis of Electrical Equipment

The electrical equipment represented in the model corresponds to typical high-power household appliances. The key factors influencing the load include the appliance’s physical dimensions, initial charging duration, and service duration. Based on findings from a national survey of electrical equipment, the operational time of these appliances follows a log-normal distribution N(17.8, 3.6), as shown in the following equation:

(1) $\begin{array}{}{f}_s\left(n_1\right)=\left(\begin{array}{c}\frac{1}{{\sigma }_s\sqrt{2\pi }}\exp \left(-\frac{\left(n_1-{\mu }_S\right)}{2{\sigma }_s^2}\right),{\mu }_s-12⩽n_1⩽24,\\ \frac{1}{{\sigma }_s\sqrt{2\pi }}\exp \left(-\frac{\left(n_1+24-{\mu }_s\right)}{2{\sigma }_s^2}\right),0⩽n_1<{\mu }_s-12,\end{array}\right)\end{array}$

where μ_s and σ_S are, respectively, the expected value and standard deviation of the initial charging time of household appliances. n ₁ is the initial charging time of the electrical device.

Simultaneously, it is assumed that the initial charging time of the electrical equipment follows a uniform distribution, as expressed in the following equation:

(2) $\begin{array}{}{f}_s\left(n_2\right)=\text{rand}\left(24,1\right),\end{array}$

where rand(24, 1) represents the random integer generated in the interval [Reference Idris, Sasongko and Kuntjoro1, Reference Walker and Kwon18]. n ₂ is the initial charging time of the electrical device.

According to the results of the national survey on electrical equipment, the use time of electrical equipment follows the lognormal distribution N(45, 0.87). Its probability density function is shown in the following equation:

(3) $\begin{array}{}{f}_D^{\prime }\left(s_1\right)=\frac{1}{s_1{\sigma }_{D_1}\sqrt{2\pi }}\exp \left(\frac{\ln \left(s_1-{\mu }_{D_1}\right)}{2{\sigma }_{D_1}}\right),\end{array}$

where μ_D1 and σ_D1 are the expected and standard deviation of the service time of the electrical equipment, respectively. S ₁ refers to the service duration of electrical devices.

Utilizing the usage guidelines and probability density function provided above, the Monte Carlo algorithm is employed to stochastically sample the initial charging time and service duration of various electrical equipment. Through an analysis of the charging profiles of electrical equipment during different time intervals and the initial state of charge (SOC) of ES batteries, it is possible to subsequently calculate the daily load demand of electrical equipment.

3. Methodology

In this paper, the enhanced multiobjective whale optimization algorithm (WOA) [Reference Nasiri and Khiyabani19] is employed in conjunction with models for ES operation efficiency and ES capacity attenuation. This study is focused on addressing ES issues, particularly concerning the efficiency of ES operations and the decline in ES capacity within the ESS. By comparing various algorithms, we investigate differences in ES capacity degradation time, the average annual income of the comprehensive energy system, and system operation efficiency. Furthermore, we delve into the optimization of the ESS using the enhanced algorithm introduced in this paper.

The paper provides a comprehensive methodology for optimizing ESSs through the utilization of the improved WOA. This methodology encompasses enhancements to the algorithm, the management of objective functions, and the modeling of ES efficiency and degradation.

3.1. Improve the WOA

The basic WOA is an optimization algorithm based on the behavior of whales to round up prey. It includes the following three search methods: bubble search, enveloping search, and random mutation. Assume that the initial population size of whales is T and the variable dimension is D. Cetacean individuals in the solution space position can be expressed as ${\overline{I}}_x=\left(i_x^1,i_x^2,\cdots ,i_x^D\right)$ . Each whale corresponds to an optimal configuration of the ESS. Upon discovering the optimal individual location, it represents the optimal solution for configuring the ESS.

3.1.1. Enveloping Search

Enveloping exploration within the whale optimization algorithm emulates the hunting actions of groups of whales, where the leading whale identifies its target and the other whales employ the lead position as a point of reference. They update their positions in a specific manner to create an encircling pattern around the prey. The equation for individual position updates is provided in the following equation:

(4) $\begin{array}{}\left(\begin{array}{c}\overline{I}\left(n+1\right)=\overline{I^{\prime }}\left(n\right)-\overline{G}\cdot \overline{D},\\ \overline{G}=2\overline{g}\cdot \overline{r}-\overline{g},\\ \overline{g}=2-\frac{2n}{T},\\ \overline{D}=\left(\overline{C}\cdot \overline{I^{\prime }}\left(n\right)-\overline{I}\left(n\right)\right),\end{array}\right)\end{array}$

where n denotes the present iteration count. Ī(n + 1) is the updated position of the individual whale. Ī߰(n) is the position of leading whale. Ḡ is a convergence factor. $\overline{D}$ is the distance between the current individual and the leading whale [Reference Liu, Guo, Jiang, Liu, Wang and Ni20]. ḡ linearly diminishes from 2 to 0 with the progression of iterations. $\overline{r}$ and $\overline{C}$ represent random values within the [0, 1] range. T indicates the maximum iteration count.

3.1.2. Bubble Search

Bubble search in the algorithm simulates the foraging behavior of whales as they spiral and release bubbles. In this process, individual particles update their positions using the spiral equation, as indicated in the following equation:

(5) $\begin{array}{}\overline{I}\left(n+1\right)=\overline{D^{\prime }}\cdot e^{hl}\cdot \cos \left(2\pi l\right)+\overline{I^{\prime }}\left(n\right),\end{array}$

where Ī(n + 1) is the updated position of the individual. Ī′(n + 1) is the leading whale position. l denotes a random number falling within the range of [−1, 1]. h serves as a constant for the spiral path. $\overline{D^{\prime }}$ stands for the distance between the current individual and the leading whale.

3.1.3. Random Variation

Individuals also introduce random changes to their positions as they search for food. This updating mechanism expands the search scope and helps avoid getting stuck in local optimal solutions. The individual position update equation is presented as follows:

(6) $\begin{array}{}\left(\begin{array}{c}\overline{I}\left(n+1\right)={\overline{I}}_{\text{rand}}\left(n\right)-\overline{G}\cdot \overline{D},\\ \overline{D}=\left(\overline{C}\cdot {\overline{I}}_{\text{rand}}\left(n\right)-\overline{I}\left(n\right)\right),\end{array}\right)\end{array}$

where Ī(n + 1) is the updated position of the individual. Ī_rand(n) is the location of a random individual whale.

The improved multiobjective WOA builds upon the basic algorithm, introducing cubic chaotic mapping and an individual information exchange mechanism. In addition, it incorporates a nonlinear convergence factor, adaptive weighting, and multiobjective function processing.

(1) Initializing the Population Using Cubic Chaotic Mapping. In the basic WOA, the initial population exhibits high randomness, often starting too far from the optimal solution. This can significantly impede the efficiency of the iterative search. To enhance population diversity and ensure even coverage of the search space, cubic chaotic mapping is employed for dynamic and uniform generation of the initial population. The cubic mapping is expressed in the following equations:

(7) $\begin{array}{}\left(\begin{array}{c}j\left(t+1\right)=4j{\left(t\right)}^3-3j\left(t\right),\\ -1<j\left(t\right)<1,j\left(t\right)\ne 0,\end{array}\right)\end{array}$

(8) $\begin{array}{}{i}_t=L_d+\left(1+j_t\right)\frac{\left(U_d-L_d\right)}{2},\end{array}$

where U_d and L_d are the upper and lower limits of the solution space, respectively. The whale population is generated according to equation (7). The location of the whale is mapped to the solution space according to equation (8) and the location of the whale is initialized.

(2) Introduce Individual Information Exchange Mechanism. In the basic WOA, the whale population conducts the search for the optimal solution by updating their positions autonomously. The individual updating mechanism primarily influences the algorithm’s iteration efficiency. To enhance the algorithm’s convergence rate, the worst solution is updated using a hybrid leapfrog algorithm. An individual information exchange mechanism is introduced to facilitate communication and learning between the weaker individuals and the best individuals. The update equation for the worst solution is provided in the following equation:

(9) $\begin{array}{}\overline{I}={\overline{I}}_m+r\cdot \left({\overline{I}}_h-{\overline{I}}_M\right),\end{array}$

where Ī is the updated individual. Ī_h is the best individual in the population. Ī_M is the worst individual. r is a random number between [0, 1].

For the updated individual Ī, if Ī is better than Ī_M, Ī is used instead of Ī_M. Otherwise, substitute Ī_h for Ī_M.

(3) Increase the Nonlinear Convergence Factor and Adaptive Weight. The search performance of the basic WOA is affected by the convergence factor G, and the value of G is determined by g. Due to the intricacy of the solution process, the conventional linear degradation approach of a value has proven to be inadequate for addressing real-world engineering needs, particularly in high-dimensional nonlinear optimization problems. As a result, this paper introduces an updated method for the traditional convergence factor, and its calculation equation is presented as follows:

(10) $\begin{array}{}g\left(n\right)=2-u_1\cdot e^{u_2\cdot n/N_{\max }}.\end{array}$

In the equation, u₁ and u₂ are nonlinear adjustment parameters and the values in this paper are 0.07 and 9/8π, respectively.

In the optimization process of the WOA, a balance must be struck between the convergence efficiency and the algorithm’s solution precision. During the initial optimization stages, the range of random search is expanded to prevent getting trapped in local optima. In the later stages of optimization, it is essential to enhance the learning capability of the optimal solution to accelerate algorithm convergence. Therefore, an adaptive weight coefficient is introduced to balance the global and local search abilities of the algorithm. After introducing this new convergence factor, equations (4)–(6) are updated as shown in the following equation:

(11) $\begin{array}{}\left(\begin{array}{c}\overline{I}\left(n+1\right)=g\left(n\right)\cdot I_h\left(n\right)-G\cdot D,\\ \overline{I}\left(n+1\right)=g\left(n\right)\cdot I_h\left(n\right)+D\cdot e^{hl}\cos \left(2\pi l\right),\\ \overline{I}\left(n+1\right)=\left(1-g\left(n\right)\right)\cdot I_{\text{rand}}\left(n\right)-G\cdot D.\end{array}\right)\end{array}$

(4) Multiobjective Function Processing. To address the optimal configuration model of the multiobjective function proposed in this paper, the basic WOA is enhanced in the following two aspects.

1. (1) A fast nondominant sort is introduced to calculate the crowding degree and the elite retention strategy.

2. (2) Choose the solution with the highest comprehensive satisfaction based on information entropy. The specific steps are as follows.

   Step 1: Normalize the objective function. The calculation method is shown in the following equation:

   (12) $\begin{array}{}{\mu }_{x,z}=\frac{f_{x,\max }-f_{x,z}}{f_{x,z}-f_{x,\min }},\end{array}$

   where f_x,z is the actual value of the k-th objective function of individual x. μ_x,z is the value after normalization of the k-th objective function of individual x. f_x,max is the maximum value of the k-th objective function of individual x. f_x,min is its minimum value.

   Step 2: Calculate the weight of each indicator. The expression is shown in the following equation:

   (13) $\begin{array}{}\left(\begin{array}{c}{M}_z=-\frac{1}{\ln \left(t\right)}\sum _{x=1}^t\left(\frac{{\mu }_{x,z}}{\sum _{x=1}^t{\mu }_{x,z}}\cdot \ln \frac{{\mu }_{x,z}}{\sum _{x=1}^t{\mu }_{x,z}}\right),\\ {\tau }_z=\frac{1-M_z}{t-\sum _{x=1}^t{M}_z},\end{array}\right)\end{array}$

   where M_z is the information entropy of the k-th attribute. τ_z is the weight of the k-th attribute.

   Step 3: Calculate the scale and fit degree, as shown in the following equation:

   (14) $\begin{array}{}\left(\begin{array}{c}{S}_{x,+}=\sqrt{\sum _{z=1}^w{\left({\tau }_z{\mu }_{x,z}-{\tau }_z{\mu }_{x,z+}\right)}^2},\\ S_{x,-}=\sqrt{\sum _{z=1}^w{\left({\tau }_z{\mu }_{x,z}-{\tau }_z{\mu }_{x,z-}\right)}^2},\\ C_x=\frac{S_{x,-}}{S_{x,-}+S_{x,+}},\end{array}\right)\end{array}$

   where S_x,+ and S_x,− are the positive ideal distance and negative ideal distance of individual x, respectively. μ_x,z+ is the maximum value of all objective functions of individual x after normalization. μ_x,z− is its minimum value. C_x is the fit degree of individual x, and the larger the value, the higher the comprehensive satisfaction of the solution.

   Step 4: Choose the solution with the highest C_x value as the compromise solution.

In the multiobjective optimal configuration model proposed in this algorithm, the variables to be optimized are the installation position and capacity of the ESS. The improved WOA simulates the spatial positions of whale groups. The initialization of the whale population’s locations corresponds to the initialization of numerous ES configurations. The whale population continuously updates its positions through enveloping search, bubble search, random variation, and individual communication until it discovers the optimal solution. This process is equivalent to comparing and optimizing multiple objective functions associated with various configuration schemes. The optimal positions of the whales represent the optimal allocation of ES. The solution steps using the enhanced WOA are outlined as follows.

1. (1) Input the power supply, load timing output parameters, and algorithm-related parameters.

2. (2) Initialize the whale population using cubic chaotic mapping. Then, calculate the objective function values for each whale in the population.

3. (3) Perform rapid nondominant sorting of the whale individuals in the algorithm.

4. (4) Utilize the enhanced WOA to update and optimize the positions of each individual whale.

5. (5) Check whether the optimization process has reached the maximum number of iterations. If not, return to step (2). If the maximum iterations are reached, output the set of frontier solutions.

6. (6) Select the compromise solution based on the normalized processing values and output the final optimal configuration scheme.

The algorithm flowchart of this paper is shown in Figure 4.

Figure 4 The algorithm flowchart of this paper.

The time complexity of the algorithm can be summarized as follows. The input parameter step (step 1) has a constant time complexity O(1). The population initialization and objective function evaluation step (step 2) have a time complexity of O(N + K ∗ N), where N is the population size and K is the number of objective functions. The fast nondominant sorting step (step 3) has a time complexity of O(N²), where N is the population size. The update and optimization of individual positions step (step 4) have a time complexity of O(N). The termination condition check step (step 5) has a constant time complexity O(1). Lastly, the compromise solution selection step (step 6) has a time complexity of O(K ∗ N), where K is the number of objective functions and N is the population size. Overall, the entire algorithm’s time complexity is O(N² + K ∗ N).

3.2. Construction of the ES Model

3.2.1. Constraints

(1) Power-Balance Constraints.

(15) $\begin{array}{}{U}_{\text{sol},n}+U_{s,n}=U_{v,n}+U_{\text{grid},n},\end{array}$

where U_sol,n is the output of the IES at time n. U_s,n is the output of the ESS in time cycle n. This paper assumes that U_s,n remains constant over time n. U_v,n indicates the total power required by the charging device during time n. U_grid,n is the power traded between the IES and the grid during n cycle.

(2) ES Constraints.

1. (a) The maximum and minimum limits of ES output are shown in the following equation:

   (16) $\begin{array}{}-U_{s,\max }⩽U_{s,n}⩽U_{s,\max },\end{array}$

   where U_s,max is the maximum output of the ESS.

2. (b) The following equation illustrates the SOC constraints for energy storage.

   (17) $\begin{array}{}{K}_{\text{soc}}^{\min ,n}⩽K_{\text{soc}}^n⩽K_{\text{soc}}^{\max ,n},\end{array}$

   where Kⁿ_soc is the SOC of the ESS during cycle n. K^min,n_soc and K^max,n_soc are the minimum and maximum SOC of the ESS in time cycle n, respectively.

   The SOC conversion of ES battery is given by the following equation:

   (18) $\begin{array}{}{K}_{\text{soc}}^{n+1}=\left(\begin{array}{cc}{K}_{\text{soc}}^n-{\eta }_n^c\frac{U_{s,n}\Delta n}{F_{ini}},& U_{s,n}⩽0,\\ K_{\text{soc}}^n-\frac{1}{{\eta }_n^{\text{disc}}}\frac{U_{s,n}\Delta n}{F_{ini}},& U_{s,n}>0,\end{array}\right)\end{array}$

   where η^c_n is the charging efficiency of the ES battery in time cycle n. η^disc_n is the discharge efficiency of the ES battery during this cycle [Reference Liu, Wang, Nie, Wang and Chen21].

3. (c) SOC constraints in the final cycle of ESS.

   To guarantee the proper functioning of the ESS in the subsequent control cycle, the SOC after the ES action in the preceding cycle needs to be adjusted to a predefined value. The calculation equation for this adjustment is provided in the following equation:

   (19) $\begin{array}{}{K}_{\text{soc}}^{N_o+1}=K_{\text{soc}}^1.\end{array}$

   N_o is the cycle number. In this article, a cycle is 1day and there are 96 cycles.

4. (d) SOC penalty at the conclusion of the ES cycle.

According to the constraint condition (19), the SOC penalty in the final cycle of ES is shown in the following equation:

(20) $\begin{array}{}{r}_{n,4}=\left(\begin{array}{cc}\left(K_{\text{soc}}^{N_0+1}-K_{\text{soc}}^1\right),& n=N_0,\\ 0,& n<N_0,\end{array}\right)\end{array}$

where r_n,4 signifies the SOC without penalty in the final cycle of ES.

In summary, the reward and penalty function is articulated as in the following equation:

(21) $\begin{array}{}{r}_n={\sigma }_1{r}_{n,1}-{\sigma }_2{r}_{n,2}-{\sigma }_3{r}_{n,3}-{\sigma }_4{r}_{n,4},\end{array}$

where r_n is the immediate return at time cycle n. σ₁, σ₂, σ₃, and σ₄ are the weight coefficients of each component of the reward and penalty, and all are positive. They are in the range [0, 20] and are chosen as σ₁ = 10, σ₂ = 13.5, σ₃ = 15, σ₄ = 18.5, respectively.

3.2.2. Operation Efficiency Model of ES

The efficiency of battery charging and discharging can be ascertained by employing the steady-state circuit equivalent model of the battery [Reference Zhang, Hou and Zhang22], as shown in Figure 5.

Figure 5 Steady-state circuit equivalent model.

V_oc is the open-circuit voltage in the figure. I_b indicates the battery current. U_out indicates the output power. R₁ is the series resistance. R₂ is a short-time response resistance associated with transient response. R₃ represents a short-time response resistor linked to transient response, as well as a long-time response resistor. All the aforementioned variables exhibit a nonlinear relationship with the battery’s SOC, as demonstrated in the following equation:

(22) $\begin{array}{}\left(\begin{array}{c}{V}_{oc}=a_0{f}^{-a_1{k}_{\text{soc}}}+a_2+a_3{K}_{\text{soc}}-a_4{K}_{\text{soc}}^2+a_5{K}_{\text{soc}}^3,\\ R_1=h_0{f}^{-h_1{K}_{\text{soc}}}+h_2+h_3{K}_{\text{soc}}-h_4{K}_{\text{soc}}^2+h_5{K}_{\text{soc}}^3,\\ R_2=c_0{f}^{-m_1{K}_{\text{soc}}}+m_2,\\ R_3=d_0{f}^{-n_1{K}_{\text{soc}}}+n_2,\\ R_t=R_1+R_2+R_3,\end{array}\right)\end{array}$

where R_t is the terminal resistance. a₀-a₅, h₀-h₅, m₀-m₂, and n₀-n₂ are all coefficients.

By solving equation (23) for a given SOC and U_out, the current through an individual cell can be determined.

(23) $\begin{array}{}{I}_b^2{R}_t-I_b{V}_{oc}+U_{\text{out}}=0.\end{array}$

The following equation conveys the charging and discharging efficacy of an individual battery.

(24) $\begin{array}{}\left(\begin{array}{c}{\eta }^c=\frac{V_{oc}}{V_{oc}-I_b{R}_t},\\ {\eta }^{\text{disc}}=\frac{V_{oc}-I_b{R}_t}{V_{oc}},\end{array}\right)\end{array}$

where η^c and η^disc represent the charging and discharging efficiency of the battery, respectively.

The energy storage system (ESS) utilized in this study comprises numerous batteries organized in both series and parallel configurations. The following equations offer the formulas for calculating the charging and discharging efficiencies of the equivalent current, respectively.

(25) $\begin{array}{}{I}_e^2{R}_t-I{V}_{oc}+U_{es}\frac{1}{T_{\text{series}}{T}_{\text{para}}}=0,\end{array}$

(26) $\begin{array}{}\left(\begin{array}{c}{\eta }^c=\frac{V_{oc}}{V_{oc}-I_e{R}_t},\\ {\eta }^{\text{disc}}=\frac{V_{oc}-I_e{R}_t}{V_{oc}},\end{array}\right)\end{array}$

where I_e represents the electric current of the ES battery. I denotes the total current considering multiple batteries arranged in both series and parallel configurations. T_series i indicates the count of battery cells arranged in series, while T_para represents the number of battery units linked in parallel. U_es is the ES battery power.

The operational efficiency was computed using equation (26). The correlation between operational efficiency and output power and SOC was modeled using a second-order polynomial, as demonstrated in the following equation:

(27) $\begin{array}{}\left(\begin{array}{c}{\eta }^c=f_0+f_1{K}_{\text{soc}}+f_2{K}_{\text{soc}}^2+f_3{U}_{es}+f_4{U}_{es}^2+f_5{U}_{es}{K}_{\text{soc}},\\ {\eta }^{1/\text{disc}}=b_0+b_1{K}_{\text{soc}}+b_2{K}_{\text{soc}}^2+b_3{U}_{es}+b_4{U}_{es}^2+b_5{U}_{es}{K}_{\text{soc}}.\end{array}\right)\end{array}$

In the equation, f₀-f₅ and b₀-b₅ denote the coefficients.

3.2.3. ES Attenuation Model

Battery capacity degradation is influenced by ambient temperature, depth of discharge, SOC, and battery runtime. It is a nonlinear process and not a simple sum of individual cycle degradation. This paper employs a semiempirical battery capacity degradation model. Over T loops, the battery-capacity degradation model is defined in the following equation:

(28) $\begin{array}{}{F}_{\text{loss}}=\left(\begin{array}{cc}\left(1-{\mu }_{\text{sei}}{f}^{-T{\gamma }_{\text{sei}}{l}_d}-\left(1-{\mu }_{\text{sei}}\right)f^{-T{l}_d}\right)F_{\text{ini}},& F_{\text{loss}}^{\prime }=0,\\ \left(1-\left(1-F_{\text{loss}}^{\prime }\right)f^{-T{l}_d}\right)F_{\text{ini}},& F_{\text{loss}}^{\prime }>0,\end{array}\right)\end{array}$

where F′_loss and F_loss denote the decay capacities before and after T loops, respectively. IfF′_loss = 0, it indicates that the battery has not undergone charging or discharging processes. μ_sei and γ_sei are the influencing factors for solid electrolyte interface (SEI) membrane formation during the use of new batteries. l_d is the capacity attenuation function of a single cycle.

The residual usable capacity F_r of the ES battery is shown in the following equation:

(29) $\begin{array}{}{F}_r=F_{\text{ini}}-F_{\text{loss}},\end{array}$

where F_ini is the maximum capacity of the ES battery before attenuation. F_loss is the discharge capacity of the ES battery.

Equations (28) and (29) represent the approach for calculating battery capacity degradation. However, obtaining parameters for each cycle can be challenging in practice because battery charge and discharge patterns are often irregular. To address this, the rain-flow counting method is utilized to compute battery loops and their associated parameters [Reference Feng, Zheng and Qiao23]. Figure 6 illustrates the process of calculating ES battery capacity degradation using this method.

Figure 6 Calculation flow of ES battery capacity attenuation.

4. Result Analysis and Discussion

This paper has established a comprehensive ESS model, encompassing the system operation model and the ES decay model. By optimizing these models, the average annual revenue status of the ESS is determined. The effectiveness of the algorithms and models proposed in this paper is demonstrated through comparisons and analyses with other approaches. The abstract underscores the importance of optimizing ESS capacity allocation to address the volatility and uncontrollability of renewable energy generation for energy stability. Given regional constraints, seasonal variations, and mid-day climate impacting photovoltaic energy systems, this paper indirectly derives the optimal capacity and power solution through the average annual return profile of the ESS.

740 arrays of photovoltaic cells are utilized in this investigation, collectively possessing an ES capacity of 14.10MW. When the capacity of 1949 storage batteries declines to 80% of their designated capacity, they are no longer deemed functional. The ESS operates at a charging and discharging efficiency rate of 90%. The overall system cost is approximately 32.8 million yuan, and the time-of-use electricity prices are detailed in Table 1.

Table 1 The time-of-use electricity prices used in the simulation case.

In this paper, three comparative algorithms are selected, and their specific details are presented in Table 2. During testing, each algorithm’s specific test conditions and parameters are described in Figure 7. For their proposed algorithm, the authors likely conducted experiments or simulations to determine the optimal parameter values. The sensitivity of these parameters to the performance of the proposed algorithm is of utmost importance in ensuring its effectiveness. Initially, the fitness function weights are set by assigning suitable weights to each objective. Altering these weights enables the algorithm to discover optimal solutions. Subsequently, termination criteria are defined. These criteria dictate when the algorithm should cease iterations. By setting a maximum number of iterations or a convergence threshold, the algorithm’s convergence is enhanced, and excessive computations are avoided. Evaluating the algorithm’s performance involves individually or collectively varying these parameters to identify their optimal values. This process aids in understanding how parameter changes impact the algorithm’s behavior, facilitating fine tuning for improved performance.

Table 2 Description of the proposed algorithm and the referenced algorithm.

Figure 7 Comparison of the four methods.

When the algorithm in this paper uses the actual data for training, only the current cycle data are obtained in the test stage and the loss of the ESS is considered. Literature [Reference Guan, Lv, Li, Su and Dong24] addresses the ESS charging and discharging approach under theoretical conditions, assuming all the factual information. In contrast, literature [Reference Walker and Kwon18] utilizes “actual data +17% deviation data within a normal distribution” during the training phase, introducing variations between the data and the factual information. In order to verify the robustness of the calculation, the method considers many uncertain factors such as comprehensive energy, electrical equipment charging, and electricity price. Literature [Reference Yin and Qin25] also assumes that the charging and discharging strategy of the ESS is solved when the actual data is known, but the loss of ES is not considered. Since the proposed algorithm converts all costs of ES, the average annual income is taken as the standard to judge the optimization effect.

Table 3 shows the benefits of ESS calculated by different algorithms. Figure 8 shows their corresponding ES capacity attenuation curves.

Table 3 Revenues of the photovoltaic-storage charging station.

Figure 8 ES capacity attenuation curve.

The results presented in Table 3 indicate that the ESS operated by the algorithm in this paper has a duration of 3,124 days, which is the longest among the compared algorithms. Moreover, it yields the highest life cycle benefit for the ESS. The average annual income from this algorithm ranks second, with only literature [Reference Guan, Lv, Li, Su and Dong24] outperforming it. However, it is worth noting that literature [Reference Guan, Lv, Li, Su and Dong24] achieves the best optimization results because it considers an ideal scenario where it already knows the actual data of all cycles in advance when optimizing. While this approach results in the highest average annual return, it necessitates having access to actual data within the optimization cycle, which might not be practical in real-world applications. The algorithm in this paper operates under more realistic conditions, considering only current data during testing, which may explain the difference in performance.

In comparison to literature [Reference Guan, Lv, Li, Su and Dong24], literature [Reference Walker and Kwon18] uses more data than the actual data. Specifically, it adds extra 17% deviation data based on the actual data, with these data following a normal distribution. This approach results in a notable reduction in the lifespan of the ESS, as well as a decrease in the system’s life cycle income and annual average income. However, when compared to the algorithm presented in this paper, the average annual return from literature [Reference Walker and Kwon18] is only slightly lower. This comparison highlights that the proposed algorithm can effectively address the charging and discharging decision problem of an ESS in an uncertain environment, offering competitive results. It suggests that the algorithm’s ability to optimize ESS performance even in the face of uncertainties makes it a valuable tool for practical applications.

The analysis of Figure 8 reveals that when not considering energy action losses, the ESS optimized by literature [Reference Walker and Kwon18] has the shortest lifespan at 1314 days. This is because literature [Reference Walker and Kwon18] does not account for energy action losses. In comparison, literature [Reference Guan, Lv, Li, Su and Dong24] uses data with an additional 17% normal distribution bias, leading to a reduction of 389 days in the lifespan of the ESS when compared to literature [Reference Feng, Zheng and Qiao23]. In contrast to literature [Reference Feng, Zheng and Qiao23], which optimizes the entire cycle of data, the algorithm presented in this paper utilizes data from the current time cycle, resulting in the longest lifespan for the system. This surpasses the performance of literature [Reference Feng, Zheng and Qiao23] by 850 days. Thus, it can be inferred that energy capacity losses, which lead to frequent charging and discharging of the ESS, accelerate the decay of its lifespan and ultimately reduce its long-term benefits. This demonstrates the importance of considering energy action losses in optimizing ESSs, as our algorithm achieves the longest system lifespan, highlighting its effectiveness.

Absolutely, we have made a succinct and accurate observation. ES capacity loss indeed leads to frequent charge and discharge loops of the ESS. This accelerated cycling of ES components can significantly reduce the overall lifespan of the system and, consequently, diminish the long-term benefits it provides. This insight underscores the importance of considering and mitigating ES capacity decay in the design and optimization of ESSs to ensure their sustained effectiveness.

The test result of energy storage SOC on a certain day is shown in Figure 9. The proposed algorithm can discharge the ESS when the electricity price is high and charge it when the electricity price is low. This shows that the improved algorithm in this paper can learn the scheduling strategy that makes the return higher. Compared to literature [Reference Walker and Kwon18], literatures [Reference Feng, Zheng and Qiao23, Reference Guan, Lv, Li, Su and Dong24] reduce the number of charging and discharging loops, resulting in lower power for these operations. In comparison to literatures [Reference Walker and Kwon18, Reference Feng, Zheng and Qiao23], our algorithm achieves the least number of charging and discharging loops for the ESS. This is because in our algorithm, we assign a higher weight coefficient to the cost of ES, which extends the lifespan of the system. During the decision-making stage, our algorithm sacrifices short-term benefits by reducing immediate gains in order to accumulate long-term benefits.

Figure 9 SOC of energy storage.

In this case, σ₂ relating to the ESS’s capacity reduction penalty is utilized to delve deeper into the impact of the weight coefficient on the outcomes of the reinforcement learning algorithm. The statistical results of sensitivity analysis with the weight coefficient σ₂ are shown in Table 4.

Table 4 Sensitivity review of the weight coefficient σ₂.

The analysis in Table 4 shows that with the increase of σ₂, the life of the ESS extends. At the same time, the life cycle returns also increase. When σ₂ is less than 5, the average annual revenue of the ESS increases. This indicates that the choice of actions is more likely to avoid the end of the ES penalty. As σ₂ increases, the way in which ES can be acted and profited from becomes more rational. However, when σ₂ is greater than 5, the average annual revenue of the ESS begins to decrease and the decrease amplitude gradually increases. This indicates that the ESS prioritizes the effects of capacity reduction from ES actions when forming judgments. As σ₂ increases, the ESS reduces the number of ES operations to delay capacity decay. As a result, prolonged operational longevity and increased overall life cycle advantages can be achieved.

The parameters utilized in the algorithm simulation within this study, like σ₂, were determined after numerous trial simulations to derive a set of comparatively superior experimental parameters. Nevertheless, this set of trial parameters might not represent the most optimal, highlighting the challenge of parameter adjustment in the whale algorithm. From Table 5, it can be observed that our algorithm requires the longest time during training, but it operates faster in the testing phase after the model is trained, with an average optimization time of 0.0031 s per cycle. The time taken to test our algorithm is 68.36 seconds, and the optimization duration per cycle is independent of the issue’s intricacy. The other three algorithms do not require pretraining, so their training time is 0 s. Since they assume a constant operating efficiency for the ESS, they all require longer testing time compared to our algorithm. Literature [Reference Feng, Zheng and Qiao23] has the longest testing time, while literature [Reference Guan, Lv, Li, Su and Dong24] has the longest average optimization time per time cycle. Due to the fact that literature [Reference Walker and Kwon18] does not consider the degradation of the ESS, its average optimization time per time cycle is less than that of literature [Reference Feng, Zheng and Qiao23] and literature [Reference Guan, Lv, Li, Su and Dong24]. However, it still requires more time compared to our algorithm. This indicates that our improved algorithm maintains a higher computational efficiency in complex problems.

Table 5 Comparative analysis of computation time among the four algorithms.

The iterative convergence curves of the proposed, literature [Reference Gharehchopogh, Nadimi-Shahraki, Barshandeh, Abdollahzadeh and Zamani26], and literature [Reference Zhang, Khayatnezhad and Ghadimi27] are shown in Figure 10. It can be seen from Figure 9 that under the same population size and the same number of iterations, the objective functions of the three algorithms decrease with the increase of iterations. However, the convergence speed of the proposed algorithm is fast, which reduces the risk of premature convergence of particles. This shows that the proposed algorithm has strong global and local search ability and faster convergence speed.

Figure 10 Fitness convergence curve.

In order to provide empirical evidence supporting the significant improvement advantages of the improved algorithm proposed in this paper over other existing algorithms, a Wilcoxon rank-sum test [Reference Gao, Zhuang, Chen and Wang28] was conducted independently to compare the results of each experiment at a significance level of 5%. The comparison results are presented in Table 6, where “p” represents the test result and “h” indicates the significance judgment result. If “p” is less than 0.05, “h” is denoted as 1, which means the significance of this algorithm is stronger than other algorithms. Conversely, if “p” is greater than 0.05, “h” is denoted as 0, indicating that the significance of this algorithm is weaker than other algorithms. In cases where “p” is displayed as N/A, it indicates that a significance test could not be performed and the significance of the proposed algorithm may be similar to that of other algorithms. Based on the statistical results presented in Table 6, it is evident that the majority of the p values are below the significance level of 0.05. This compellingly indicates that there exist significant differences between the proposed algorithm and other algorithms under comparison. Notably, the performance of the proposed algorithm consistently outperforms the comparison algorithms across various functions. These findings strongly support the effectiveness and superiority of the proposed approach, lending further credibility to its potential for improved performance and applicability.

Table 6 Results of the Wilcoxon rank-sum test.