Optimization of wind farm operation with a noise constraint

. This article presents a method for performing noise constrained optimization of wind farms by changing the operational modes of the individual wind turbines. The optimization is performed by use of the TopFarm framework and the PyWake wind farm modeling as well as two sound propagation models: the ISO 9613-2 model and the Parabolic Equation model, WindSTAR. The two sound propagation models introduce different levels of complexity to the optimization problem with the WindSTAR model taking a broader range of parameters, like the acoustic ground impedance, the complex terrain 5 elevation and the flow field from the noise source to the receptor, into account. Thus, as the WindSTAR model introduces a higher complexity of the sound propagation computations, it likewise introduces a higher computational time. Wind farm optimization using each of the two sound propagation models is therefore performed in different atmospheric conditions and for different source/receptor setups, and compared through this study in order to evaluate the advantage of using a more complex sound propagation model. The article focuses on artificial wind farms in flat terrain as well as arbitrarily chosen dwellings at 10 which the noise constraints are applied. By this, the study presents the potential of an optimization algorithm focusing on the sound propagation and wind farm operation trade-off.


Introduction
As the demand for onshore wind farms increases, the social acceptance of wind turbines becomes a larger challenge.One of the main factors contributing to neighbour annoyance is the aerodynamic noise from the wind turbine blades.Previous social studies have shown how neighbours to wind farms experience annoyance and sleep disturbances caused by the noise emitted from wind farms (Michaud et al. , 2016;Poulsen et al. , 2019Poulsen et al. , , 2018)).Thus, in order to successfully continue the expansion of onshore wind energy by either constructing new wind farms or by repowering existing ones, methods to ensure a low noise level are needed.Several wind turbine developers introduce multiple operational modes in the turbine design with the aim of switching to a more noise reducing operation by for example decreasing the rotational speed of the rotor.However, modifying the rotational speed to reduce the emitted noise causes a curtailed power output of the wind turbine.A method for effectively choosing at which operational mode each of the wind turbines should operate during varying atmospheric conditions is therefore needed.By performing noise constrained optimization of wind farm operations, it is possible to maximize the overall power production of a wind farm while still keeping the noise received at each neighbour under a defined limit.This is done by letting each wind turbine in a wind farm individually switch to the optimal discrete operational mode.Previously, noise constrained optimization has been performed through layout design of onshore wind farms (Tingey et al. , 2017;Wu et al. , 2020;Sorkhabi et al. , 2016;Mittal et al. , 2017;Cao et al. , 2020).Furthermore, optimization of discrete design variables has previously been done through other layout optimization problems (Riva et al. , 2020;Feng et al. , 2017).However, the work in this article presents a novel way of considering the discrete operational modes in noise constrained optimization.The optimization of the operation of an existing layout can be needed in the case of repowering of a wind farm or in wind farms that are already heavily de-rated in order to limit the emitted sound.Currently, the ISO 9613-2 sound propagation model (ISO 9613-2 , 1997) is extensively used to determine the location and operation of onshore wind turbines.The model is adapted to the regulations of the specific country, by considering varying values i.e. for the acoustic hardness of the ground.Some examples of regulations in four European countries are summarized in (Nieuwenhuizen et al. , 2015).In addition, each country has specified noise limits that may vary from day to night or from area to area.In Denmark, the ISO 9613-2 model is used such that the 'worst case scenario' noise is modeled.The ground type parameter is set thus to 0, representative of a hard, reflecting ground surface.The Danish regulations further set noise limits for the wind speeds U 10m = 6 m/s and 8 m/s at 10 m height above the ground.In noise sensitive areas the limits are for example defined as 37 dB(A) for U 10m = 6 m/s and 39 dB(A) for U 10m = 8 m/s (Nieuwenhuizen et al. , 2015).The noise limits are thus in general very vaguely defined by only considering the wind speed and the site of interest.This definition is naturally originating from the limitations of the ISO 9613-2 model and the uncertainty of the measured meteorological conditions at the site.However, using a higher fidelity model, where the wind direction and a more detailed flow field from the wind turbine to the receptor as well as the complexity of the terrain elevation can all be considered, can yield the possibility of designing the layout or the operation strategy of a wind farm based on more parameters than simply the wind speed at two specific weather scenarios.Thus, the attenuation of the sound has further been shown to change with the stability of the atmosphere (Barlas et al. , 2018).
Previous studies have further shown, how the wind direction or the upward/downward refraction of the atmosphere can have an immense effect on the propagation of sound (Lee et al. , 2016;Bolin et al. , 2020;Barlas et al. , 2018;Evans et al. , 2012).
An upward refracting atmosphere can especially cause significant acoustic shadow zones in the far field of the wind turbine, and thereby result in highly reduced sound levels.In order to take phenomena like this, and in general more details about the flow and the terrain, into account, the WindSTAR model based on the Parabolic Equation (PE) method (Barlas et al. , 2017;Barlas , 2017;Cao et al. , 2022) is along with the ISO 9613-2 model used for optimization in this study.Both models are coupled to the Topfarm optimization framework (Pedersen et al. , 2021;Réthoré et al. , 2014), and the PyWake framework (Pedersen et al. , 2019) is used for the wind farm modeling.The WindSTAR model has previously been validated against field measurements of sound propagation (Nyborg et al. , 2022;Cao et al. , 2022), and compared to other sound propagation models in order to further verify the model.These studies showed overall good results from the WindSTAR model.By comparing the optimization performed by using WindSTAR to the optimization with the ISO 9613-2 model, the aim of this paper is to show the advantage of taking more parameters into account when considering the propagation of sound from wind turbines.This study focuses on the optimization of operational modes of the wind turbines, but could be transferred to layout optimization with noise constraints by considering the Annual Energy Production (AEP) of the wind farm instead of the power production in the specific flow case.

Sound propagation models
The sound pressure level, L p , of each receptor surrounding the wind farm in question sets the constraints of the optimization problem.The sound pressure level at a receptor near a source of sound is generally derived by the source strength or sound power level and the propagation of the sound through the atmosphere from the source to the receptor.The source strength of the wind turbines considered in this study is determined by the manufacturer of the wind turbine.Thus, the octave band depending sound power levels, L W , are obtained by the predefined operational modes of the wind turbine of interest.Each operational mode of the wind turbine is designed to reduce the overall emitted noise by a few dB by slowing down the rotational speed of the turbine rotor or by pitching the blades.This implies that as the operational mode is switched, the L W spectrum of the wind turbine changes.As shortly mentioned, this study involves two sound propagation models of different complexities: The ISO 9613-2 model (ISO 9613-2 , 1997) and the WindSTAR model (Barlas et al. , 2017;Barlas , 2017;Shen et al. , 2019).Following the ISO 9613-2 model, the general equation for the sound at a nearby receptor is given by: where a uniform directivity is assumed leading to D C = 0 while the attenuation parameter, A is expressed as: Only the attenuation caused by the geometrical spreading, A div , the atmospheric absorption, A atm , and the ground effects, Generally, the frequency, f , dependent sound pressure level, L p (f ), can be written as (Salomons , 2001):  ISO 9613-1 , 1993).The atmospheric absorption depends on the temperature and the relative humidity and is seen to become more dominant at longer distances and higher frequencies.The third right hand side term is the geometrical spreading of a spherical wave with S 0 given at a reference distance of d 0 = 1 m.Independent of the model used, the geometrical spreading of the sound is the major contributing factor to the attenuation of the sound.However, the sound pressure level relative to the free field sound pressure level, ∆L p , can contribute to the L p being pushed over the defined constraint at a receptor.
∆L p includes any effects from atmospheric refraction and terrain elevations.Alternatively, L p can be defined by the complex sound pressure where p f ree is the propagation of a reference point source in a free field and the complex sound pressure, p c , can be expressed by Helmholtz wave equation (Salomons , 2001).The Helmholtz equation is solved in WindSTAR through a Parabolic Equation (PE) method by use of the Crank Nicholson (Gilbert et al. , 1989;West et al. , 1992) approach and by introducing a coordinate shift at ground elevation changes, the model has been adapted to propagation over complex terrain.This method is commonly referred to as the Generalized Terrain Parabolic Equation or GTPE in short (Sack et al. , 1995).The GTPE method replaces the moving atmosphere with a hypothetical motionless atmosphere with an effective speed of sound, expressed as c ef f = c 0 + u where c 0 is the adiabatic speed of sound and u is the wind speed field projected into the plane of propagation.In addition, the GTPE model is approximated to a 2D model by assuming independence of the direction of propagation from the source.
The model is a one way propagating model, meaning that back-scattering of sound is neglected.The GTPE method is one of many existing PE models of which each model introduces an individual method for solving the system of PEs.Where some approaches like the Greens Function PE (GFPE) model (Gilbert et al. , 1993;Salomons , 1998) allow for a large grid step size in the r−direction of the computational domain, the GTPE requires a grid resolution of ∆r = ∆z = λ/8, where λ is the wave length of the considered frequency (Salomons , 2001).Hence, the resolution of the grid becomes more and more refined as the frequency increases.This wave length dependency of the grid spacing further introduces numerical issues at too high frequencies, f = 4 kHz and f = 8 kHz.The attenuation at these frequencies is therefore obtained by the ISO 9613-2 model regardless of the sound propagation model used for the remaining frequencies.Moreover, it is experienced that the memory allocation becomes too excessive at longer distances and frequencies of f = 1 kHz and f = 2 kHz.Thus, for these frequencies at distances from the bottom of the turbine tower to the receptor exceeding 3.5 km, the ISO 9613-2 model replaces the WindSTAR model as well.
The bottom boundary condition of the GTPE is defined by the acoustic impedance at the ground surface computed by the model by (Attenborough , 1985) and characterized by the flow resistivity, σ, while an artificial absorbing layer with a thickness of 50λ is assumed at the top boundary (Salomons , 2001).The height of the computational domain is set to span 500 where z hub is the height of the wind turbine hub above the ground and R is the rotor radius.This decision is made based on studies performed with point sources distributed at the rotor (Cotté , 2019;Barlas et al. , 2017) and on the source positional study by (Oerlemans et al. , 2007).In previous work with the WindSTAR model 36 point sources have been used (Cao et al. , 2022).The computational time of running 36 individual computations for each octave band frequency is however too excessive for optimization purposes, and the use of the 3 distributed point sources has previously shown good comparison with field measurements (Nyborg et al. , 2022).The 3 point sources are assumed to be uncorrelated and the average L p (f k , d ij ) of the 3 point sources from the ith wind turbine to the jth receptor can be obtained by assuming equal source strength and uniform directivity From which the attenuation of the sound will henceforth be referred to as the transmission loss, T L. All computations are done in octave band frequencies, according to the ISO 9613-2 model since the considered L W spectra are provided in this form.The overall integrated L p,j at the jth receptor is thus obtained by L p,j = 10 log 10 where k is the octave band frequency index.Although accounting for the turbulence in the atmosphere can have a significant influence on the L p,j at a receptor, this effect is not included in the optimization.Hence, only steady computations of the specific flow case are considered in order to avoid excessive computational times.The turbulence introduces increased scattering of the sound, which can result in larger sound pressure levels in regions subject to shadow zones caused by upward refraction or by complex terrain (Gilbert et al. , 1990;Bolin et al. , 2020).By not including the turbulence in the computations, a higher uncertainty of the estimated L p,j is therefore expected in the case that the jth receptor is positioned in a shadow zone.
Due to the different complexities of the two sound propagation models, the computational time varies significantly as well.
While the ISO 9613-2 model can be evaluated on a laptop, the amount of physics included in the WindSTAR model require a cluster for the computations.

Topfarm optimization framework
The optimization framework used is the Topfarm framework developed at the Technical University of Denmark (DTU) (Pedersen et al. , 2021).Topfarm was developed as a package in Python with the intention of performing economical optimization of wind farms.The famework uses the OpenMDAO package for the optimization (Gray et al. , 2019) and has its own implemented among others previously been used for layout optimization by introducing load constraints (Riva et al. , 2020).In the work done in this article, a random search optimization algorithm is used (Feng et al. , 2015), which has been adapted to discrete design variable problems (Feng et al. , 2017).The adapted random search algorithm for discrete design variable problems is modified to optimize the discrete operational modes of each discrete wind turbine.The optimization thereby switches the modes of each wind turbine in question during a defined number of iterations.The number of wind turbines in question at each iteration and the corresponding operational modes are randomly chosen by the algorithm.The choice is however to a certain extend made through heuristics by choosing the operational modes based on the modes chosen in the previous iterations.When using the random search algorithm, it should be kept in mind that there is no guarantee of the algorithm finding the global optimal solution especially as the number of design variables increases.Furthermore, it should be noted that the computational time and scalability of the random search algorithm when increasing the number of design variables may not be appropriate for the optimization of larger wind farms.In these cases a gradient based optimization method may be more appropriate (Martins et al. , 2021), which is part of the future work with the presented framework.However, the random search algorithm is easy to apply to an optimization problem considering discrete design variables and is deemed feasible for the purpose of this work, since it aims to demonstrate the potential of the type of optimization presented.

Wind farm flow modeling in PyWake
The wind farm modeling is done in the PyWake framework (Pedersen et al. , 2019) in order to apply the different implemented engineering wake models to the problem.The 2-dimensional flow fields obtained from the PyWake framework are interpolated to the grid points in the computational domain used in WindSTAR spanning from each wind turbine to each receptor.In this way, PyWake is loosely coupled to WindSTAR and wrapped with the Topfarm optimization framework.Furthermore, the power output of the considered wind farm is computed by applying the power-and C T curve of the wind turbine type.The wind farm is modeled in an iterative downstream manner, neglecting any blockage effects of the wind turbine but providing fast computations in return.In the work presented, the Gaussian wake deficit model developed by (Bastankah et al. , 2014) is used.WindSTAR has previously been coupled with the engineering wake model by Qian (Qian et al. , 2018;Cao et al. , 2022;Nyborg et al. , 2022).Besides the velocity deficit, the Qian wake model provides a turbulence intensity deficit which can be useful for noise source modeling.However, the source strength is provided through the operational modes and the Qian wake model is further not yet available in the PyWake framework, but only loosely coupled with WindSTAR through an implemented Fortran version.The implementation of the Qian wake model to PyWake was out of scope in the presented work.
Engineering wake models naturally yield a more simplified flow field than Computational Fluid Dynamics (CFD) methods such as Reynolds Average Navier-Stokes (RANS) (van der Laan et al. , 2015) or Large Eddy Simulations (LES) (Jimenez et al. , 2007).However, the engineering wake model provides a fast estimate of the wake field which can be more appropriate for optimization purposes.If a wind farm is placed in very complex terrain, more advanced modeling such as RANS computations can be implemented in order to estimate the speed up effects in the background flow field.In this case, the flow field including wake effects is obtained by superposition of the background flow field from the RANS computations and the velocity deficits from the engineering wake model.

Optimization flow
Figure 1 illustrates the flow chart of the algorithm designed for the noise constrained optimization problem when using the WindSTAR sound propagation model.As seen, the modeling in the problem is divided into two parts: the wind farm wake modeling and the sound propagation modeling.Initially, the models are provided with a flow case including the wind direction, θ, the free field wind speed, U 0 , as well as the temperature, T , and relative humidity, ϕ, needed for calculating the atmospheric absorption.Furthermore, site information like the terrain elevation, ground impedance and positions of the wind turbines and receptors are given.Lastly, the initial operational mode, m 0,i , the lower, m l , and upper, m u , bounds of the operational modes as well as the noise constraint at each receptor, L p,lim,j , are provided to the optimizer.It is noted that for some cases it may be necessary to shut down a wind turbine completely.Such a mode is however not included in the presented work.The updated operational mode, m i , of each wind turbine, i, is parsed into both model parts (red.the wind farm model and the sound model) in every iteration of the optimization, which in return parse the updated total power of the wind farm, i P i (U 0 , θ, m i ), and the updated sound pressure level at each receptor, L p,j (U 0 , θ, m i ).Thus, the updated operational mode modifies the C T of the wind turbines and thereby the wakes through the wind farm, while the sound power level is likewise dependent on the  In short, the optimization problem when using either of the two sound propagation models can be mathematically described as:     In the results presented in Figure 5, it is observed how the range dependent T L at each octave band frequency, f k , is altered as the operational mode of the wind turbine is switched.However, the effect of the change in operational mode is mostly observed in the far field of the wind turbine at distances longer than d = 2 km.Moreover, the effects become more apparent for higher frequencies of f = 1 kHz and f = 2 kHz with significant differences appearing after d = 1500 m for f = 2 kHz and after d = 2000 m for f = 1 kHz.In Figure 6 the effect on the propagation from a wind turbine in the wake of another turbine subject to changes in operational mode can be observed.The changes in the T L is generally seen to be less significant than the ones observed for the first setup in Figure 5.The most distinct differences are observed when switching the front wind turbines from mode 3 to mode 6.Still, the major effects are apparent at further distances and higher frequencies.When observing the T L at U hub = 6 m/s and 14 m/s in Figures A1-A4, the sensitivity of the change in operational modes is less significant at U hub = 6 m/s.At U 0 = 14 m/s, the T L shows changes in the far field similar to Figures 5 and 6.It is observed that especially for some frequencies at U hub = 6 m/s, the results obtained for mode 0 and mode 3 are similar, making mode 0 hardly noticeable 280 in Figure A1 and A3.Since the sensitivity to the operational modes in the presented cases is generally observed for longer distances at higher frequencies, it is considered negligible compared to the high transmission losses expected.The sound propagation modeling in WindSTAR can therefore be excluded from the iterative function calls during the optimization and instead performed separately prior to the optimization by using the initial modes, m 0,i , of the wind turbines.The T L obtained from the initial run of WindSTAR is therefore used as a transfer function in the optimizations presented in this article.receptors are positioned in arbitrarily chosen locations, but in a way such that some of the receptors will either be in the upwind or downwind positions of wind turbines in the farm.By choosing these positions, some of the distinct differences between the computed noise by WindSTAR and by ISO 9613-2 caused by refraction in the atmosphere are expected to be captured (Barlas et al. , 2018).The hub height wind speed is kept at U 0 = 10 m/s, while two different wind directions are used, θ = 0 • and θ = 225 • .In this way, the effect of all the wind turbines being in the free flow field compared to the majority of the wind turbines being in the wake can be analysed.The temperature is T = 15 • C and the relative humidity is ϕ = 80 %.For the wind profile in the free field, a logarithmic profile for neutral conditions is used where the roughness length is z 0 = 0.1 m, the von Karman constant is κ = 0. defined for the noise sensitive areas in Denmark are used for all receptors in all optimization cases presented.For the wind profile chosen, the wind speed at 10 m height is U 10m ≈ 6 m/s.The constraint is thus set to L p,lim = 37 dB(A) (Nieuwenhuizen et al. , 2015).Lastly, the initial mode of every wind turbine is set to m 0,i = 0, starting the optimization from the least noise reducing mode.As a result of the study done in Section 4.2, the optimization with the WindSTAR model is performed by conducting the initial sound propagation computations for the specified flow cases and using the results as a transfer function for all wind farm layouts in this article.

Row of 7 wind turbines
As mentioned, the optimization of the wind farm operation is done for both a row of 7 SWT-2.The convergence of the total power output, the operational modes of each wind turbine, m i , and the L p,j at each receptor during the optimization performed for the flow case of θ = 0 • are shown in Figure 8 for the 7 SWT-2.3-93turbines.For all presented wind farm layouts and flow cases the convergence is shown for both the optimization using the ISO 9613-2 model and for the one using the WindSTAR model outside of the iterative function calls.A total number of 1000 iterations have been chosen for the optimization, but it is observed that convergence is reached after around 100 iterations for the ISO 9613-2 model and after 50 iterations for the WindSTAR model.Thus, the optimization is observed to converge relatively fast.It is observed how the ISO 9613-2 model is generally predicting higher L p,j values than the WindSTAR model at the initial operational mode, m 0,i .This could be the cause of the larger number of iterations needed to reach convergence, and a lower optimized power output of the wind farm.It is observed how the optimization with the ISO 9613-2 model first overcompensates for the violated constraint, and eventually increases the L p,j at all receptors to approach the L p,lim,j and thus increasing the power output of the wind farm.The WindSTAR model predicts significantly lower L p,j for all receptors causing only a few of the turbines to switch to a noise reducing mode.
Scatter plots of the operational modes of each wind turbine along with the corresponding L p,j before, during and after the optimization are given in Figure 9.The scatter plot at the time during the optimization represents the iteration at which the estimated power output is at its minimum.The receptors at which the noise constraint is violated are colored red, while the ones at which the constraint is satisfied are blue.It is noticed how L p,j estimated by ISO 9613-2 model is estimated to violate the constraint at the initial mode, m 0,i .When observing the WindSTAR results for m 0,i , two receptors are estimated to experience L p,j values lower than the constraint.Especially one receptor positioned in the upwind of all wind turbines is exposed to a significantly reduced L p,j as a result of upward refraction.It is however noted, that the uncertainty of WindSTAR in the shadow zone might be significantly higher due to the fact that turbulence is omitted.Thus, a higher L p,j could be expected at these positions due to scattering of sound into the shadow zone.In general, it should for all optimization cases be kept in mind, that the optimization are done for a single wind direction and wind speed which in this case causes two receptors to be directly upwind.Thus, normally a small variation in i.e. the wind direction would be expected for each considered flow case.The effect   During the optimization it is observed that the noise reducing modes are distributed to all wind turbines when using the ISO it is observed that the turbines in the outer positions of the row have the highest curtailment while the turbines in the centre are less noise curtailed.This is even though the receptors closest to the centre of the row initially are exposed to the highest L p,j .
Convergence plots are similarly presented for the optimization of the row of 7 SWT-DD-142 turbines in Figure 10 at θ = 0 • .
In general slightly higher L p,j are experienced for all receptors due to the increased turbine size.It is observed that during the ISO 9613-2 optimization, the operational modes are gradually modified causing receptor 1, 2 and 4 to quickly reach L p,j values below the constraint.However, although the L p,j of receptor 3 initially is at the same level as the L p,j at receptor 1 and 2, it is not as significantly reduced with the change of operational modes.Thus, the power output is seen to reach a certain level at which the constraints are satisfied and after this point not being able to optimize any further.The large curtailment in order to bring the L p,j of receptor 3 below 37 dB(A) further causes the remaining receptors to experience an L p,j significantly below the constraint.For the WindSTAR optimization the L p,j at all receptors is on the other hand kept close to the constraint.
Thus, although the L p,j at receptor 3 is just below 37 dB(A), in the further iterations the optimizer still manages to find a more optimal solution that increases the power output.
The scatter plot in Figure 11 emphasizes the observations done in the convergence plots in Figure 10.Thus, it is apparent that a higher L p,j is generally computed at each receptor.The L p,j at receptor 2 which for the SWT-2.3-93turbine type was estimated by the WindSTAR model to be slightly below 37 dB(A), is now violating the noise constraint.Moreover, receptor 1 with significantly lower L p,j in Figure 9, which was expected to be caused by a shadow zone in the upward refracting atmosphere captured by WindSTAR, is observed to receive an L p,j similar to the L p,j estimated at the remaining receptors.
This can be a result of the increased hub height when going from the SWT-2.3-93 to the SWT-DD-142 turbine, yielding an increased height of the source positions.The higher source positions is expected to thereby cause the sound waves to travel further before attenuating due to the upward refraction (Bolin et al. , 2020).Furthermore, the 3 point sources are distributed over a larger area due to the increased rotor diameter.Thus, averaging over points that are further apart could further cause the effect from the shadow zone to be reduced.For the optimization using the ISO 9613-2 model it is observed how 3 of the wind turbines are switched to the most noise reducing mode, m i = 6, and only one wind turbine is kept at m i = 0.The power output is thereby heavily curtailed in order to satisfy the noise constraint.The corresponding optimization with the WindSTAR model results in a higher power output due to the estimated lower L p,j .
The optimization of the row of 7 turbines is further performed for a flow case with a wind direction of θ = 225 • .The resulting flow field in the x/y-plane at hub height is presented in Figure 12 for the SWT-DD-142 turbines.The choice of wind direction results in the majority of the wind turbines being positioned in a wake field.This will cause reduced effective wind speeds, presumably leading to a lower L p,j at the receptors.Furthermore, none of the receptors are positioned directly in the wake or directly upwind of a turbine.Hence, they are all positioned in a free flow field.
Comparing the convergence for the SWT-2.3-93turbine in Figure 13 with the convergence for the SWT-DD-142 turbine in  mized solution.It is noticed how the L p,j at receptors 3 and 4 positioned downstream along the row of turbines is significantly lower when θ = 225 • .This is caused by the reduced effective wind speed, U ef f , of the nearest wind turbines since they are now in a wake position.This reduction is more apparent when observing the WindSTAR optimization in which the L p,j at receptor 3 and 4 is well below the defined noise constraint.Receptor 1 and 2, positioned more upstream of the row of wind turbines, are however still exposed to higher L p,j , which for the ISO 9613-2 optimization results in the observed mode reduction.For the WindSTAR optimization, only the L p,j at receptor 1 is violating the noise constraint leading to only two of the wind turbines operating at noise reducing modes.
A different behaviour of the optimizer is noticed when considering the row of SWT-DD-142 turbines in Figure 14.Thus, it is observed that as the turbine operation is modified to a noise reducing mode, the power output increases.This is caused by the reduced C T at the noise reducing modes leading to a higher U ef f at the rotor positioned in the wake.Thereby, although the front turbine will produce a decreased amount of power, the overall output of the wind farm will be improved.This optimization also yields a significant decrease in the L p,j at each receptor, automatically keeping it below the noise constraint.This tendency of the optimization is only obtained for the larger turbine type, which is expected to be caused by the larger differences in the power-and C T curve observed for the SWT-DD-142 turbine at U hub = 10 m/s in Figure 3 than for the SWT-2.3-93turbine at U hub ≈ 9.3 m/s in Figure 2. Hence, the C T curves for the operational modes defined for the SWT-2.3-93may not be sufficiently large and cause a higher overall power output of the wind farm.uses more iterations before reaching convergence due to the increased number of design variables.As has been discussed for the row of 7 wind turbines, the WindSTAR model generally estimates a lower L p,j at all receptors than the ISO 9613-2 model.
It is further noticed, that even though the size of the wind farm has increased significantly, the L p,j is still similar to the L p,j estimated at the row of turbines.Hence, the noise characteristics of the nearest wind turbines seem to have a larger impact on the received L p,j than the total number of noise sources does.The upwind positions of some of the receptors are further seen to not significantly affect the L p,j estimated by WindSTAR, which is expected to be due to the contribution from the remaining wind turbines nearby.The higher L p,j estimated by the ISO 9613-2 model causes the wind turbines to be generally heavily curtailed.It is observed how at the iteration evaluating the minimum power output, the heavily noise reducing modes are distributed to all turbines in the wind farm.At the end of the optimization, the turbines positioned at the edges of the wind farm are heavier curtailed, while the centre turbines are modified to lower operational modes.Although the L p,j is estimated by WindSTAR to violate the noise constraint at almost all receptors at the initial operational mode m 0,i = 0, the WindSTAR optimization still manages to reach a power output very close to the initial power output.The computed L p,j at receptor 4 positioned upwind of the wind farm is just below the noise constraint of 37 dB(A), resulting in that the upper right row of turbines closest to receptor 4 proceeds to operate at m 0,i .

Further discussion
In the optimizations performed in the presented work, the optimization method using a random search algorithm has been applied (Feng et al. , 2015(Feng et al. , , 2017)).The method is easy to apply to an optimization problem of discrete design variables and by being a global search algorithm, it is more certain to find the global minimum/maximum.However, as briefly mentioned, the method has the disadvantage of getting less efficient as the size of the wind farm in question, the number of receptors or the number of possible operational modes increases.Therefore, it can be beneficial to use a gradient-based optimization method for the defined problem (Martins et al. , 2021).The gradient-based approach requires that the functions in the optimization can be assumed to be continuous such that the gradient of the ISO 9613-2 model and the WindSTAR model, respectively, as well as 435 the overall power output of the wind farm obtained from PyWake with respect to the discrete operational modes can be derived.
By doing so, the optimization method using the WindSTAR model can be applied to larger problems and to a broader range of flow cases in order to obtain an estimate of the optimized AEP.It should however be noted that the optimizations presented in this article are considered a 'proof of concept' of the developed approach to optimize the wind farm operation based on the advanced sound propagation modeling.The use of the random search algorithm is therefore concluded to be a feasible choice 440 for this purpose.This had been further emphasized through this article by the fast convergence of the presented optimization studies.The use of any of the two sound propagation models introduces a certain uncertainty to the predicted sound pressure level at each receptor.First of all, both the WindSTAR model and the ISO 9613-2 model compute the sound propagation based on simplified flow fields using a logarithmic inflow profile and an engineering wake model.Thus, these simplifications introduce uncertainties already in the flow field modeling which is expected to propagate as uncertainties in the sound propagation modeling.It should however be noted that the use of the logarithmic inflow profile is deemed acceptable for the flat terrain in the studied wind farm cases.For a complex terrain wind farm, higher fidelity flow modeling like RANS should be considered in order to obtain the speed ups in the flow field.In addition, the turbulence effects in the atmosphere are neglected due to the high computational costs.This will in some scenarios, i.e. when considering receptors in the upwind position of a wind turbine, lead to higher uncertainties due to the omitted scattering of sound.The turbulence effects in the sound propagation modeling of WindSTAR could be included by i.e. developing a surrogate model based on a limited amount of model evaluations (Martins et al. , 2021).This was however considered out of the scope of the work done in this article.In general the ISO 9613-2 model is observed to estimate higher sound pressure levels at the different receptors compared to the WindSTAR model.Although this suggests that the ISO 9613-2 model is more conservative, it on the other hand gives a higher insurance that the noise constraints are not violated.Thus, the lower estimated sound pressure level of WindSTAR may lead to that the noise constraints in reality are not satisfied at the obtained optimal mode.This could i.e. be accounted for by adding the uncertainty of the WindSTAR model to the integrated sound pressure levels prior to the optimization.However, in general the higher fidelity model gives a better prediction of the noise at each receptor and allows for a broader exploration of the flow parameters and their influence on the L p,j .

Conclusions
Through the work of this article a new approach for performing optimization of wind farm operation was presented.The optimization considers noise constraints at nearby receptors of an onshore wind farm.By the use of the ISO 9613-2 and WindSTAR sound propagation models as well as the Topfarm optimization framework and PyWake flow model the overall power output is optimized in a specific flow case while assuring that the sound pressure level satisfies the given noise constraints.This is done by individually changing the defined operational modes of each wind turbine in the wind farm.The approach was tested on a smaller wind farm of 7 wind turbines and 4 receptors which showed a fast convergence for both sound propagation models and a significant gain in power output when using the WindSTAR model over the ISO 9613-2 model.Especially for cases in which one or more receptors are in the upwind positions of the wind farm, the use of the WindSTAR model in the optimization results in lower estimated sound pressure levels at the receptors and a higher overall power output of the wind farm.While being a more advanced sound propagation model, it is also evident that the use of the WindSTAR model requires longer computational times.It was therefore tested whether the sensitivity of the WindSTAR model to the operational modes is negligible, such that the WindSTAR computations can be performed once prior to the optimization and later used as a transfer function during the iterations in the optimization.It was shown that variations in the sound attenuation are most apparent for far distances where the sound pressure levels are already low.These variations were therefore omitted and WindSTAR was used as a transfer function.As an analysis for future work with the presented framework, the potential and uncertainty in replacing WindSTAR computations at the cases of high frequencies and long distances with ISO 9613-2 computations will be investigated.This is already done for frequencies and distances where the memory of the computations becomes too excessive.However, there is a potential value in implementing this replacement at shorter distances for f = 1 kHz and f = 2 kHz and thereby ideally reducing the computational time even further.Finally, the optimization framework has been tested on an artificial onshore wind farm of the size of 4x5 SWT-DD-142 4.1 MW wind turbines and 4 nearby receptors.Although being a larger wind farm, both sound propagation models show that the sound pressure levels at each receptor do not necessarily increase, implying that the noise characteristics of the nearest wind turbines are of higher importance than the number of turbines in the considered wind farm.
As it has been discussed, the use of a random search algorithm for the optimization does not guarantee a global optimum.

485
In order to fully exploit the capabilities of the framework and to further approach a globally optimal solution, a gradient based approach should be implemented.This requires that the gradient of the sound pressure level at each receptors with respect to the operational modes of each wind turbine is derived.Thus, this is considered the next step in the development of the framework.
https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.The paper is organized as follows.Section 2 shortly describes the WindSTAR and ISO 9613-2 models used for sound propagation as well as the Topfarm framework used for optimization and PyWake used for the wind farm flow modeling.Section 3 presents the optimization problem and the flowchart of the different models in question.Section 4 defines a few selected cases used for the tests of the optimization framework, and further presents a sensitivity study of the C T variation in the WindSTAR model.The results of the performed optimization are presented and discussed in Section 5, while Section 6 lists the conclusions of the work done.

A
gr , are included in the version of the ISO 9613-2 model implemented.The attenuation caused by barriers in the propagation path, A bar , and the attenuation caused by any miscellaneous effects, A misc , are neglected in the presented work.
) https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.where α(f )d is the atmospheric absorption of the sound along the distance d between the source and the receptor.Both the ISO 9613-2 and the WindSTAR model include the atmospheric absorption, α(f )d , in the computations by following the procedure of (

m
from the bottom to the top boundary.The propagation of sound from a wind turbine is normally considered as the sound https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.propagating from a point source positioned at hub height.Thus, a single point source representation is used for the ISO 9613-2 computations.The individual computations in the WindSTAR model are as well assuming a single point source at the specified location.However, in order to represent the wind turbine rotor, 3 point sources are positioned at z = z hub and z = z hub ±85%R https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.cost model, which estimates i.e. the Levelized Cost of Energy (LCoE) and the AEP of the wind farm in question.Topfarm has https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.
operational mode causing a change in the integrated sound pressure level, L p .The propagation of the sound itself is altered by the operational mode through the updated (u ij , v ij , w ij ) flow field parsed from the wind farm wake model to the sound propagation model.This step is only applicable for the WindSTAR model, since the ISO 9613-2 model does not take the entire flow field into account.It is later examined in section 4.2 whether the sensitivity of the L p,j provided by the WindSTAR model to the operational modes is negligible such that the sound propagation model alternatively can be computed only once prior to the optimization and used as a transfer function in each iteration.In this way, the computational time can be reduced considerably.When using the more simple ISO 9613-2 sound propagation model for the optimization, the flow chart reduces.Thus, only the effective wind speed at each wind turbine in the wind farm is parsed from the wind farm flow model to the sound model in order to obtain the L W of the wind turbines.Furthermore, the flow case information and terrain elevation are only parsed as inputs to the wind farm wake model in PyWake, and not to the sound propagation model, since the ISO 9613-2 model does not take these parameters into account.Since the results from the ISO 9613-2 model are not depending on the flow field in the wind farm, the operational mode only changes the L W,i of each wind turbine.Thus, the propagation of the sound provided by the ISO 9613-2 model remains unchanged in each iteration.

Figure 1 .
Figure 1.Flow chart of the optimization framework structure 7) https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.Where n wt is the total number of wind turbines in the wind farms and n rpt is the number of receptors.The objective of the optimization is to maximize the total power output i P i (U 0 , θ, m i ) at a given flow case, where U 0 is the free wind speed at hub height and θ is the wind direction.The design variables are given by the operational modes of each wind turbine in the wind farm, m i , which are subject to a lower and upper bound determined by the design of the wind turbine in question.A set of constraints are given, by which the overall L p at each receptor integrated from each wind turbine must stay under the given limit in dB(A).4Test of optimization4.1 Site and wind turbine typesFor the optimizations done in this article, a wind farm in flat terrain is considered.For the layout of the wind turbines, the Lillgrund wind farm is used as a reference site.Although being an offshore wind farm, Lillgrund provides a flat terrain case consisting of a wind turbine type with various noise reducing operational modes.The Lillgrund wind farm has a size of 48 wind turbines, but only parts of the wind farm have been used for the tests performed in this work.Thus, the tests consider one row of the wind farm consisting of 7 wind turbines and the North-East corner of the wind farm consisting of a layout of 4x5 wind turbines.In order to have dwellings in a near distance of the wind farm at which the noise constraints should be fulfilled, artificial receptors are arbitrarily placed around the wind farm with a distance to the nearest wind turbine no closer than four times the total height of the wind turbine type.The original wind turbines in the Lillgrund wind farm are of the type Siemens SWT-2.3-93.A number of defined operational modes are provided for this wind turbine type with information about the L W spectra and the corresponding power-and C T curves.The operational modes of the SWT-2.3-93are however only given for hub height wind speeds up to 11 m/s.Thus, additional optimizations using the larger Siemens SWT-DD-142 4.1 MW wind turbine are performed as well.The operational modes of the SWT-DD-142 turbine are defined for hub height wind speeds up to 20 m/s which allows for more exploration of the sensitivity of the L p , while the higher L W values of the larger turbine type introduce a larger need for optimization.The specifications of the two wind turbines are listed in Table1.The power-and C T curves as well as the L W spectra at U hub = 10 m/s are shown in Figure2for the SWT-2.3-93turbine and in Figure3for the SWT-DD-142 turbine.For both turbine types a reduction in the L W at each octave band frequency as well as the corresponding power and C T is observed during the discrete steps from operational mode 0 to 6.It is further observed that the operational modes of the two turbine types introduce a similar reduction in L W while the power reduction is more distinct for the SWT-DD-142.The larger rotor diameter of the SWT-DD-142 requires a rescaling of the distance between the turbines in the chosen layouts.In the original layout using the SWT-2.3-93, the distance between turbines is 4.3D in the direction from South-West to North-East and 3.3D in the direction from North-West to South-East.These distances are therefore used to scale the wind farm layout to fit the rotor diameter of D = 142 m of the SWT-DD-142 turbine.

Table 1 .Figure 2 .
Figure 2. SWT-2.3-93:The operational mode depending power-and CT curves and the sound power level, LW , curve at U hub = 10 m/s.From dark to light blue: the least noise reducing mode, m = 0 to the most noise reducing mode, m = 6.

Figure 3 .
Figure 3. SWT-DD-142: The operational mode depending power-and CT curves and the sound power level, LW , curve at U hub = 10 m/s.From dark to light blue: the least noise reducing mode, m = 0 to the most noise reducing mode, m = 6.

Figure 4 .
Figure 4. Left: Setup 1 with the downwind attenuation of sound from a single wind turbine at different operational modes.Right: Setup 2 with the downwind attenuation of sound from a single wind turbine in the wake of wind turbine at different operational modes.The distance between the two wind turbines is d ≈ 3.8D.

Figure 5 .
Figure 5. Setup 1: CT sensitivity of WindSTAR obtained transmission loss, T L, from a single wind turbine subject to changing operational modes at each octave band frequency at a free field wind speed of U hub = 10 m/s at hub height.

Figure 6 .
Figure 6.Setup 2: CT sensitivity of WindSTAR obtained transmission loss, T L, from a single wind turbine positioned in the wake of a wind turbine subject to changing operational modes at each octave band frequency at a free field wind speed of U hub = 10 m/s at hub height.
cases and constraintsFor the initial tests of the optimization, the row of 7 wind turbines either of the type SWT-2.3-93 or SWT-DD-142 are used.As mentioned, the distances between the turbines are scaled to fit the rotor diameter of the turbine type.As mentioned, four https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.
4 and the friction velocity is set to u * = 0.57 m/s.The choice of parameters yield a hub height wind speed of U hub = 10 m/s for the SWT-DD-142 turbine with z hub = 109 m and U hub = 9.3 m/s for the SWT-2.3-93turbine with z hub = 68.5 m.The ground flow resistivity used for computing the acoustic impedance is kept at σ = 2 • 10 4 k Pa s m −2 (Wagner et al. , 1996) in the WindSTAR model, while the ground factor in the ISO 9613-2 model is kept at G = 0 (ISO 9613-2 , 1997).The chosen values for the ground parameters in both models are representative to a hard, reflecting ground.In the final optimization of the presented work, the described layout of 4x5 wind turbines is considered.In this layout the larger SWT-DD-142 turbine is used.In a similar way as for the row of 7 wind turbines, 4 receptors at different arbitrarily chosen positions near the wind farm are considered.The chosen layout results in 4 • 20 = 80 individual WindSTAR computations.For the ISO 9613-2 model the 80 individual computations are performed in every iteration of the optimization.The noise constraint 3-93 turbines and a row of 7 SWT-DD-142 turbines.The flow fields through the row of wind turbines for the chosen flow cases with the two different wind directions θ = 0 • and θ = 225 • are observed for the SWT-DD-142 wind turbine in Figures 7 and 12. Figures of the flow fields for the SWT-2.3-93turbines are not included here, since the layout is very similar to the row of SWT-DD-142 turbines.It is noticed, that the four receptors will be in either the upwind, downwind or crosswind positions of the wind turbines depending https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.on the wind direction.For θ = 0 • , two receptors are positioned directly in the wake of a wind turbine, while the two remaining receptors are positioned in the free flow field directly upwind of a wind turbine.For θ = 225 • all receptors are positioned in the free flow field either dominantly downwind or upwind of the wind turbines.

Figure 7 .
Figure 7.The flow field at hub height through the row of 7 SWT-DD-142 wind turbines at the wind direction θ = 0 • .The distances between the turbines are scaled to be approximately 4.3D.

Figure 8 .
Figure 8.The operational mode, mi, for each wind turbine (left), the overall power, P , of the row of 7 SWT-2.3-93wind turbines (middle) and the integrated sound pressure level, Lp,j, at each receptor (right) during the optimization at θ = 0 • and U hub = 9.3 m/s.The noise constraint is set to L p,lim,j = 37 dB(A) and represented by the dashed line in the right figure.

345
of varying the wind direction in the WindSTAR model is shortly investigated for θ = 0 • ± 15 • in Table2.Observing receptor 3 and 4, the L p,j is further reduced when considering θ = −15 • since the positions of the receptors become even more upwind relative to the closest wind turbine.At θ = +15 • the position of receptor 4 relative to the closest turbine is more crosswind resulting in a distinct increase in L p,j .On the other hand, receptors 1 and 2 experience smaller variations in the L p,j with the https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License.change in θ.Thus, the more abrupt changes in sound propagation appear when receptors are in the upwind position.It should be noted that these variations only appear in the WindSTAR results since the propagation of the ISO 9613-2 model is not sensitive to changes in the wind direction when the wind turbines are positioned in the free field.

Figure 9 .
Figure 9. Scatter plots of the operational modes of the 7 SWT-2.3-93wind turbines before (left), at the minimum power iteration during (middle) and after (right) the optimization at θ = 0 • and U hub = 9.3 m/s.The noise constraint is set to L p,lim,j = 37 dB(A).
9613-2 model, while the WindSTAR optimization only modifies 2-3 turbines.For the optimum, all wind turbines are switched https://doi.org/10.5194/wes-2022-80Preprint.Discussion started: 1 September 2022 c Author(s) 2022.CC BY 4.0 License. to a noise reducing mode in the ISO 9613-2 optimization while the operation of only two turbines, close to the receptors initially being subject to constraint violations, is modified in the WindSTAR optimization.For the ISO 9613-2 optimal mode

Figure 14 ,
Figure14, both at θ = 225 • , some noticeable differences can be observed.Similar to the flow case of θ = 0 • , the optimization using the SWT-2.3-93turbine and the ISO 9613-2 model yields reduced modes for a majority of the wind turbines in the opti-

Figure 10 .
Figure 10.The operational mode, mi, for each wind turbine (left), the overall power of the 7 SWT-DD-142 wind turbines (middle) and the integrated sound pressure level, Lp,j, at each receptor (right) during the optimization at θ = 0 • and U hub = 10 m/s.The noise constraint is set to L p,lim,j = 37 dB(A) and represented by the dashed line in the right figure.

Figure 11 .
Figure 11.Scatter plots of the operational modes of the 7 SWT-DD-142 wind turbines before (left), at the minimum power iteration during (middle) and after (right) the optimization at θ = 0 • and U0 = 10 m/s.The noise constraint is set to L p,lim,j = 37 dB(A).
larger wind farm layout of 4x5 wind turbines is tested by using the SWT-DD-142 turbine.The larger wind turbine type is chosen due to the more distinct differences in the defined operational modes seen in Figure3.The flow field including the wind turbine wakes through the wind farm are shown in Figure 15 with θ = 0 • .Similar to the row of 7 wind turbines analysed thus far, 4 receptors are arbitrarily positioned around the wind farm.As can be noticed in Figure 15, the receptors are positioned downwind, upwind or crosswind of the wind turbines.The convergence of the ISO 9613-2 and WindSTAR optimization, respectively, of the 4x5 wind farm is presented in Figure 16 and scatter plots are shown in Figure 17.It is observed that for both the ISO 9613-2 and the WindSTAR model, the optimizer

Figure 12 .
Figure 12.The flow field at hub height through the row of 7 SWT-DD-142 wind turbines at the wind direction θ = 225 • .The distances between the turbines are scaled to be approximately 4.3D.

Figure 13 .
Figure 13.The operational mode, mi, for each wind turbine (left), the overall power of the 7 SWT-2.3-93wind turbines (middle) and the integrated sound pressure level, Lp, at each receptor (right) during the optimization at θ = 225 • and U hub = 9.3 m/s.The noise constraint is set to L p,lim = 37 dB(A) and represented by the dashed line in the right figure.

Figure 14 .
Figure 14.The operational mode, mi, for each wind turbine (left), the overall power of the 7 SWT-DD-142 wind turbines (middle) and the integrated sound pressure level, U hub = 10 m/s, at each receptor (right) during the optimization at θ = 225 • and U0 = 10 m/s.The noise constraint is set to L p,lim,j = 37 dB(A) and represented by the dashed line in the right figure.

Figure 15 .
Figure 15.Example of flow at hub height through the 4x5 SWT-DD-142 wind turbines at the wind direction θ = 0 • .The distances between the turbines are scaled to be approximately 4.3D and 3.3D.

Figure 16 .
Figure 16.The operational mode, mi, for each wind turbine (left), the overall power of the 4x5 SWT-DD-142 wind turbines (middle) and the integrated sound pressure level, Lp,j, at each receptor (right) during the optimization at θ = 225 • and U hub = 10 m/s.The noise constraint is set to L p,lim,j = 37 dB(A) and represented by the dashed line in the right figure.

Figure 17 .
Figure 17.Scatter plots of the operational modes of the 4x5 SWT-DD-142 wind turbines before (left), at the minimum power iteration during (middle) and after (right) the optimization at θ = 0 • and U0 = 10 m/s.The noise constraint is set to L p,lim = 37 dB(A).

Figure A2 .
Figure A2.Setup 1: CT sensitivity of WindSTAR obtained transmission loss, T L, from a single wind turbine subject to changing operational modes at each octave band frequency at a free field wind speed of U0 = 6 m/s.

Figure A3 .
Figure A3.Setup 2: CT sensitivity of WindSTAR obtained transmission loss, T L, from a single wind turbine positioned in the wake of a wind turbine subject to changing operational modes at each octave band frequency at a free field wind speed of U0 = 14 m/s.

Figure A4 .
Figure A4.Setup 2: CT sensitivity of WindSTAR obtained transmission loss, T L, from a single wind turbine positioned in the wake of a wind turbine subject to changing operational modes at each octave band frequency at a free field wind speed of U0 = 14 m/s.

Table 2 .
Overall Lp,j estimated by WindSTAR at each receptor for θ = 0 • ± 15 • for the row of 7 SWT-2.3-93wind turbines operating at