Multi-dimensional holomorphic embedding method (MDHEM) is a new method for solving the energy flow of integrated energy systems (IES). This method boasts two key advantages: freedom from initial value dependence and non-singularity of the Jacobian matrix. However

it faces notable challenges in probabilistic energy flow (PEF) calculations

particularly the requirement for high-order power series expansions and the resulting heavy computational burden. In order to address the above issues

this paper proposes a multi-germ solution and multi-dimensional holomorphic embedding method (MG-MDHEM) to calculate the PEF of IES via the semi-analytical solution characteristics of HEM. Firstly

the MDHEM-based energy flow model of IES is formulated. Then

the fluctuation range of the multidimensional stochastic inputs in IES is divided into multiple subintervals

and physical solutions of MDHEM in each subinterval are set up. In each subinterval

the low-order semi-analytical solution expressions for the variables can be recursively obtained. Further

for the stochastic samples obtained by Monte Carlo simulation

the samples are preprocessed and substituted into the low-order semi-analytical solutions of the variables in the subintervals

and the Pade approximation is used to solve the low-order semi-analytical formulations using parallel multi-dimensional difference methods to achieve the energy flow result of IES in each sample

and the probability distribution of PEF for the IES is hereby obtained using the energy flow calculation results of Monte Carlo simulation. Finally

the accuracy and effectiveness of the proposed method are proved by the energy flow calculation results of the E14-G6 and E118-G20 test systems.