## 1 Introduction  
As one of the critical operation tasks of power systems, economic dispatch (ED) is usually employed to minimize the total operating cost via optimally calculating the active power outputs of all the generators to balance the active power demand under various operating constraints [1]. However, the conventional ED only considers the electricity energy dispatch of the thermal units, which cannot satisfy the operation requirement of an integrated energy system (IES) with various energy types and sources. In light of this issue, the ED problem of an IES has attracted extensive investigation. By incorporating customer aggregators' flexible energy demand into the centralized energy dispatch model, a two-level optimization problem [2] was formed in the heat and power IES. In [3], a regional integrated energy system including wind turbines (WT), photovoltaics (PV), gas turbines and battery energy storage was introduced to minimize the operating costs of the system and determine the optimal coordination between the various energy sources. In order to fully exploit the economic and environmental advantages of the system, considering the difficulties of information collection from subareas, a novel decentralized optimal multi-energy flow (OMEF) for large-scale IESs in a carbon trading market was proposed in [4]. The improved differential evolution algorithm [5] was designed to obtain the minimum operation cost of the IES, while considering battery lifetime loss. In [6], a novel multi-agent bargaining learning (MABL) model was employed for the distributed energy hub economic dispatch (EHED) of multiple energy carrier systems (MECS). Among these studies, the combined heat and power ED is the most frequent research interest. Similarly, it usually aims to minimize the total operating cost while satisfying the electricity and heat energy balance constraints and other operating constraints. In order to take the wind and PV energy into account, their stochastic and intermittent power outputs were introduced in the combined heat and power ED via construction of the probability density functions of the wind speed and solar irradiance [7]. Particularly, this new ED is also called the optimal energy dispatch (OED) of the combined heat and power-thermal-wind-photovoltaic systems in this paper. In essence, it is a highly complex optimization with highly non-smooth, non-linear, and non-convex features, due to the uncertainty of wind and PV energy outputs, the valve-point loading effect of thermal generators, and the non-convex feasible operating region constraint of combined heat and power (CHP) units. Generally speaking, the presented OED can be addressed by two types of methods. The first type is a classical optimization algorithm, e.g., the quadratic programming method [8], interior point method [9], and the Newton method [10]. However, most of them are gradient-based optimization methods, which are highly dependent on a specific mathematical model. If the optimization is non-linear with discontinuous functions and multiple local optimums, this will easily trap into a low-quality local optimum, especially for the OED. On the other hand, it can be also solved by metaheuristic optimization algorithms, such as genetic algorithms (GA) [11], particle swarm optimization (PSO) [12], differential evolution (DE) algorithms [13] and the grey wolf optimizer (GWO) [14], which are highly independent from the specific mathematical models, and are much easier to apply for OED compared to the first type of method. Motivated by this advantage, this paper proposes a novel algorithm named multi-searcher optimization (MSO) for the OED of combined heat and power-thermal-wind- photovoltaic systems, which has the following novelties:  
* Traditional ED only considers the electricity energy dispatch of power systems, but does not consider the heat energy dispatch of thermal systems.