The Optimal Selection of Renewable Energy Systems Based on MILP for Two Zones in Mexico

: This paper presents a series of enhancements to a previously proposed mixed-integer linear programming (MILP) model for investment decisions and operational planning in distributed generation (DG) systems. The main contribution of this study consists of integrating a wind generation system and multiple loads at different buses in a network. The model considers dynamic weather data, energy prices, costs related to photovoltaic and wind systems, storage systems, operational and maintenance costs, and other pertinent factors, such as efficiencies, geographical locations, resource availability, and different load profiles. The simulation results obtained through implementation in Julia’s programming language illustrate that the MILP formulation maximizes the net present value, and four configurations for hybrid power generation systems in Mexico are analyzed. The objective is to enable profitability assessment for investments in large-capacity DG systems in two strategic zones of Mexico. The results show that the configurations in the NE zone, especially in Tamaulipas, are the most cost-effective. Case 1 stands out for its highest net present value and shortest payback time, while Case 2 offers the highest energy savings. In addition, Cases 3 and 4, which incorporate storage systems, exhibit the longest payback periods and the lowest savings, indicating less favorable economic performance compared with Cases 1 and 2. Moreover, the sales of two case studies, one without a storage system and the other with a storage system, are shown. The model also incorporates instruments for buying or selling energy in the wholesale electricity market, including variables that depict the injected energy into the electrical grid. This comprehensive approach provides a detailed overview of optimal energy management.


Introduction
Several studies have focused on determining the optimal sizing of distributed generation (DG) systems to meet demand requirements.Some of these studies present a techno-socio-economic analysis to incorporate the DG systems in rural zones [1].In [2], an economic evaluation of batteries as a load-shifting technique in a photovoltaic generation system for a building is presented.The study considers photovoltaic (PV) system generation, electricity consumption, and tariffs.The research concludes that the Carnot battery studied is not cost-effective for their specific case while highlighting that batteries are becoming more financially competitive in current scenarios.The authors in [3] present the appropriate combination of a PV system with a battery storage system (Bat) as an option that maximizes savings and provides flexibility and supply availability in the Italian energy market.On the other hand, the study in [4] proposes an optimal model for inserting energy storage systems in electric transmission networks using mixed-integer linear programming.The objective is to reduce operating costs and unsupplied energy.It highlights reducing dependence on fossil fuels, savings in generation costs, and stability in energy supply.Therefore, incorporating energy storage systems into the power grid offers several benefits.These systems optimize fuel requirements, reduce environmental impact, and are highly responsive to load changes, which dampens demand fluctuations.In addition, they allow for greater integration of distributed and renewable resources and provide backup during peak hours or supply other needs, as shown in [5].However, adding a Bat increases the investment value and requires a larger PV system capacity.Then, consumers must adjust consumption and battery storage according to market conditions.According to [6], the authors develop a model considering intermittent renewable power plants with low marginal cost based on stochastic and nontradable inputs.Batteries have been demonstrated to smooth prices over time periods; however, they exhibit diminishing marginal benefits.Consequently, the incorporation of Bat is limited even as battery prices decrease.In this scenario, relying only on a hybrid system composed of only one intermittent generation source is not a viable option, even with a Bat.Additionally, when considering the economic risks of an electric market, [7] analyzes the impact of integrating a large-scale Bat in two power grids-France and Germany.The research aims to determine the increase in consumption costs and concludes that large-scale Bat can have significant repercussions on the operation of conventional generation and influence the adoption of renewable sources [8].Thus, large-scale Bat is crucial for transitioning to renewable generation.The optimal combination is established as wind (Wd) and PV generation [9], as Wd energy development is increasing [10].However, the random nature of wind speeds poses challenges in dispatching energy efficiently [11].
Most papers employ traditional budgeting techniques, such as the analysis of cash flows and net present value (N PV ) [12,13].Renewable sources are easily accessible in remote and rural areas, where they can meet energy demands [14].Additionally, hybrid systems not only provide support to remote communities but also present potential solutions for sustainable buildings [15,16].Furthermore, the industrial sector, which consumes the most energy worldwide, progressively incorporates renewable energy sources.Efforts to efficiently manage peak demands led several nations to implement Time-Of-Use (TOU) fees.These prices incentivize manufacturing businesses to shift their power usage from peak to off-peak hours [17].As a result, industrial organizations are facing a growing necessity to adopt energy-efficient operations and maintenance programs.Amid fluctuating energy demands, maintenance planning and production scheduling become intertwined with patterns of electricity use [18].At the same time, the approaches presented in [19] imply that customers could choose self-consumption over grid electricity when PV and wind (Wd) self-consumption systems become economically viable.However, significant obstacles to electrical stability are associated with the mass deployment of decentralized solar PV [20].Therefore, the challenge is to develop cost-effective and systematic solutions to enable the smooth integration of these variable-generating systems.
The authors propose models based on mixed-integer linear programming (MILP) and utilize the net present value (N PV ) as the evaluation metric [21].Other works, such as [22], incorporate the concept of locational marginal price and annual increases in electricity demand to assess the profitability of photovoltaic (PV) and battery (Bat) power systems.Moreover, the potential of offshore wind resources in the Gulf of Mexico and their integration into the power system are highlighted in [23,24].
This research aims to enhance the current literature by introducing additional constraints to consider an offshore wind power system, building upon the model presented in [21,22].
To address the gaps identified in the literature and contribute to the ongoing research, this study aims to achieve answers to the following research questions using MILP:

•
How does the inclusion of locational marginal price (LMP) and annual percentage increases in electricity demand affect the economic profitability of large-capacity photovoltaic (PV) and battery (Bat) systems?• What are the impacts of incorporating offshore wind energy resources in the Gulf of Mexico into the Mexican power system?
• How does adding a second load influence the economic feasibility and performance of hybrid PV, wind (Wd), and battery (Bat) systems within the Mexican electricity market?
Moreover, the economic profitability of incorporating a second load into the model is assessed.The main contributions of this study to the existing literature include:

•
The objective function of the model is adjusted to include the locational marginal price (LMP), and a parameter representing the annual percentage increase in electricity demand is introduced to evaluate the profitability of large-capacity PV and Bat systems.

•
The findings reported in [23], which evaluate the offshore wind resource potential in the Gulf of Mexico, were incorporated as input into the simulation process.

•
The configuration of offshore wind energy systems (Wd) and their integration into the Mexican power system are analyzed, introducing the necessary constraints to consider a Wd system in the model presented in [22].

•
The economic profitability is evaluated through the inclusion of a second load into the modeling.Moreover, unlike the model introduced in [21], which solely considered one load and a residential tariff, our modeling approach can incorporate loads at multiple network nodes, like nodesP and nodesD, within the Mexican electricity market.
This paper is organized as follows: Section 2.1 presents the mathematical model based on MILP.Sections 2.2 and 2.3 describe and prepare the data required for the model.Section 3 analyzes the results, showing the behavior of the hybrid PV, Wd, and Bat system determined by the model.Finally, the conclusions derived from this research are presented in Section 4.

Mathematical Model
A mixed-integer linear programming (MILP) model is proposed for sizing a PV-Bat system in the study published [21,22].However, unlike the references mentioned earlier, the present paper focuses on the MILP formulation of two specific nodes, denoted as EP1 p and EP2 p in period (p : 1, 2, . . ., P).The problem involves decision-making across multiple periods and incorporates a wind (Wd) system with four candidate wind turbine models.To accommodate the inclusion of the Wd system and its interaction with the second load, modifications were made to constraints (4), ( 5), ( 6), ( 8), (12) and (20).The obtained results serve as inputs for another linear programming model aimed at maximizing energy sales in response to variations up in the local marginal price (LMP) at the selected nodes.
Figure 1 presents the configuration of a PV-Wd-Bat system that powers two distinct loads.The PV system generates direct current (DC) power, while the Bat system utilizes DC power with the assistance of a DC-DC converter acting as a voltage regulator before converting it into alternating current (AC) using a DC-AC inverter.Load 1 P Load 1 p is supplied by the PV η PV in = 0.97 and Bat η B in = 0.97 inverters [21,25].The efficiencies of these components remain constant and are independent of external factors.The DC power generated by the PV system can directly charge the batteries through a DC charge controller, eliminating the need for AC conversion.The charge and discharge efficiencies of the batteries are assumed to be 90% each (η C and η D ) [18].On the other hand, the Wd system supplies AC power, and the wind power is rectified before charging the Bat system to ensure a stable power supply from all four sources (Grid, PV, Wd, and Bat).+ - Configuration of a photovoltaic-wind-battery generation system, considering two loads, connected to the grid.

Model Constraints
The model incorporates I discrete PV systems (i : 1, 2, 3, . . ., I) with a corresponding capital expenditure (CAPEX) CE PV and area A i .These PV systems have an efficiency of η PV = 20% [17,26].Additionally, J Bat systems are included in the model, each with its capital expenditure CE B and capacity (j : 1, 2, 3, . . .J), as well as their respective storage capacities.The model also considers K Wd systems (k : 1, 2, 3, . . ., K) with power generation capacity G AG .To represent the selection of PV, Bat, and Wd system candidates, binary variables are utilized for i, j, and k, respectively.
If the binary variable is equal to 1, the candidate system is selected.To limit the number of systems, PV i , Bat j , and Wd k , N PV , N Bat , and N Wd are determined, respectively: The total area where the PV system is intended to be installed, A max can be delimited and is given by If a PV system (PV i ) is selected, it provides power in the form of direct current (DC).In each period p, this power can be utilized by the loads or transmitted to the grid through a DC-AC inverter.Alternatively, in situations where there is excess generation and low demand, the power can be directed toward charging the Bat system.This applies to all 1 ≤ i ≤ I and 1 ≤ p ≤ P.
The decision variables P PG ip , P PL 1 ip , and P PL 2 ip represent the power transmitted from PV i to the grid (G) and to loads P Load 1 p and P Load 2 p , respectively.Similarly, P PB ijp denotes the DC power transmitted from PV i to Bat j .On the other hand, the Wd system (Wd k ) supplies AC energy during each period p through the variables P WG kp , P W L 1 kp , and P W L 2 kp Therefore, the energy supplied by the Wd system is analogous to that of PV and is expressed for all periods p and all Wd systems Wd k .This relationship holds for 1 ≤ k ≤ K and 1 ≤ p ≤ P.
The load P Load 1 p can receive power supply from four distinct sources: PV i , Bat j , Wd k , and the grid G.This relationship is expressed as follows: for all 1 ≤ p ≤ P.
In the same way as the previous equations, an equation for the second load is formulated, with a reliability of 100%.Let P grid 2 p represent the energy injected into the grid to provide energy to Load 2 .
So, the energy supplied by the generation systems PV i and Wd k is then used to provide energy to Load 2 .
Adding Bat storage devices is essential for improving the supply reliability of a system where there are a lot of intermittent sources.Batteries can stabilize prices over time, however this price-smoothing impact results in declining returns on their marginal benefit.Consequently, even with falling battery prices, the expected rate of battery adoption is limited [6].This restriction makes the integration of a Bat system a viable option.For any 1 ≤ j ≤ J and 1 ≤ p ≤ P, the energy flow balance of battery B jp is thus determined by the following constraint at any given period p, where P BG jp represents the energy transmitted from battery j to the grid during period p.
Another variable introduced is the State Of Charge (SOC) of the battery; unlike [21,22], we propose a formula for the SOC that does not require a higher computational expense and is formulated as The SOC must be controlled within specified upper (SOC Max ) and lower (SOC min ) limits during operation.
However, the battery cannot charge and discharge simultaneously.To incorporate this constraint, we introduce the binary variable y B jp , which takes a value of 1 when the battery is charged by the PV, Wd, or grid.The following constraints are imposed for all 1 ≤ j ≤ J and 1 ≤ p ≤ P to ensure the validity of this constraint.
The parameter M represents a large constant, with a value of 10,000 utilized in the application of the Big M method [27].Furthermore, the battery is subject to charge and discharge limits specified by the values (CR j , DR j ) for each period p.This constraint holds for all 1 ≤ j ≤ J and 1 ≤ p ≤ P.
To estimate the optimum quantity of energy, PV, Wd, and Bat systems are unable to export power to the grid while simultaneously supplying power to Load 1 during period p.This condition can be expressed through the following equations.The variable y G takes a value of 1 when energy is imported from the grid and 0 otherwise.These constraints hold for all 1 ≤ p ≤ P. P

Net Present Value Optimization Program
In the context of renewable energy investment projects, net present value (N PV ) evaluates economic profitability through the difference between assets (energy sold to the grid) and liabilities (any cost and purchase of energy from the grid).Therefore the "N PV " is defined as: where the following are used: • Assets-It is a resource with a value that generates a future benefit.

•
Liabilities-Any resource that consumes part of the asset.
Then, with the variables and constraints stipulated above, the calculation of the objective function for the Net Present Value "N PV " over the planning horizon H is given by max The equation takes into account the safety factor F s to impose a limit on the energy supplied to the grid, preventing congestion.It also incorporates the interest rate r, the annual operation and maintenance costs FOM PV i , FOM Wd k , and FOM Bat j relative to the capital costs CE i,j,k , as well as the annual percentage increase in energy ie as a factor of demand growth.These aspects are considered within the MILP model assuming constant inverter efficiency.

Linear Sales Optimization Program
Once the optimal value N PV for the PV-Wd-Bat systems is determined, the maximum revenue in Mexican currency can be calculated by optimizing the energy flows in each period p.The linear programming formulation for energy sales is similar to the candidate N PV , but without the presence of binary variables and with adjustments made based on Equation (22).The objective is to maximize the revenue generated from energy sales throughout the planning horizon, thus the objective function is defined by Equation (23).In summary, the model for maximizing energy sales can be expressed as follows: Therefore, the model is complemented by Equations ( 5)-( 20), aiming to maximize energy sales in each period p as specified in Equation (23).

Data from Northeastern and Eastern Mexico
Two specific regions in Mexico, namely Tamaulipas, and Oaxaca states, are chosen due to their favorable conditions for the generation of PV and Wd energies.The necessary data for the analysis is sourced from "CENACE" [28], including the demand data for the Northeast zone (NE) and the Eastern zone (EA) [29], which corresponds to the year 2022.Additionally, solar irradiation and wind speed data are available for these regions.

Demand in the Northeast and East
The demand data for the Northeastern (NE) zone is available at an hourly granularity, as illustrated in Figure 2. The dataset consists of 8760 data points, representing one year of data.For the analysis, we consider that 15% of the total NE zone demand is allocated to Load 1 , which corresponds to an average value of approximately 1017 MW.Load 2 corresponds to a portion of the demand in the Eastern zone (EA).It accounts for only 5% of the total demand in that zone, as depicted in Figure 3.The average value of Load 2 is estimated to be approximately 325 MW.

Locational Marginal Price
The locational marginal price (LMP) depends on various factors, which include transmission line congestion and system losses.These factors result from the characteristics of the power grid.Historical data from the real-time market (MTR) are utilized to capture fluctuations in electricity prices, with particular attention to the LMP at two specific nodes.One of these nodes is "06PAE-400," located in the NE zone [30], as illustrated in Figure 4. On the other hand, we consider a node located in the EA zone, node "02JUD-155" [30], which is shown in Figure 5.

Global Horizontal Solar Irradiation and PV Generation
Solar Global Horizontal Irradiance (GHI) represents the total solar energy reaching a horizontal surface on Earth, combining direct radiation from the sun and diffuse radiation from the atmosphere [31].It is measured in Watts per square meter (W/m 2 ) and is a crucial indicator for assessing solar energy potential in specific locations.GHI data obtained from the Modern-Era Retrospective analysis for Research and Applications, Version 2 (MERRA-2) dataset [32].For the NE zone, the GHI data is collected near San Fernando, Tamaulipas, while for the EA zone, from Salina Cruz, Oaxaca.The average GHI values are shown in Figure 6.To estimate the generation of a PV system, Equation (3) relies on the value of GHI p as well as the area utilized by the PV system.Typically, the required is around 10,000 m 2 /MW [33,34].An alternative approach is to consider the occupied by existing PV plants in the country, which provides a more approximate estimation.Table 1 presents the areas occupied by various PV plants.The data exhibit a clear linear trend, as Figure 7 depicts.By applying the linear regression method, we derive an equation that describes the required area of a PV system based on its capacity.This relationship is expressed by Equation (24), where x represents capacity in MW and y represents area in km 2 .
Figure 7. Linear regression of approximate areas of PV systems installed in Mexico.

Wind Speed and Wind Generation "Offshore"
Wind speed is a measure of the velocity at which air flows in a particular direction.It is commonly quantified in units of kilometers per hour (km/h) or meters per second (m/s).Several factors can influence the magnitude of wind speed, including atmospheric pressure, temperature, altitude, and the topography of the surrounding terrain.Additionally, obstacles along the wind's path, such as buildings, trees, or other structures, can affect its speed [35].
The wind power generation of the Wd system is determined using the findings presented in [23].The wind speeds at specific geographical locations are taken into account: (25 • N, 97 • 30 ′ W) for the NE zone (refer to Figure 8) and ( 16• N, 94 • 375 ′ W) for the EA zone (refer to Figure 9).These wind speeds measured at a height of 10 m and are then extrapolated to a height of 100 m using the following relationship [36]:  The relationship between wind speeds at different heights is described by Equation ( 25), as discussed in [37].The variables v 1 and v 2 represent the wind speeds at heights h 1 and h 2 , respectively, while γ denotes the vertical wind shear coefficient [38].For this specific case, a value of γ = 0.11 is used, as reported in [23].The calculated wind speeds at 100 m height are stored in a vector denoted as V w .These values are then interpolated based on the wind turbine generation curves [39].The offshore wind turbine models considered in the analysis are presented in Table 2, and their corresponding generation curves can be observed in Figure 10.The generation output for each wind turbine model is obtained accordingly.These values are utilized in Equation ( 4) for the parameter G Ag k .
Energy Capacity of PV Plants, Wd, and Bat Systems The capacities used in the candidate models for N PV are detailed in Table 3. Estimating the number of wind turbines N AG is a complex process that involves various considerations [40].It requires analyzing several factors to ensure an efficient and cost-effective design.One of the initial steps is to assess the characteristics of the selected marine zone, including wind patterns, water depth, seabed topography, and weather conditions [41].These data are crucial for evaluating the energy potential and determining the optimal turbine locations [42,43].To simplify this process, existing wind farms can serve as reliable references, and their databases are considered.The number of turbines and their spatial boundaries from these databases are provided in [44] and summarized in [44,45].Hence, the average number of wind turbines from existing operational wind farms is estimated to be approximately N Ag = 50.The capital expenditure (CAPEX) for each system candidate is approximately MXN 20,442,352

System Specifications and Programming Environment
The N PV and Sale MILP were solved using the CPLEX algorithm in the "Julia" programming language.The computations were performed on a computer with the following specifications:

Results and Discussions
The decrease in the prices of PV, Wd, and Bat storage systems in recent years [46] has increased the interest in using grid-connected PV, Wd, and Bat hybrid systems, from minimizing the electricity consumption bill to integrating them into a power system.A model that determines the profitability of an investment in renewable energies, such as PV, Wd, and Bat, is considered to obtain an accurate and financially viable evaluation of the project, as well as the value of the assets over time, which means that the future cash flows generated by the project are discounted to present value.In this paper, four test simulations were carried out, which are described below: In this scenario, we evaluate a flexible selection of PV-Wd-Bat systems, focusing on optimizing the systems located in the NE zone.Load 1 represents 15% of the demand of the NE zone, while Load 2 is 5% of the EA zone demand.The estimated value for N PV , as calculated by Equation ( 22), is approximately USD 1.6791 × 10 11 .This calculation considers a planning horizon of 30 years, an interest rate of 8.5%, a grid usage value of USD 2.7857/MW, and an annual increase in demand of 6%.The resulting selected systems are summarized in Table 4.
In this particular case, a Bat system was not selected.However, the PV system with a capacity of 1000 MW (i = 4) and the Wd systems (k = 8, 10, 13) with capacities of 3.6 × 50, 4 × 50, and 5 × 50, respectively, were chosen.Consequently, an offshore wind farm of 630 MW is optimally selected, comprising 150 wind turbines.The computation time for this selection process was 6811.09s.The NE zone exhibited an average demand of Load 1 : 1017.21MW, resulting in an average savings of 424.64 MW.Similarly, the EA zone had an average of Load 2 : 325.42 MW, leading to savings of 264.60 MW.

Case 2: Hybrid System with Load 1 in EA Zone and Load 2 in NE Zone
In the following case, we propose to analyze a different scenario where Load 1 and Load 2 are exchanged, along with the relocation of the PV system in the vicinity of Salina Cruz and the Wd system at coordinates (16 • N , 94 • 375 ′ W).With these adjustments, a flexible selection of PV-Wd-Bat systems is located in the EA zone, where Load 2 represents 15% of the NE zone demand, and Load 1 represents 5% of the EA zone demand.The value of the objective function (22) for this case is USD 1.3110 × 10 11 .The systems selected by the optimization process are presented in Table 5.
The combination includes PV systems (i = 1, 4) with a total capacity of 1150 MW and the Wd systems (k = 8, 10, 13, 16), totaling approximately 630 MW with 50 wind turbines, and a running time of 6085.53 s.In this case, it is possible to appreciate the flexibility of the model to determine combinations of systems according to the demand requirements, thus having an average demand for Load 1 : 325.42 MW of the EA zone, resulting in an average saving of 128.56 MW.For Load 2 : 1017.21MW of the NE zone with an average saving of 685.39 MW.

Case 3: Hybrid System with Load 1 in NE and Load 2 in EA Zone, Forced
Considering that the Bat systems were not selected in any of the previous cases, the constraint ( 1) is modified as follows: Therefore, Equation ( 26) modifies the feasible region to seek an optimal combination of PV-Wd-Bat systems, i.e., the model is constrained to find the optimal system.A value is assigned to the parameters N PV = N Bat = N Wd = 1.With this modification, the results considering the selection of one system of each technology, with the loads and GD systems as in Case 1, the optimal selected systems are shown in Table 6.
In this case, the optimal combination yielded an objective function value of USD 1.4609 × 10 11 , with the PV system with a capacity of 1200 MW (i = 6) and the Wd system (k = 10) model "Siemens SWT-4.0-130" with a capacity of 4 × 50 MW, i.e., 200 MW and a Bat system (j = 1) of 1000 MWh, with a computation time of 27,208.08s.The average demand for Load 1 is 1017.21MW with a savings of 388.84 MW.For Load 2 , the average demand was 325.42 MW, yielding savings of 208.01 MW.

Case 4:
Hybrid System with Load 2 in EA Zone and Load 1 in NE Zone, Forced Similar to the previous case, the model is adjusted to select a single candidate from the possible systems using Load 1 15% of the demand of the NE zone and Load 2 5% of the demand of the EA zone to determine the optimal selection of systems located in the EA zone, resulting in the following optimal candidates presented in Table 7.
This simulation determines that the optimal combination comprises a PV system with a capacity of 1500 MW (i = 9) and the Wd system (k = 10) model "Siemens SWT-4.0-130" with a capacity of 4 × 50 MW, i.e., 200 MW and a Bat system (j = 1) of 1000 MWh, with a computation time of 7695.59 s.Moreover, the average demand for Load 1 is 325.42MW, resulting in savings of 502.11MW, while for Load 2 , the average demand was 1017.21MW, yielding savings of 126.97 MW.

Economic Analysis of the Four Cases
This section aims to analyze each case using the energy sales program, characterized by the objective function of Equation (23).A distinct analysis from that in [3] is undertaken to compare and determine the optimal location among the four cases, which would be the best location for the generation systems; Table 8 shows that the cases with higher profitability, that is, higher net present value, are Case 1 and Case 3, which are the cases of the generation systems located in the NE zone.Additionally, Case 2 shows greater energy savings, demonstrating the highest benefit from the energy sales program.Finally, Case 4 presents the lowest N PV compared with the rest of the cases.To compute the payback time, we use the ratio between the CAPEX and the value of the objective function representing the maximum energy sales program.Consequently, the best location to install the hybrid system is in the NE zone of Mexico (the state of Tamaulipas), according to Case 1.

Graphs of the Daily Average Operation of Case 1
In many regions, solar and wind energy production peaks do not coincide, allowing for greater and more efficient utilization of available natural resources.So, the combination of PV and Wd systems can provide a constant and reliable energy source, with each source compensating for the other's low production.This complementarity is of great interest as hybrid systems can reduce long-term operating and maintenance costs by spreading the supply between different types of technologies.Adjusting the capacities of the previously selected candidates in Case 1 to have a single hybrid system PV = 1000 MW and one of Wd = 630 MW without a storage system, the objective function (23) has an energy sales value of approximately S = USD 6.7038 × 10 9 .The results obtained represent the averaged behavior over 24 h.In this case, Figure 11 shows the behavior of the energy flows supplied by the 4 different sources to the local load, or Load 1 .On the other hand, Figure 12 graphically illustrates the areas of energy contributions to Load 1 , where the grid contributes the most energy with 58.25% of the load, followed by the PV system supplying 32.11%, and finally the Wd system contributing only 9.63% of the total Load 1 .The daily average of energy purchase and sale to Load 1 is presented in Table 9.Similarly, Figure 13 shows the behavior of the flows reaching Load 2 through the grid.In Figure 14, the energy contributions supplied to Load 2 are depicted.The Wd system contributes 52.5%, the grid delivers 18.7%, and finally, 29.69% of the energy is supplied by the PV system.
On the other hand, the average energy purchase and sale to Load 2 is shown in Table 10.

Graphs of the Daily Average Operation of Case 3
The candidates have been adjusted to those shown in Table 6, a PV system = 1200 MW, a Bat storage system = 1000 MWh, and a Wd system = 200 MW.The objective function has an approximate sales value of S = USD 5.3480 × 10 9 .Additionally, the average energy in Load 1 and Load 2 are shown in Figures 15 and 16   This new PV, Wd, and Bat system configuration provides a different average energy sale than Table 9.The average sale is shown in Table 11.In systems that integrate a PV and Wd system, a storage system can help manage and provide an alternate power supply to the PV and Wd systems when they are not generating.Therefore, unlike Figure 12, in Figure 15 the integration of the batteries contributes 1.71% to the Load 1 supply, thus reducing the contribution of the energy supplied by the Wd system to a value of approximately 1.08%, as well as 34.42% of the PV system, achieving a 37.21% reduction in the consumption of the electrical system.To illustrate the behavior of the Bat storage system, Figure 17 shows the behavior of the average daily loading and unloading in the storage system, where it can be seen that the system loading occurs on average during the first 16 h of the day and the unloading occurs at 20 h.Another key consideration is the state of charge (SOC).This value is significant as it directly influences the efficiency of the Bat system since operating at optimal levels of charge and discharge benefits the system's life and allows for better planning and energy management.For ease of visualization, Figure 18 only depicts the behavior of the SOC during the first 1000 periods (h).On the other hand, Figure 19 shows the daily average behavior of the SOC.This is highly beneficial, as an adequate model and optimal behavior of the SOC significantly improve the performance and profitability of the entire hybrid system.

Sensitivity Analysis on Load, Capex, FOM, and Interest Variation
In the context of mixed-integer linear programming (MILP), the complexity of sensitivity analysis arises primarily due to the discrete nature of integer variables; this characteristic is called "Non-convexity".Due to the requirement to modify the integer variables, slight variations in the coefficients of the variables or constraint constants may cause sudden changes in the best solution.Hence, the nonconvex property makes performing sensitivity analysis by standard methods challenging.As a result, using simulations to conduct an experimental sensitivity analysis is the best course of action because it enables changing the parameters and objectively analyzing the results.A sensitivity analysis of the optimal results is carried out for Case 1, varying by 25% the demands at the upper and lower limits, then the load factor for Load 2 is FL 2 = (0.025, 0.0375, 0.05, 0.0625, 0.075) and load factor for Load 1 is FL 1 = (0.075, 0.1125, 0.15, 0.1875, 0.225).Similar techniques have been utilized in other studies [48][49][50]   According to the sensitivity of CAPEX, a similar methodology is used as the demand variation sensitivity analysis mentioned above.In this case, the analysis is feasible since a variable is varied without becoming a nonlinear component of the objective function.The substitution is as follows: as the FOM value of the three technologies, PV, Bat, and WD, a value that is directly proportional to the CAPEX is used, i.e., where α has a value of 0.5% for the PV systems, 1.5% for the Bat storage system and 2.61% for the WD system.Therefore, the objective function of Equation ( 23) can be reformulated as follows: In this way, the CAPEX is included in the objective function of Equation ( 28), with X representing the CAPEX variation variable for each technology, which is assigned values of X = −2 × 10 6 , −1 × 10 6 , 0, +1 × 10 6 , +2 × 10 6 .Therefore, as the value of X increases, the value of the objective function S * decreases accordingly.Figure 22 displays the simulation results.Finally, to examine how the objective function in Equation ( 23) behaves when the interest rate changes, we define S(r) as a function of the interest rate r given by expression (29) as follows: Once again, it has been experimentally verified by simulating this case with the following values for r = (3/100, 4/100, 6/100, 8.5/100, 10/100, 16/100), the Figure 23a demonstrates how the objective function S decreases linearly for H = 1, i.e., for one year.In contrast, for H = 30, i.e., 30 years, S decreases exponentially as shown in Figure 23b.It is worth mentioning that the sensitivity analysis is consistent with the findings reported in the literature across different settings [51].

Conclusions
This paper underscores the need to address the MILP problem in the power sector, utilizing an economic approach that reflects the present value of future energy usage if a hybrid system is invested in, based on reliable and valid historical data.The problem is framed as a significant future decision-making strategy using Net Present Value.This study develops a model to optimize the combination of PV, Wd, and Bat systems for power generation in two specific zones in Mexico.The limitations and recommendations for future work on this model are summarized as follows: Limitations and Future Works

•
The model assumes a constant inverter efficiency that does not account for variability under different conditions.Developing a model that incorporates the dynamics of inverter efficiency is advised.

•
The inclusion of uncertainty of irradiance, wind, and electricity demand data has yet to be addressed, prompting a move towards a stochastic programming approach.

•
The interaction of the hybrid system with the electric power system network infrastructure is optimal for power flow in DC or AC.

•
Incorporating more details in the modeling leads to Mmixed-integer nonlinear programming (MINLP) problems.Therefore, it is advisable to utilize relaxation and decomposition techniques to manage complexity and enhance the model's decisionmaking feasibility.
The model considers various factors such as demand, solar irradiation, wind speed, and financial parameters to determine the optimal system configurations.The results obtained from the optimization model demonstrate its effectiveness in achieving the desired objective and the optimal combination of systems for each case.In Case 1, the model recommended a PV system with a capacity of 1000 MW, coupled with offshore Wd systems with capacities of 180 MW and 200 MW.The selected configuration resulted in significant savings, with an average reduction of 424.64 MW in Load 1 and 264.60 MW in Load 2 .In Case 2, where the locations of Load 1 , Load 2 , PV system, and Wd system were reversed, the model suggested a PV system near Salina Cruz and Wd systems at specific coordinates.The objective function value for Case 2 was slightly higher than in Case 1, indicating a marginally different optimal solution.Besides, by varying the demand values, the sensitivity analysis reveals that the model's recommendations adjusted accordingly, showing its adaptability and robustness.The profitability of the optimal configuration was evaluated by considering factors such as capital costs, operation and maintenance costs, and financial parameters.On the other hand, Case 3 includes a storage system in which the energy stored in the Bat system can be used during periods of high demand or when energy costs are higher, thus maximizing energy sales.However, the inclusion of this Bat system achieves a relevant supply in Load 1 since it reduces the energy supplied by the Wd system in Load 1 , and this will be used by Load 2 , showing a new energy management.The analysis confirmed the economic viability and potential return on investment of the optimal system combinations.Overall, this study showcases the functionality and effectiveness of the proposed model in determining the optimal combination of PV, Wd, and Bat systems for power generation.The model's ability to adapt to different demand scenarios and provide economically viable solutions makes it a valuable tool for decision-making in the renewable energy sector.This is of great help since Mexico's energy policy seeks to diversify the energy matrix, increase the participation of renewable energies, and reduce dependence on fossil fuels [33].The results of this analysis indicate that the state of Tamaulipas in the NE zone is particularly suitable for the installation of hybrid energy systems, aligning with the policies of taking advantage of natural resources efficiently.Additionally, Case 1 is not only the most profitable but also shows a shorter payback time, incentivizing investment in renewables.This supports sustainability goals and the attraction of private investment in the energy sector.

Figure 8 .
Figure 8. Wind speed at 10 m from NE zone.

Figure 9 .
Figure 9. Wind speed at 10 m from the EA zone.

3. 1 .
Case 1: Hybrid System with Load 1 in NE and Load 2 in EA Zone

Figure 12 .
Figure 12.Energy supplied to Load 1 and its respective percentages of Case 1.

Figure 14 .
Figure 14.Energy supplied to Load 2 and its respective percentages Case 1. .

Figure 15 .
Figure 15.Energy supplied to Load 1 and its respective percentages Case 3.

Figure 16 .
Figure 16.Energy supplied to Load 2 and its respective percentages Case 3.

Figure 17 .
Figure 17.Energy balance in Bat system.

Figure 18 .
Figure 18.State of Charge during the first 1000 periods (h).

Figure 19 .
Figure 19.Daily average State Of Charge.
. The objective function (23) was more impacted by the variation in Load 2 , as illustrated in Figure 20 when increasing the demand in Load 2 , decreases the objective function (23), contrasts with the situation shown in Figure 21 in which the Sale is not affected by a variation of almost 50% in Load 1 .

Figure 20 .
Figure 20.Sensitivity analysis of the Sale in Case 1, varying Load 2 .

Figure 21 .
Figure 21.Sensitivity analysis of the Sale in Case 1, varying Load 1 .

Figure 23 .
Figure 23.Combined sensitivity analysis of the interest rate: (a) sensitivity analysis in one year; (b) sensitivity analysis in thirty years.

Table 1 .
Approximate areas of installed PV systems in Mexico.
.63/47]for PV systems,MXN 16,187,500.00/MW for Bat systems, and approximately MXN 34,799,980.34/MW for offshore Wd systems.The Fixed Operation and Maintenance costs (FOM) are estimated based on a percentage of the corresponding CAPEX values.Specifically, 0.5%, 1.5%, and 2.61% of the CAPEX for PV, Bat, and Wd systems, respectively, are considered[46,47].

Table 3 .
Capacities of PV power plants, Wd, and battery systems.

Table 8 .
Comparison of the four Cases.

Table 9 .
Daily average energy sales supplied to Load 1 Case 1.

Table 10 .
Daily average Sale of energy supplied to Load 2 Case 1.

Table 11 .
Daily average sale of energy supplied to Load 1 and Load 2 in Case 3.

j
Fixed operation and maintenance costs of PV i system FOM BatFixed operation and maintenance costs of Bat j system FOM Wd

k
Fixed operation and maintenance costs of Wd k system G Upper limit of SOC for Bat j systems SOC m

j
Lower limit of SOC for Bat j systems V wWind speed set to 100m CR j Maximum charge of Bat j system DR jMaximum discharge of Bat j system EP1 p Electricity price in period p at node 1 EP2 p Electricity price in period p at node 2 N usage Network usage fee in the electrical system Energy from the grid to Load 1 in period p PGL 2Energy from the grid to Load 2 in p Energy injection in p from PV i system to Load 2 P W L 2 Energy injection in p from Wd k system to Load 2 P PB ijp Energy injection during period p from PV i system to Bat j system p p p p