Minimizing Harmonic Distortion Impact at Distribution System with Considering Large-Scale EV Load Behaviour Using Modified Lightning Search Algorithm and Pareto-Fuzzy Approach

This research is focusing on optimal placement and sizing of multiple variable passive filter (VPF) to mitigate harmonic distortion due to charging station (CS) at 449 bus distribution network.There are 132 units of CS which are scheduled based on user behaviour within 24 hours, with the interval of 15minutes. By considering the varying ofCS patterns andharmonic impact,Modified Lightning Search Algorithm (MLSA) is used to find 22 units of VPF coordination, so that less harmonics will be injected from 415V bus to the medium voltage network and power loss is also reduced. Power system harmonic flow, VPF, CS, battery, and the analysis will be modelled in MATLAB/m-file platform. High Performance Computing (HPC) is used to make simulation faster. Pareto-Fuzzy technique is used to obtain sizing of VPF from all nondominated solutions. From the result, the optimal placements and sizes of VPF are able to reduce the maximum THD for voltage and current and also the total apparent losses up to 39.14%, 52.5%, and 2.96%, respectively. Therefore, it can be concluded that the MLSA is suitable method to mitigate harmonic and it is beneficial in minimizing the impact of aggressive CS installation at distribution network.


Introduction
The vision to have less carbon dioxide (CO2) emissions and less dependency on natural resources has encouraged Electric Vehicle (EV) to become an important option compared to conventional vehicle.Based on report in 2015, around 56.4% global crude oil is used for transportation sector [1].Furthermore, the unstable price for crude oil has also influenced EV to become suitable alternative.The demand on EV, indirectly, has increased the number of CS installation in the distribution system [2].However, power losses will increase when a large number of unplanned CS are installed.Furthermore, it will also introduce harmonic distortion due to the power electronic devices that convert alternating current (AC) to direct current (DC) at CS [3,4].This harmonic distortion will cause negative impact such as increment heating losses, shorter insulation lifespan, increased temperature and insulation stress, decreased power factor, and lower efficiency [5,6].
Currently, most of the researchers focus on strategies to reduce loss by coordinating the charging of EVs.Alonso et al. [7] have developed an optimization algorithm to coordinate the charging of EVs using Genetic Algorithm (GA).The charging schedule is based on optimal load pattern for EV charging-based reliability which takes into consideration the thermal line limits, the load on transformers, voltage limits, and parking availability patterns.The proposed method is able to reduce the cost for power system new investment.Furthermore, Bharati and Paudyal [8] have proposed framework to coordinate EVs in the distribution system using bilevel hierarchical vehicle-grid (VG) optimization.This framework requires information exchange among EV aggregators and grid controller and then the EV charging can be scheduled based on the minimum losses to the distribution network.Next, Masoum and Nabavi [9] have also proposed onlineovernight PEV coordination using metaheuristic technique which is particle swarm optimization.This method introduces high priority and low priority customer that indicate  the willingness of customer to be charged.High priority customer will get faster charging services with high tariff and the low priority customer will only charge their vehicle after power loss control in the distribution system is done.Even though many approaches have been introduced to minimize the EV impact, the focus is more towards minimizing the power loss.However, beside power loss, the harmonic injection to the distribution system when many CS are connected must also be taken into consideration in ensuring that the network can operate optimally.
In general, there will be two types of harmonic that will be mitigated in distribution network, which are voltage and current harmonic.Voltage harmonic is measured at the bus, while current harmonic is measured at the lines and cables.Passive filter is the common component used to eliminate these harmonics.Passive filters can provide a low-impedance path to harmonic currents and hence prevent it from flowing into the systems.There are many papers that utilized passive filter to mitigate harmonics in power system using filter devices [10][11][12][13][14] and the most common ones are single tuned filters.Single tuned filters are designed to eliminate or reduce single frequency from the system depending on the design resistance, capacitance, and inductance value [15].Sakar et al. [16] have suggested a hosting capacity determination for a distorted distribution system due to Photovoltaic connection.In their research, passive filter is used to increase the harmonic-constrained hosting capacity which then improves the voltage, power factor, and filtering the harmonics.Author in [17] proposed an Asymmetric Synchronous Reference Frame Control scheme and harmonic voltage compensator with a Static Var Compensator connected to the grid to mitigate the harmonic.The author used harmonic voltage compensator to reduce voltage THD.Even though many researchers have used passive filter to mitigate the harmonic, very few researchers considered both voltage and current harmonics in their analysis especially for CS.
Therefore, this research proposes a coordination of 22 VPF when 132 unit of CS is adopted in a 449 bus radial distribution system.Section 2 presents the modelling of load profile, CS, battery, passive filter, harmonic, CS scheduling, and radial load flow used in this study.The introduction of MLSA technique will also be covered in this section.The methodology of the proposed coordination will be discussed clearly in Section 3. Result and discussion for 22 units of filter coordination in 449 bus are presented in Section 4. Lastly, conclusion of this research is in Section 5.

Problem Modelling and Formulation
A typical 449 bus radial distribution system is used in this research to determine the optimal placement and sizing of VPF.Furthermore, forward/backward sweep method is chosen due to the accuracy of this load flow analysis in solving distribution load flow [18][19][20][21].The harmonic pattern for individual EV charger is modelled based on actual single phase charger impact as measured in [22].most typical and simplest one is by using current injection analysis method [26].The electrical parameters that affected the harmonic existence are the line impedances, load variation impedances, and filter impedances values.Equations ( 1) and (2) represent the line impedance formulation and the impedance for single-tuned filter at harmonic ℎ, respectively, which indirectly caters for resonance impact.Equations (3), (4), and (5) show the impedance formulation for the load used in setting up the harmonic admittance matrix.The harmonic flow is calculated using (6), which consists of harmonic admittance matrix and harmonic current injection as per (7) and (8).SOC state will be updated in every 15 minutes (Δ) and can be calculated using ( 9), while the current produced by CS can be calculated using (10).
where  oc, is open circuit voltage for th node (V),   is rated battery ampere hour for the th PEV (Ah),   is battery equivalent internal resistance for the th node (ohm), CR (Δ  ,) is maximum charging rate for the th PEV (A), SOC(Δ  , ) is state of charge of the th PEV at th time slot (%), SOC(Δ +1 , ) is state of charge of the th PEV at next th time slot (%), (Δ  , ) is charging current for the th PEV at current time slot (A), and  CS (Δ  , ) is consumed power for the th PEV (kW).
Figure 3 shows the open circuit voltage (OCV) for lithium iron phosphate (LFP) batteries based on the SOC of battery [27] that can also be presented using (11).Thus, the relation between SOC and  OCV will influence the total current injected during charging process.
Since there is conversion from AC to DC during charging process, the harmonic impact to the network is considered in this study.The harmonic injection is basically based on SOC state, which depended on battery SOC status and is being modelled based on practical data from [18].Since there are 132 CS units installed within the 22 low voltage residential areas (6 units for every residential area) it might have significant  impact on THD  and THD  value when all CS operated simultaneously.All CS and battery data in the residential area are presented as per Table 2.

Modelling of Load Profile and Charging Station Operation.
The load profile for 449 bus system is design based on typical load consumption in residential area during normal day.The peak demand in the system is occurring at the afternoon due to the several usages of electrical appliances, while the lowest demand is at the night when majority of people are asleep.However, the new electrical load CS is normally being charged at the night time.In this study, it is assumed that harmonic produced in the distribution system is coming from CS. Figure 4 shows the loading produced by all CS in the distribution system based on design of CS operation as in Table 3. CS operation is designed based  on random behaviour and will penetrate at the night hour.
Figure 5 shows the summation of original load profile with CS load.Even though most of the CS is operating during low loaded condition at distribution system, it will cause higher impact on the total harmonic distortion as well as power losses.

Assumptions and Constraints.
There are few assumptions and constraints are considered in this research as shown below: (i) CS can be plugged in/plugged out at any time according to the customer's request.Customers will input their requested plug-out time and requested final SOCs at the time of plug-in.Once SOC reaches requested SOC, CS will be switched to a standby mode.
(ii) The time slot is 15 min (Δ) or equal to 96 slots for one day.
(iii) The aggregator has access to CS information including their bus locations, harmonic distortion information, charger types, battery sizes, plug-in time, and plug-out time.
(iv) CS are controllable and have variable charging functions.During the charging process, each CS is assumed as a variable active load.The power that is used to charge the battery is based on calculation that is done by the aggregator.(v) The requested time for each CS must be greater than the minimum charging time required to charge the battery.
Since 22 variable filters will be placed in the network, the total unknown variables will become 44, which are 22 locations and 22 optimal sizes.All these parameters will have their own constraints that need to be considered.In general, the parameters can be divided into 2 categories which are as follows: (i) Filter location: the filter will be placed in low voltage bus to avoid any harmonic injection to upper feeder.There are 22 locations in 449 bus radial distribution system.Furthermore, only one VPF will be placed at every low voltage 415 buses system.The constraints are as follows: (ii) Filter reactive value: each individual variable filter reactive value is limited to 2.6.Modelling of Passive Filter.The main function of a passive filter is to sink the harmonic current that flows in the system based on a selected frequency.The filter impedance will become very low to allow the harmonic to sink.Single tuned filter is the most popular type of filter which is used widely in dealing with harmonic pollution especially in the industrial area [28,29].In this research, four units of single tuned filters are considered as one set of filter which can eliminate four frequencies as shown in Figure 6.Equation in [26] is used to calculate capacitor, inductance, and resistance components, respectively, as per (14).Capacitor is calculated based on injected reactive power () and voltage () at that bus; meanwhile, inductance and resistance are based on the chosen harmonics () that need to be reduced.Four sets of filter will be used to eliminate 3rd, 5th, 7th, and 9th harmonic order in the network.
2.7.Multiobjective MLSA with Pareto Optimization.In this research, the optimization technique is used to identify optimized locations and sizes for 22 sets of filters in 449 bus radial distribution system.The Lightning Search Algorithm (LSA) is a metaheuristic optimization method inspired by lightning propagation from cloud to ground [30].The LSA process basically starts with generation of random population before fitness calculation.The worst step leader will be eliminated every 10th iteration before the direction of the step leader is updated.Next, step leader movement is updated based on direction, shape parameter, and scale parameter.Lastly, forking phenomenon is taken place as per (7) for 1% of the better step leader.Since the problem faced in this research is too complicated with many local minima, modification on the existing LSA is needed.In this research, the MLSA is proposed with better convergence for this application compared to existing LSA.MLSA basically is the improvement of LSA method which is proven beneficial for this research.
The placement and sizing of variable filter will be based on five parameters, which are maximum THD  , maximum THD  , THD  summation, average THD  , and apparent power losses for overall system.The three parameters that involve THD  are to ensure that the THD  at all buses is within acceptable range.These parameters are then normalized to get the most accurate solution as defined by ( 15)- (19).The weight summation with normalized fitness function is shown in (20) after several manual trials are done in order to get the best optimal result.Last but not least, basic equations to calculate THD  and THD  are as shown by ( 21) and (22).
Fit = 0.1fit 1 + 0.2fit 2 + 0.2fit 3 + 0.2fit 4 + 0.3fit 5 , (20) In order to reduce computational time in this research, bus location is gathered from MLSA with weight summation approach; meanwhile, Pareto optimization is used to determine size of VPF.Normally, Pareto technique will produce a "nondominated" solution which has more than a single solution [31].Thus, fuzzy satisfying approach is adopted in this research to get the best solution from the "nondominated" solutions.Equation ( 23) is the weight membership function used in fuzzy approach for all "nondominated" solutions in the fuzzy set.The fit max in the formula represents the existing fitness value before implementing the filters, while  represents the number of solutions listed in the "nondominated" list.The highest   value indicates the best solution for this research.

Methodology
The variable filter is used to minimize the impact of harmonics due to large-scale deployment of CS in the big distribution system.The numbers of filter are based on how severe the harmonic distortion and power losses on the existing system.By placing 132 CS units (6 units for every low voltage residential area) at specific locations, MLSA will be used to determine optimal placement and size for 22 sets of filters (1 set for every low voltage residential area).The THD  , THD  , and  loss will be recorded based on specify CS and filter locations.The overall process to evaluate system performance is as shown in Figure 7.
There are four modifications on MLSA that are made to improve LSA method.The first modification is on channel time; the original channel time is not suitable when dealing with the problems that have many local minima.The second modification is the updating approach; a new updating approach able to check the forward or backward direction is proposed.The third modification is on the scale parameter; in LSA, it is based on exponential distribution, which will cause the step leader movement to be active at 20% of the early iterations as per (24).Thus, in MLSA, the Laplacian distribution equation is used to increase the movement of the step leader during 30% until 70% iterations.Last but not least, the learning factors which make the movement of the step leader more active is implemented in the fourth modification as per (25).where   is the scale parameter at th iteration (1/3) −(−50/15) ,  is learning factor (2.0 in this research),    new is new position in lead projectile at th iteration,    is position in lead projectile at th iteration,   is shape parameter at th iteration,   is the position at th iteration, and  best is the best individual during minimum.
Figure 8 shows the flowchart for MLSA optimization process using weightage summation approach for this research.
In this process, simulation involves 500 iterations and 50 populations.Next, the bus location results will be used in the next stage by the MLSA-Pareto-Fuzzy combination technique.Figure 9 shows the flowchart for MLSA optimization process with Pareto-Fuzzy technique.In this stage, Pareto technique will be adopted in finding the best solution, while fuzzy approach will be used at the end of process to find the best solution among nondominated solutions.

Result and Discussion
The simulation on this research was done based on the scenario where 132 CS units are installed in low voltage 415 V buses for 24 hours (96 states, 1 state = 15 minutes).The schedule for CS operation can be referred to in Table 3.The THD  and THD  for all medium voltage buses and lines are shown in Figures 11 and 12, respectively, for sampling time at 30 and 63.Sampling time 63 is chosen due to the worst harmonic recorded in 24 hours, while sampling time 30 is chosen due to the lowest harmonic injection.From the result, without VPF, maximum THD  for sampling time 30 is recorded at bus 15 with 0.8750% and maximum THD  is recorded at lines between buses 26 and 27 (line number 27) with the value 0.6780%.For sampling time 63, maximum THD  is 0.9032% at bus 15, while maximum THD  Next, 22 sets of VPFs are placed at low voltage 415 V buses using MLSA optimization with multiobjective function.The optimal location and size for all filters are tabulated in Table 4, while the parameter value is shown at Table 5.From the results in Table 5, there is a significant reduction for maximum THD  from 1.258004% to 0.74478% and maximum THD  from 2.019816% to 1.16256% after MLSA get the optimal placement for VPF.However, the total apparent losses in the network are increased.Thus, a Pareto-Fuzzy approach is used to cater this issue at the next step.
The results from Pareto consist of multiple set solutions where the possible solution is the best solution as long as one of their objectives functions is dominant compared to other solutions.Fuzzy approach will be used to find the best among solutions in Pareto set.In this study, the solution consists of THD  , THD  (minimum, average, and summation), and power loss.From the analysis, the total Pareto result for sampling time at 63 in MLSA is 137 sets.The fuzzy approach is used in order to get the best single solution.The highest value of weight membership function in fuzzy approach will be chosen as the best solution in this research.Tables 6  and 7 show the five best solutions with the highest   for sampling times 63 and 30.From the result, the best sizing for all filters is shown in Table 8 at sampling time 63.The size of filter to cater harmonic is different between times 30 and 63 which indicates that the size of passive filter is very crucial in this study.The THD  and THD  for buses and lines, after the filter implementation as per MLSA with Pareto-Fuzzy best solution, are shown in Figures 10 and 11.
From the results, it shows significant reduction of THD  at all medium voltage buses, while majority of THD  also shows reduction.From the best solution, the THD  and THD V for the 449 bus distribution system have been reduced effectively after filter installations.The results for all 96 slots are shown in Figures 12-16.Based on Figures 12 and 13, maximum THD  and THD  show reduction for all states which indicate that the system is better compared to previous value, while apparent losses at Figure 16 show the capability of this method to also reduce losses together with harmonic.Next, Figures 14 and 15 also show reduction which assist THD  improvement in overall system indirectly.From the result, implementation of VPF with changes every 15 minutes is able to reduce harmonic for 24 hours.Other than that, it is obvious that the implementation of MLSA with Pareto-Fuzzy is able to identify the best solution that provides the greatest improvement to the overall system.
Based on the collective result gathered from all 96 simulations that represent 24 hours, minimum and maximum values of VPF at dedicated bus are acquired.Table 9 shows the minimum and maximum VPF required to solve harmonic problem after EV is installed.

Conclusion
Optimal placement and sizing of VPF in the distribution system can avoid harmonic injection to the medium voltage  network as well as reducing apparent power losses.This research has shown that MLSA with multiobjective function and Pareto-Fuzzy are able to determine appropriate filter location and size to reduce harmonic distortion and apparent losses in 449 bus system with 132 units of CS.Furthermore, the placement of twenty-two sets of four unit single-tuned filters is able to reduce the four harmonic orders.Based on simulation run for every 15 minutes, the proposed technique is able to improve maximum THD  , maximum THD  , and  loss up to 39.14%, 52.5%, and 2.96%, respectively.The study is very important for future distribution grid that might have many CS in the network.For future work, the study can be expanded to other types of passive filter able to have better impact on the system such as C-type filter and others.

Figure 1 :
Figure 1: Single-line diagram for 449 bus radial distribution system.

Figure 3 :Figure 4 :
Figure 3: Relationship between open circuit voltage and SOC for LFP battery.

Figure 5 :
Figure 5: Load profile before and after CS operation.

Figure 8 :
Figure 8: Flowchart for MLSA method to get minimum fitness with weightage summation approach.

Figure 9 :
Figure 9: Flowchart for MLSA method to get minimum fitness with Pareto-Fuzzy technique.

Figure 10 :Figure 11 :
Figure 10: THD  at medium voltage for times 30 and 63 with and without passive filter.

Figure 15 :
Figure 15: THD  summation for 24 hours with and without passive filter.

Table 1 :
Line data for 415 V buses.

Table 2 :
CS and battery data.

Table 4 :
Bus location and filter size after 500 iterations.

Table 5 :
Result using MLSA with multiobjective function.

Table 6 :
The five best solutions based on the highest   using fuzzy approach for 137 sets of solutions (sampling time = 63).

Table 7 :
The five best solutions based on the highest   using fuzzy approach for 10 sets of solutions (sampling time = 30).Figure13: THD  maximum for 24 hours with and without passive filter.THD  average for 24 hours with and without passive filter.

Table 8 :
Filter size based on the best solution after fuzzy at sampling times 30 and 63.

Table 9 :
Minimum and maximum VPF size based on 96 simulations.