## 3 Modeling Method Used in This Work
The mathematical problem statement is as follows:  
Find the minimum of the objective function  
\[Z=\sum_{i\in I}(F_{i}+\mu_{i})k_{i}\overline{N}_{i}+\sum_{j\in J_{2}}(F_{j}+ \mu_{j})k_{j}\overline{Q}_{j}+\frac{H}{T}\sum_{t=1}^{T}\sum_{j\in I_{1}}p_{j} \sum_{i\in I}N_{it}\alpha_{ij}(s_{it})\rightarrow\min \tag{4}\]  
under the primary energy consumption constraints;  
\[\sum_{i\in I}N_{it}\alpha_{ij}(s_{it})\leq P_{jt},\quad j\in J_{1} \tag{5}\]the final energy consumption constraints;  
\[\sum_{i\in I}N_{it}\beta_{ij}(s_{it})\geq L_{it},\quad j\in J_{3} \tag{6}\]  
the secondary energy balances;  
\[\sum_{k\in I}N_{kt}\beta_{kj}(s_{kt})\Delta t-\sum_{k\in I}N_{kt}\alpha_{kj}(s_{ kt})\Delta t-X_{it}+Y_{it}-\Delta Q_{it}=0,\quad j\in J_{2} \tag{7}\]  
\[Q_{it}=Q_{j,t-1}+\eta_{j}X_{it}-Y_{jt}\,\ j\in J_{2} \tag{8}\]  
variables constraints  
\[0\leq N_{it}\leq\overline{N}_{i},\quad i\in I \tag{9}\]  
\[X_{jt}\geq 0,\ \ Y_{jt}\geq 0,\ \ \Delta Q_{jt}\geq 0,\ \ \ \ j\in J_{2} \tag{10}\]  
\[0\leq Q_{jt}\leq\overline{Q}_{j},\quad j\in J_{2} \tag{11}\]  
subject to the initial condition  
\[Q_{j0}=\frac{1}{T}\sum_{t=1}^{T}Q_{jt},\ \ j\in J_{2} \tag{12}\]  
for \(t=1,2,\,\ldots,\,T\) and \(s_{it}\in S_{i}\). Here  
\[F_{i}=\frac{\ln(1+d)}{1-\left(1+d\right)^{-\Delta T_{i}}} \tag{13}\]  
is the capital recovery factor at the annual discount rate \(d\); \(H\) is the duration of the year in selected time units (for example, \(H=8760\) h/year); \(p_{j}\) and \(P_{jt}\) are the price and the maximum available volume of the primary energy resource, respectively; \(L_{jt}\) is the final energy demand; \(X_{jt}\) and \(Y_{jt}\) are the secondary energy flows into and out of the energy storages, respectively; \(\Delta Q_{jt}\) is the secondary energy losses. The decision variables are the installed \(\overline{N}_{i}\) and current \(N_{it}\) capacities of conversion energy installations, secondary energy reserves \(Q_{jt}\) in the energy storages and their capacity \(\overline{Q}_{j}\), energy flows \(X_{jt}\) and \(Y_{jt}\), losses \(\Delta Q_{jt}\), and controlled operating modes \(s_{it}\). Stochastic modes of energy sources using renewable energy (for example, solar) cannot be optimized because they are determined by external environmental conditions and have preset values for the model. The objective function (4) is the total (per year) discounted costs for the creation and operation of the system. For a given amount of energy supplied to final energy consumers, the minimum of \(Z\) corresponds to the maximum of the net present value. The term \(\Delta Q_{jt}\) in Equation (7) accounts for possible (at certain moments) \"excess\" production of secondary energy. For example, it is sometimes more beneficial to employ a dump load to consume the excess energy of RES using a free resource (solar or wind) than to increase the battery capacity or reduce the installed capacity of the wind turbine. The initial condition (12) imitates a random selection of the beginning of the simulation period on the time axis (to estimate a random value of the initially available energy resource, we use its most probable value, i.e., its arithmetic mean).