## 2 Matrix-Based Model Formulation for MES
### District Gas Network
#### 2.3.1 Matrix-Based Universal DGN Model  
Prior to the derivation of the mathematical model, several assumptions for the steady-state gas flow calculation are made as the foundational basis of this work [30]: (1) no storage effect, leakage, or any form of energy loss exists during gas transmission; (2) the pipeline is fully rough, and the friction factor remains constant throughout the pipeline; (3) a perfect mixture is assumed after the alternative gas source is injected; (4) gas volume flow rates mentioned are measured under standard temperature and pressure (STP) conditions. According to [PERSON]'s equation, the gas flow within the pipeline can be generalized in matrix form as Equation (18), whose detailed parameters are defined by Equation (19):  
\[f_{\text{P}}=S_{\text{g}}\left(\mathbf{K_{\text{P}}}^{-1}\cdot|\mathbf{\Pi_{s }}-\mathbf{\Pi_{e}}|\right)^{\frac{1}{\lambda}} \tag{18}\]  
\[\left\{ \begin{array}{l}S_{\text{g}}=\text{diag}\Big{[}(\mathbf{\Pi_{s}}- \mathbf{\Pi_{e}})\cdot|\mathbf{\Pi_{s}}-\mathbf{\Pi_{e}}|^{-1}\Big{]}\\ K_{\text{P}}=\text{diag}(\mathbf{M_{P}}\cdot\mathbf{Z_{P}})\overline{-A_{\text{g }\text{g}\text{d}}}\rho^{\mu}\\ M_{\text{P}}=C_{\text{P}}T_{\text{g}}\boldsymbol{\varphi_{\text{P}}}\cdot \mathbf{L_{P}}\cdot\mathbf{D_{P}}^{\kappa}\\ A_{\text{gd}}=A_{\text{g}}S_{\text{g}}\end{array} \right. \tag{19}\]where \(f_{\rm P}\) is the gas pipeline flow vector; \(S_{\rm g}\) is the directional factor vector of gas pipelines; \(K_{\rm P}\) is the pipeline coefficient vector; \(\Pi_{\rm g}/\Pi_{\rm e}\) is the squared pressure vector of the start/end node; \(M_{\rm P}\) represents the constant part in \(K_{\rm P}\); \(Z_{\rm P}\) is the pipeline compressibility factor vector; \(\rho\) is the specific gravity vector; \(\varphi_{\rm P}/L_{\rm P}/D_{\rm P}\) is the friction factor/length/diameter vector of the pipeline; \(C_{\rm P}/T_{\rm g}/\chi/\kappa/\mu\) is the constant parameter/gas temperature/flow exponent/diameter exponent/specific gravity exponent; \(A_{\rm g}\) is the node-pipeline incidence matrix of the DGN; and \(A_{\rm gd}\) is the directional node-pipeline incidence matrix. Compressors are installed to regulate the pressure level of the DGN. The adiabatic horsepower \(H_{\rm C}\) of the compressor can be calculated by  
\[\left\{ \begin{array}{l}\mbox{$H_{\rm C}$}=\mbox{$B_{\rm C}$}\mbox{$\cdot$} \mbox{$f_{\rm C}$}\mbox{$\cdot$}\left[\left(\mbox{$\pi_{\rm o}$}\mbox{$\cdot$} \mbox{$\pi_{\rm i}$}^{-1}\right)^{\theta_{\rm C}}-1\right]\\ \mbox{$B_{\rm C}$}=\mbox{$C_{\rm C}$}\mbox{$\cdot$}\mbox{$P_{\rm st}$}\mbox{$T _{\rm st}$}^{-1}\mbox{$T_{\rm g,in}$}\mbox{$Z_{\rm C,in}$}\mbox{$\cdot$}\mbox{$ \theta_{\rm C}$}\mbox{$\cdot$}\\ \mbox{$\theta_{\rm C}$}=\left(\mbox{$\omega$}-1\right)\mbox{$\cdot$}\mbox{$ \omega$}^{-1}\end{array} \right. \tag{20}\]  
where \(\mbox{$B_{\rm C}$}/\mbox{$f_{\rm C}$}/\mbox{$\theta_{\rm C}$}\) is the power coefficient/compressor gas flow/compression exponent vector of the compressor; \(\mbox{$\pi_{\rm i}$}/\mbox{$\pi_{\rm o}$}\) is the inlet/outlet pressure vector of the compressor; \(C_{\rm C}\) is the compressor characteristic constant; \(P_{\rm st}/T_{\rm st}\) is the standard temperature/pressure; \(T_{\rm g,in}\) is the gas temperature at the compressor inlet; \(\mbox{$Z_{\rm C,in}$}\) is the gas compressibility factor vector at the compressor inlet; and \(\omega\) is the gas polytropic coefficient vector.