Impact of wind profiles on ground-generation airborne wind energy system performance

This study investigates the performance of pumping-mode ground-generation airborne wind energy systems (AWESs) by determining cyclical, feasible, power-optimal flight trajectories based on realistic vertical wind velocity profiles. These 10 min profiles, derived from mesoscale weather simulations at an offshore and an onshore site in Europe, are incorporated into an optimal control model that maximizes average cycle power by optimizing the trajectory. To reduce the computational cost, representative wind conditions are determined based on k-means clustering. The results describe the influence of wind speed magnitude and profile shape on the power, tether tension, tether reeling speed, and kite trajectory during a pumping cycle. The effect of mesoscale-simulated wind profiles on power curves is illustrated by comparing them to logarithmic wind profiles. Offshore, the results are in good agreement, while onshore power curves differ due to more frequent non-monotonic wind conditions. Results are references against a simplified quasi-steady-state model and wind turbine model. This study investigates how power curves based on mesoscale-simulated wind profiles are affected by the choice of reference height. Our data show that optimal operating heights are generally below 400 m with most AWESs operating at around 200 m.

1 Introduction

Airborne wind energy systems (AWESs) aspire to harvest stronger and less turbulent winds at mid-altitude, here defined as heights between 100 and 1000 m, presumably beyond what is achievable with conventional wind turbines (WTs). The prospects of higher energy yield combined with reduced capital cost motivate the development of this new class of renewable energy technology (Lunney et al., 2017; Fagiano and Milanese, 2012). Unlike conventional WTs, which in recent decades have converged to a single concept with three blades and a conical tower, several different AWES concepts and designs are being investigated by numerous companies and research institutes (Cherubini et al., 2015; Vermillion et al., 2021; Fagiano et al., 2022). These kite-inspired systems consist of three main components: one or more tethered aircraft or kites, one or more ground stations, and one or more tethers to connect the flying components to the ground. This study focuses on the two-phase cyclic, ground-generation concept, also referred to as pumping mode (Luchsinger, 2013). During the reel-out or production phase, the kite pulls a tether from a drum on the ground, which is connected to a generator, thereby producing electricity. This is followed by the reel-in phase during which the kite returns to its initial position and reduces its aerodynamic forces in order to de-power. There are several ways to reduce the aerodynamic forces on a kite, such as adjusting its angle of attack or flying it out of the wind window. Various other concepts such as fly-gen, aerostat, and rotary lift are not considered in this study (Cherubini et al., 2015). Since this technology is still in a relatively early stage of development, validation and comparison of power production estimates is difficult. Several studies (van der Vlugt et al., 2013; Schelbergen et al., 2020; Vander Lind, 2013; Ranneberg et al., 2018) compared computed power curves with experimental performance data. At present, there is no standardized power curve definition or reference design that would allow for comparison between different concepts and conventional wind turbines. However, the goal of Eijkelhof and Schmehl (2022) and Eijkelhof et al. (2020) was to create a reference design for a multi-megawatt AWES. It is not the goal of this study to determine a general power curve but rather to investigate the power variation in a specific design derived from realistic wind profiles.

Recent consensus among the scientific community defined a power curve as the maximum average cycle power, which is the combination of consecutive reel-out and reel-in phases, as a function of wind speed at pattern height, which is the time-averaged height during the reel-out power-producing phase (Airborne Wind Europe, 2021). Together with the site-specific wind resource, wind park planners and manufacturers can use power curves to estimate annual energy production (AEP), which can be combined with a cost model to determine the levelized cost of electricity (LCOE) and financial viability (Malz et al., 2020a). Unlike conventional WTs, where the wind speed probability distribution at hub height is used to determine AEP, AWES continuously change their operating height, making it difficult to determine AEP with this approach. AWES performance highly depends on the shape and magnitude of the wind speed profile over the operating height range. Using simple wind profile approximations, such as logarithmic or exponential wind speed profiles, can provide an estimate of long-term average conditions. However, these approximations cannot account for the wide range of profile shapes that occur over short periods of time or changes that occur on a daily or seasonal basis (Emeis, 2013). This can reduce the accuracy of the predicted power output. However, such wind speed profiles can be employed to estimate average performance and are the standard in most AWES power estimation studies.

In their study, van der Vlugt et al. (2013) described TU Delft's 20 kW inflatable wing technology demonstrator and compared a statistically derived power curve to results from a theoretical performance model. The wind and power models used in this study were taken from Fechner and Schmehl (2013). The wind model is based on a standard exponential wind speed profile approximation, while the power model uses a multi-phase quantitative susceptibility mapping (QSM). A follow-up study (van der Vlugt et al., 2019) added more detail to specific cycle trajectories. Heilmann and Houle (2013) used exponential wind speed profiles and a standard Rayleigh distribution to estimate performance and cost. Their power curve is modeled based on a QSM by Luchsinger (2013) using the averaged flight path height of the kite as the wind speed reference height. An LCOE between EUR 40 and 110 per megawatt hour was estimated for different annual average wind speeds. Faggiani and Schmehl (2018) used a similar pumping-mode QSM to estimate power output based on wind speed at the operating height of the kite. They developed a cost model to estimate the achievable LCOE of an entire kite wind farm. Their analyses showed that the cost of energy decreases and the quality of the electrical power increases with increasing number of kites. Ranneberg et al. (2018) studied the performance of a soft-kite pumping-mode AWES by determining its power curves at various reference heights for different logarithmic wind profiles. The study found that the yield variation for these logarithmic wind profiles was quite small. Additionally, the yield for a specific site was estimated using detailed wind speed profiles from COSMO-DE, and the results were found to be consistent with more detailed simulations of the EnerKíte EK30 prototype (EnerKíte GmbH, 2022). Leuthold et al. (2018) investigated power-optimal trajectories of a ground-generation multikite configuration for a range of logarithmic wind speed profiles. Three distinct operational regions were identified: Region I where power is used to maintain altitude, Region II where power harvesting increases up to the designed wind speed, and Region III where power extraction is intentionally limited due to the physical constraints of the system. Licitra et al. (2019) estimated the performance and power curve of a fixed-kite ground-generation AWES by generating power-optimal trajectories using a power law approximation of the wind speed profile. The results were validated against data from Ampyx Power AP2 (Licitra, 2018; Malz et al., 2019; Ampyx Power BV, 2020). Eijkelhof and Schmehl (2022) found that mass had a detrimental effect on power-optimal trajectories for a large-scale single-kite fixed-wing AWES. To determine power curves, the authors used normalized average offshore wind speed profiles from the Ijmuiden measurement tower. Sommerfeld et al. (2022) used the same methodology and wind data as in this study to examine the effects of size scaling and improvements in aerodynamic efficiency on a single-kite fixed-wing reference system. The authors found that it is likely better to deploy multiple smaller-scale devices rather than a single large-scale system because of negative mass-scaling effects. De Schutter et al. (2019) analyzed the performance of utility-scale, stacked multikite systems using logarithmic wind speed profiles as boundary conditions for a nonlinear optimization problem. The authors used the same optimization framework as in the present investigation. They found that this multikite strategy could make power generation largely independent of wing size. Malz et al. (2020b) efficiently estimated the performance of a single-kite drag-mode AWES for large wind data sets by combining an optimal control performance model with smart initialization and machine learning. Wind speed profiles from the MERRA-2 reanalysis model (Gelaro et al., 2017) were clustered into characteristic profile shapes and interpolated using Lagrangian polynomials (Sect. 4.4). The authors showed that ordering the wind parameters by the wind speed at the average operational height (300 m) significantly reduced the computation time. Aull et al. (2020) explored the design and sizing of fly-gen rigid-kite systems based on a steady-state model with simple aerodynamic and mass-scaling approximations. At each scale, the relationships between size, efficiency, power output, and cost were determined. The wind resource was described by an exponential wind shear model with a Weibull distribution. The main conclusion was that physics and economics favor a larger number of small units. Bechtle et al. (2019) used ERA5 reanalysis data to assess wind resources at high altitudes throughout Europe. They described the available wind energy without considering a specific conversion mechanism and included a description of wind speed and probability at various heights. The effect of variable height harvesting was demonstrated for a location in the English Channel. Schelbergen et al. (2020) proposed a clustering procedure to obtain wind statistics from the Dutch Offshore Wind Atlas (DOWA) data set. Principal component analysis and k-means clustering were used to determine representative wind speed profile shapes. To estimate the AEP of a small-scale pumping AWES located at Cabauw in the center of the Netherlands, several power curves were derived for each wind speed profile shape using a flexible-kite, pumping-mode QSM developed by van der Vlugt et al. (2019). Faggiani and Schmehl (2018) studied the economic impact of various design aspects of wind parks, including the spatial stacking of systems, the number of units, the size of kites, and phase-shifted operation. The performance of the system was estimated using a QSM developed by Schmehl et al. (2013) and van der Vlugt et al. (2019), assuming a range of wind speeds at the operating height of the kite. The AEP LCOE could be assessed by combining a detailed cost model with an assumed Weibull probability distribution. The study found that increasing the number of kites had several scale effects, such as decreasing the cost of energy and increasing the quality of electrical power. Wind speed profiles are governed by weather phenomena and environmental and location-dependent conditions on a multitude of temporal and spatial scales. The preferred means of determining wind conditions for wind energy converters are long-term, high-resolution measurements. At mid-altitudes, these measurements can only be obtained through long-range remote sensing methods such as lidar (light detection and ranging) or sodar (sonic detection and ranging). Measuring wind conditions at mid-altitudes is costly and difficult due to reduced data availability (Sommerfeld et al., 2019a). Additionally, publicly available measurements are scarce. Therefore, wind data in this study are derived from Weather Research and Forecasting (WRF) model mesoscale simulations (Skamarock et al., 2008). However, the described trajectory optimization methodology can be applied to any wind data set such as wind atlas data or measurements. Numerical mesoscale weather prediction models such as WRF, which is well known for conventional WT siting applications (Salvação and Guedes Soares, 2018; Dörenkämper et al., 2020), are used to estimate wind conditions on timescales of a few minutes to years. Sommerfeld et al. (2019b) compared the simulated onshore data used in this study, located in northern Germany near the city of Pritzwalk, to lidar measurements and found a good agreement between both data sets. Data from the FINO3 research platform in the North Sea can be used as a reference for the simulated offshore conditions in this study. The present study investigates the performance of AWES subject to 10 min average wind data, which is the standard for conventional WTs, while the New European Wind Atlas (NEWA) only provides 30 min average data (Witha et al., 2019). We use these higher-resolution wind data because the higher temporal, spatial, and vertical resolution reduces averaging and allows for the investigation of more realistic wind conditions.

This paper's main contribution is the examination of how realistic onshore and offshore wind profiles, compared to a standard log profile, affect the power-optimal performance of AWESs, as well as how the choice of reference height impacts the power curve, particularly given the wide range of wind speed profile shapes. This study is a continuation of previous analyses of lidar measurements (Sommerfeld et al., 2019a) and WRF simulations (Sommerfeld et al., 2019b) at the onshore location. To demonstrate the validity and applicability of the data, several wind characteristics are described. These include annual wind speed probability distributions up to 1000 m, annual wind direction statistics, and wind speed profile shapes. The data are categorized using k-means clustering (Lloyd, 1982; Hartigan and Wong, 1979) which classifies the wind data at each location into groups of similar wind speeds and vertical profile shapes. Three representative 10 min wind velocity vectors U are sampled from each of these 20 clusters (total 60 vectors out of 52 560 wind data points per location) and serve as boundary conditions for the AWEBox trajectory optimization (De Schutter et al., 2020). The profiles with the 5th, 50th, and 95th percentiles of wind speed at an a priori guess of the pattern trajectory height $z_{\mathrm{pth}}\approx \mathrm{100}\mathrm{m}\le z\le \mathrm{400}\mathrm{m}$ (Airborne Wind Europe, 2021) are selected because they encompass the most probable wind conditions within each cluster while excluding non-representative extremes. Section 6.3 verifies that choice and compares the impact of reference height on the power curve. This drastically reduces the computational cost as only a few selected profiles are used to represent the entire spectrum with sufficient resolution to investigate the variation in power.

The AWEBox optimization model allows for the investigation of dynamic performance parameters, such as aircraft trajectories, tether tension, tether reeling speed, and power which highly depend on the wind conditions. The aircraft model is based on the well-investigated and published Ampyx Power AP2 prototype (Licitra, 2018; Malz et al., 2019; Ampyx Power BV, 2020), adjusted to a projected surface area of A=20 m² to generate results for more realistic and probable devices. The maximized power curves, estimated based on average cycle powers and wind speed at reference height ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}=U\left(z_{\mathrm{ref}}\right)$, are compared to performance predictions using a simple AWES QSM and a steady-state WT model. The variation in performance caused by realistic wind data is referenced against predictions based on simple logarithmic wind speed profiles. The a priori guess of $\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m}$ is confirmed by comparing the impact of different reference heights on the power curves.

The structure of this research is as follows. Section 2 introduces the mesoscale WRF model, analyzes the offshore and onshore wind resource, introduces the k-means clustering algorithm, and summarizes the results of the clustered wind vectors. Section 3 introduces the dynamic AWES model, which includes an aircraft and tether model, as well as ground station constraints. Section 4 describes the AWEBox optimization toolbox and summarizes the aircraft parameters, system constraints, and initial conditions. This is followed by a description of the WT and AWES reference models in Sect. 5. The results presented in Sect. 6 include flight trajectories and time series data for various performance parameters, as well as a statistical analysis of the tether length, operating altitude, and power curve estimations. Section 7 summarizes the findings and concludes with an outlook and motivation for future work.

2 Wind conditions

Section 2.1 introduces the model and setup of the onshore and offshore mesoscale WRF simulations. Section 2.2 analyzes wind statistics to give an insight into the wind regime at both locations. Clustering, which is introduced in Sect. 2.3, is used to identify groups of similar vertical wind profiles and to select representative profiles from these groups. This significantly reduces the computational cost as only a few selected profiles are necessary to represent the wind regime. Section 2.4 and 2.5 describe the resulting clusters and their statistical correlation with temporal and meteorological phenomena.

2.1 Mesoscale simulations

This study compares the performance of airborne wind energy systems (AWESs) at two locations in Europe, one onshore and one offshore (Fig. 1). The wind conditions for the chosen years are assumed to be representative of wind conditions at these locations. However, the wind data have not been compared to long-term wind atlas data and have not been corrected using long-term simulations. The onshore data represent wind conditions at the Pritzwalk Sommersberg airport (lat: 53^∘10^′47.00^′′ N, long: 12^∘11^′20.98^′′ E) in northern Germany and comprise 12 months of WRF simulation data between September 2015 and September 2016. The area surrounding the airport consists mainly of flat agricultural land with the town of Pritzwalk in the south, making it a suitable location for wind energy generation (Sommerfeld et al., 2019a, b). The FINO3 research platform in the North Sea (lat: 55^∘11.7^′ N, long: 7^∘9.5^′ E) was chosen as a representative offshore location due to its proximity to several offshore wind farms and the number of comprehensive reference measurements (Peña et al., 2015). The offshore simulation covers the time frame between September 2013 and September 2014.

Figure 1 Map of northern Germany with the representative onshore (Pritzwalk) and offshore (FINO3) locations highlighted by black dots.

The mesoscale simulations used the Weather Research and Forecasting (WRF) model (Powers et al., 2023). The onshore simulation was performed with version 3.6.1 (Skamarock et al., 2008) prior to the 2018 release of WRF version 4.0.2 (Skamarock et al., 2021), in which offshore simulations were computed. The setup of the model was adapted and constantly improved for wind energy applications by the authors of the present article for various projects and applications in recent years (Dörenkämper et al., 2015, 2017, 2020; Hahmann et al., 2020; Sommerfeld et al., 2019b). The focus of this study is not on the detailed comparison between mesoscale models but on AWES performance subject to realistic onshore and offshore wind conditions. Both WRF models provide adequate temporal and spatial resolution for preliminary performance assessment, despite the differences in setup and time frame.

Table 1 Key setup parameters of the onshore and offshore mesoscale WRF simulations to generate the wind data used in this study.

Each simulation consists of three nested domains around their respective location (black dot in Fig. 1). The innermost domain (D03) has the finest resolution and is nested within the middle domain (D02), which is nested within the outermost domain (D01) with the coarsest resolution. The simulations use one-way nesting where the outer domains define the boundary conditions for the inner domains. Atmospheric boundary conditions are defined by ERA-Interim (Dee et al., 2011) for the onshore location and by ERA5 (Hersbach and Dick, 2016) reanalysis data for the offshore location. Sea surface parameters for the offshore location, such as sea surface temperature and sea ice analysis, are based on OSTIA (Donlon et al., 2012). These data sets have proven to provide good results for wind-energy-relevant heights and sites (Olauson, 2018; Hahmann et al., 2020). Four-dimensional data assimilation (FDDA), also known as analysis nudging, nudges the simulation of the outer domain towards the atmospheric boundary conditions throughout the simulation time to reduce numerical drifting and provide smoother boundary conditions. Both simulations use the Mellor–Yamada–Nakanishi–Niino (MYNN) 2.5 (Nakanishi and Niino, 2009) level scheme for the planetary boundary layer (PBL) physics. The onshore simulation was carried out in a single 12-month simulation run from 1 September 2015 to 31 August 2016. The offshore simulation period covered a 410 d period from 30 August 2013 to 14 October 2014 and was divided into 41 simulations of 10 d each, with an additional 24 h of spin-up time per run. Spin-up is the period in which the model produces results that may not be reliable due to initialization using coarser global atmospheric reanalysis data. WRF calculates the vertical coordinate using a hybrid hydrostatic pressure coordinate, which is a function of surface and atmospheric pressure (Skamarock et al., 2021). Data at each vertical, terrain-following pressure coordinate (sigma level) are converted to geometric heights using the postprocessing methodology described by Dörenkämper et al. (2020). Table 1 summarizes the key model settings used in this study. All simulations were performed on the EDDY high-performance computing clusters at the University of Oldenburg (Carl von Ossietzky Universität Oldenburg, 2018).

2.2 Wind regime

Figure 2 depicts the wind roses of the computed annual wind conditions at 100 m (Fig. 2a and b) and 500 m (Fig. 2c and d) height onshore (left) and offshore (right). The dominant wind direction at both locations is southwest, turning clockwise with increasing altitude.

Figure 2 Wind roses of annual wind direction and speed statistics at Pritzwalk (onshore) and FINO3 (offshore) for heights of 100 and 500 m during the simulated years.

Directional variability decreases, and the wind speed U, which is the magnitude of the wind vector U, increases with height, following the expected trends in the Northern Hemisphere (Arya and Holton, 2001; Stull, 1988). The average onshore wind direction turns about 14^∘ between 100 and 500 m, whereas average offshore wind direction only veers approximately 5^∘. The offshore wind direction turns approximately 10 additional degrees above 500 m, resulting in roughly the same westerly wind direction at around 1000 m. Due to the prevailing unstable conditions offshore accompanied by strong vertical mixing, the investigated heights show less veer than onshore. Wind shear at the offshore location is lower compared to the onshore location due to lower surface roughness.

Figure 3 shows the annual horizontal wind speed probability distributions at each height level for both locations. These distributions give insight into the wind statistics at specific heights but not into the statistics of the wind profile shapes. The nonlinear color gradient allows for the representation of the entire relative probability range. Onshore (Fig. 3a) U are relatively low and have a fairly narrow deviation below 300 m due to dominant surface effects. Above this height the distribution broadens, but a high probability of low wind speeds remains for the entire height range. The distributions show bimodal characteristics caused by different atmospheric stratification. Low wind speeds are commonly associated with unstable and high wind speeds with neutral or stable atmospheric conditions.

Figure 3 Comparison of WRF-simulated annual wind speed U probability distributions at each height level between onshore (a) and offshore (b) up to 1000 m. The nonlinear color scheme represents the high probability of low altitude, particularly onshore, while still differentiating the lower, wide-spread frequencies at higher altitudes.

Such multimodal distributions at higher altitudes are better described by the sum of two or more probability distributions, as standard Weibull or Rayleigh distributions cannot capture this phenomenon (Sommerfeld et al., 2019a). Offshore U (Fig. 3b) display a wider distribution at all heights, as they are less affected by surface effects. Similarly to onshore, the offshore frequency distribution also shows a high probability of lower U (between 5–10 m s⁻¹) at all heights. Higher U at lower altitudes benefit conventional WT and weaken the argument for offshore AWESs, as one of their benefits would be to harness energy from the stronger winds at higher altitudes. However, other reasons for placing AWES offshore are the safety and land use regulations and the potential cost benefits of a smaller support structure (offshorewind.biz, 2018; Lunney et al., 2017; Ellis and Ferraro, 2016).

The Obukhov length ℒ (Obukhov, 1971; Sempreviva and Gryning, 1996),

$$\begin{array}{}\text{(1)}& \mathcal{L}={\displaystyle \frac{-u_{*}^{\mathrm{3}}{\mathit{\theta }}_{\mathrm{v}}}{kg}}\left({\displaystyle \frac{\mathrm{1}}{Q_{\mathrm{S}}}}+{\displaystyle \frac{\mathrm{0.61}}{Q_{\mathrm{L}}\mathit{\theta }}}\right),\end{array}$$

commonly characterizes the near-surface atmospheric stability, which highly affects the shape of the wind speed profile U, which is the magnitude of the wind velocity profile U. We extend the concept to mid-altitudes between 100 and 1000 m. The Obukhov length is a function of the simulated friction velocity u_*, the virtual potential temperature θ_v, the potential temperature θ, the kinematic virtual sensible surface heat flux Q_S, the kinematic virtual latent heat flux Q_L, the von Kármán constant k, and the gravitational acceleration g. Various stability classifications based on the Obukhov length have been defined for different wind energy sites. Table 2 summarizes the Obukhov length bin widths used by Floors et al. (2011) and the frequency of occurrence of each stability class onshore and offshore, consistent with Sommerfeld et al. (2019b).

Table 2 Stability classes based on Obukhov length ℒ (bins from Floors et al., 2011) and associated annual probability at Pritzwalk (onshore) and FINO3 (offshore), based on WRF simulations.

Neutral stratification occurs approximately 20 % of the year at both locations. The lower heat capacity of the land surface leads to a faster heat transfer and a quicker surface cool-off which favors the development of stable stratification (≈17 % onshore vs. ≈6 % offshore). The offshore location has a higher probability of unstable conditions, which is likely caused by a warmer ocean surface compared to the air above (Archer et al., 2016).

Figure 4 Representative wind velocity profile U^′ in its original direction and rotated so that the primary wind direction of the resulting profile U points in the x direction. Panel (a) shows the top view, panel (b) shows the front view, panel (c) shows the side view, and panel (d) shows the isometric view. The primary wind direction is defined by the average wind vectors $\stackrel{\mathrm{‾}}{U^{\prime }}\left(z_{\mathrm{ref}}\right)$ (black) and $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ (blue) between 100 and 400 m (dashed lines). As an example, the wind vector at z=150 m (red) (a, d) is decomposed into u and v (b, c).

Both unstable and stable conditions can lead to non-logarithmic and non-monotonic U profiles. Unstable conditions are often accompanied by almost uniform U profiles due to increased mixing, whereas low-level jets (LLJs) can develop during the nocturnal stable onshore boundary layer (Banta, 2008). Both locations have a high chance of unassigned conditions (labeled as “other”) which are mostly associated with low wind speeds.

2.3 Clustering of wind conditions

An accepted method to describe the near-surface atmosphere is atmospheric stability, commonly quantified by the Obukhov length (Obukhov, 1971; Sempreviva and Gryning, 1996), which exclusively uses surface data (Sect. 2.2 and Eq. 1). Previous studies (Sommerfeld et al., 2019a, b) showed that Obukhov-length-classified wind speed profiles U diverge with height, especially during neutral and stable conditions. This indicates vertically heterogeneous atmospheric stability and suggests that surface-based stability categorization is insufficient for higher altitudes. Unlike classifying wind regimes by atmospheric stability, which requires additional temperature and heat flux data, clustering only uses wind data at multiple heights to group profiles by similarity. This results in more cohesive profile groups (Schelbergen et al., 2020). Therefore, clustering can also be applied to wind-only measurements such as lidar.

The k-means clustering algorithm (Pedregosa et al., 2011) used in this study is chosen for its ease of use and scalability due to the high dimensionality of the data set. Many other algorithms produce similar results, but a comparison between clustering algorithms is beyond the scope of this research. The wind velocity profiles U^′ (Fig. 4, black) are rotated such that the main wind component $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)=\stackrel{\mathrm{‾}}{U}(\mathrm{100}\mathrm{m}\le z\le \mathrm{400}\mathrm{m})$ points in the positive x direction U (blue) in order to remove directional dependencies. The velocity components at each height level are decomposed into u in the main wind direction and into v perpendicular to it (red). The wind speed profile U is not shown. This results in more homogeneous clusters and simplifies the comparison of wind data and AWEBox results. It is analogous to assuming omnidirectional operation.

Figure 5 The k-means clustering inertia over of number of clusters k (a) for 1 year of onshore (blue) and offshore (orange) wind velocity profiles U up to 1000 m. The onshore (b) and offshore (d) silhouette coefficients express the distance to neighboring clusters and are color coded according to average $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ between 100 and 400 m, the same as in Figs. 6, 9, 10, and 11. The dashed red line represents the silhouette score, which is the average silhouette coefficient. Silhouette score (c) over number of clusters k for both locations. The number of clusters k=10 has been chosen for presentation purposes only. The analyses in the “Results and discussion” section use k=20 clusters.

Figure 6 Onshore (a, c) and offshore (b, d) average annual wind speed profiles (centroids) resulting from k-means clustering for k=10 (a, b). All WRF-simulated wind speed profiles U are depicted in gray. The centroids are sorted, labeled, and color coded in ascending order of average $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ between heights of 100 and 400 m. The corresponding cluster frequency of occurrence f for each cluster C is shown in (c) and (d) below.

The wind velocity data up to 1000 m are comprised of data points at 30 height levels and in two directions. The clustering algorithm assigns each data point to one of the k clusters represented by their respective cluster mean, also called the centroid. These centroids are chosen to minimize the sum of the Euclidean distances to every data point within each cluster. This cost function is also referred to as inertia or within-cluster sum of squares. As such, the centroids are usually not actual data points, but rather the clusters' average, and will at best coincide with a data point by chance. The cluster labels are the result of random initialization and have no mathematical meaning. We therefore sort and label the clusters by average $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ between 100–400 m for the following analyses in Sect. 2.4. The variable k refers to the fixed and predefined number of clusters. The choice of k significantly affects the accuracy of the wind resource description, as well as the computational cost. The choice of k is informed by the elbow method, named after the characteristic line chart that resembles an arm, and the silhouette score. The “elbow” (the point of inflection of the curve) is a good indication that the underlying model fits well for the corresponding number of clusters. k can be chosen at a point where the inertia reduction becomes marginally small or decreases linearly (Pedregosa et al., 2011).

Absolute values of inertia (Fig. 5a) are not a normalized metric and scale with the size of the considered data set. On the other hand, the silhouette coefficients (Fig. 5b and d) are normalized between −1 (worst) and 1 (best). They indicate the membership of a data point to its cluster in comparison to other clusters, i.e., the proximity of each data point in one cluster to data points in neighboring clusters (Pedregosa et al., 2011). A negative value suggests that a data point could be assigned to the wrong cluster. The silhouette score, depicted by a dashed red vertical line, is the average of all silhouette coefficients for a fixed number of clusters k. For visualization purposes, k=10 clusters have been chosen. Each cluster is sorted and color coded according to average $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ between 100 and 400 m and color, the same as in Fig. 6. Performing these silhouette score analyses for multiple k results in the trend shown in Fig. 5c. A k value of 20 seems to be a good choice for the available data sets because the decrease in inertia for a higher number of clusters is only moderate, suggesting that the additional computational cost may not be worthwhile. Similarly, the silhouette score remains almost constant for higher numbers of clusters. Therefore, k=20 has been chosen for the analyses in Sect. 6.

Figure 7 Vertical onshore U categorized using the k-means clustering algorithm. The analyses in Sect. 6 employ k=20 clusters. Here k=10 is chosen for visualization purposes. The average profile (centroid) is shown in blue, and the profiles associated with this cluster are shown in gray. Clusters 1 to 10 (a–j) are sorted and labeled in ascending order of average centroid $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ between 100 and 400 m. The corresponding cluster frequency f for each cluster C is shown in Fig. 6. The optimization-toolbox-implemented U are highlighted with a red line.

2.4 Analyses of clustered profiles

For visualization purposes, the following subsections describe the wind conditions at both locations using only k=10 clusters. The analyses in the “Results and discussion” section (Sect. 6) use k=20 clusters.

Figure 6a and b show the average of the clustered wind speed profiles U, also referred to as centroids. Their colors correspond to the average $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ between heights of 100 and 400 m. All WRF-simulated U are depicted in gray. The cluster probabilities (Fig. 6c and d) are sorted by average centroid speed within the considered height range, represented by their colors and labels (C=1–10).

As expected, offshore U (Fig. 6b) at low altitudes are higher and wind shear is lower than onshore (Fig. 6a). In general, offshore centroids are more monotonic, as they do not exhibit a distinct U peak (i.e., LLJs), and achieve higher maximum wind speeds than onshore. The U profiles within each cluster cover a relatively small range, suggesting consistent clusters. Figure 7 (onshore) and Fig. 8 show the distribution of U within each of the clusters. At both locations, the first two clusters (Figs. 7a, b and 8a, b) exhibit very low wind shear with a low and almost constant U above 200 m. These low wind speed clusters amount to approximately 25 % onshore (Fig. 6c) and 20 % offshore (Fig. 6d), as can be seen in the corresponding cluster frequency of occurrence f. A standard logarithmic wind profile does not accurately describe such almost constant profiles which could lead to an overestimation of U at higher altitudes. AWESs must either be capable of functioning at low U or be able to safely land and take off autonomously. Onshore clusters 4 and 5 (Figs. 7d, e and 8d, e) seem to mostly consist of non-monotonic profiles, as these centroids show a distinct LLJ nose at about 200 and 300 m. The offshore centroids of clusters 7 and 8 (Figs. 7g, h and 8g, h) also show a slight wind shear inversion at higher altitudes.

Figure 8 Vertical offshore wind speed profiles categorized using the k-means clustering algorithm. The analyses in Sect. 6 employ k=20 clusters. Here k=10 is chosen for visualization purposes. The average profile (centroid) is shown in blue, and the profiles associated with this cluster are shown in gray. Clusters 1 to 10 (a–j) are sorted and labeled in ascending order of average centroid $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ between 100 and 400 m. The corresponding cluster frequency f for each cluster C is shown in Fig. 6. The optimization-toolbox-implemented U are highlighted with a red line.

Clusters C=1 (a) to C=10 (j) are sorted by the average centroid (blue line) wind speed between $\stackrel{\mathrm{‾}}{U}(z_{\mathrm{ref}}=\mathrm{100}$–400 m). The red lines indicate the profiles associated with the 5th, 50th, and 95th percentiles of $\stackrel{\mathrm{‾}}{U}\left(z_{\mathrm{ref}}\right)$ within each cluster. To reduce computational cost, only these profiles are later implemented into the AWEBox optimization toolbox. We selected these profiles because they are less likely to be irregular outliers than the extrema of the cluster while still representing the variation within their respective cluster. These profiles illustrate the variations within their respective cluster and are not average profiles like the cluster centroids or scaled or semi-empirical approximations such as the logarithmic wind profile. Evidently, the magnitude of the wind speed plays a dominant role in the clustering process. A clearer wind profile shape distinction could be achieved by normalizing the data before clustering it (Molina-García et al., 2019; Schelbergen et al., 2020).

Figure 9 Monthly frequency of k-means clustered onshore (top) and offshore (bottom) wind velocity profiles U for a representative k=10. All clusters are sorted and color coded according to their U(z_ref=100–400 m). The corresponding centroid associated with each cluster can be found in Fig. 6.

2.5 Analysis of clustered statistics

This subsection examines the relationship between the clusters and monthly, diurnal, and atmospheric stability. These analyses reveal patterns that give insight into the wind regime and the resulting changes in AWES performance. Subsequent sections examine wind data from k=20 clusters, while here only k=10 clusters are chosen for presentation purposes. Clusters are sorted in ascending order of average centroid wind speed $\stackrel{\mathrm{‾}}{U}(z_{\mathrm{ref}}=\mathrm{100}$–400 m) and color coded as in Figs. 5 and 6.

Both locations exhibit a clear annual pattern as shown in Fig. 9. High wind speeds are more common during winter, while low wind speeds are more prevalent in summer. This is likely due to the seasonal difference in surface heating and the resulting differences in atmospheric mixing. The two onshore and offshore clusters associated with the highest wind speed are almost exclusively present during November to February.

Offshore data indicate minimal diurnal variation as shown in Fig. 10, with only a slight increase in the frequency of lower-wind-speed clusters during the day. Onshore clusters, on the other hand, are more dependent on the diurnal cycle with a higher likelihood of low-speed clusters after sunrise. The frequency of onshore cluster 4, which includes a LLJ nose (Fig. 6), decreases to almost zero during the day and rises at night, supporting the notion that this cluster is linked to nocturnal LLJs.

Figure 10 Diurnal frequency of k-means clustered onshore (top) and offshore (bottom) wind velocity profiles U for a representative k=10. All clusters are sorted and color coded according to their U(z_ref=100–400 m). The corresponding centroid associated with each cluster can be found in Fig. 6.

Figure 11 Atmospheric stability (U: unstable; NU: nearly unstable; N: neutral; NS: nearly stable; S: stable; VS: very stable) distribution of k-means clustered onshore (top) and offshore (bottom) wind velocity profiles U for a representative k=10. The associated stability classes are based on Obukhov length (Table 2). All clusters are sorted and color coded according to their U(z_ref=100–400 m). The corresponding centroid associated with each cluster can be found in Fig. 6.

The wind velocity clusters (Fig. 11) show a correlation with atmospheric stability. Low wind speed clusters make up about 20 % to 30 % of the annual wind resource. These clusters exhibit Obukhov lengths close to zero (probably caused by very low friction velocity u_*) and are classified as “other” according to Floors et al. (2011) (Table 2). Unstable (U) and near-unstable (NU) conditions are associated slightly higher wind speeds. The highest wind speeds develop during neutral (N) and near-stable (NS) conditions. It should be acknowledged that strong winds driven by large pressure gradients can lead to neutral stratification. LLJ profiles associated with onshore cluster 4 are most likely to develop during stable (S) and very stable (VS) conditions.

In summary, k-means clustering can effectively group wind velocity profiles with similar characteristics up to high altitudes. These clusters are correlated with seasonal and diurnal changes, as well as atmospheric stability. The magnitude of the wind velocity profiles appears to have a greater impact on the resulting clusters than the shape of the profile. The algorithm is able to identify non-monotonic profile shapes that are less common, including profiles with LLJs. Normalizing the profiles before clustering can provide more information about the different shapes of the vertical profiles, but this was not done in this study.

3 Dynamic AWES model

This section introduces the dynamic AWES model used in the AWEBox trajectory optimization toolbox (De Schutter et al., 2020). Section 3.1 provides a summary of the system configuration. The aerodynamic model is presented in Sect. 3.2, while the aircraft mass model is introduced in Sect. 3.3.

3.1 Model configuration

The rigid-body model considers a 6-degrees-of-freedom (DOF) fixed-wing aircraft which is connected to the ground via a straight tether. The introduction of the tether reduces the DOFs to 5 (Terink et al., 2011). The model uses precomputed second-order polynomials to describe the aerodynamic coefficients (Sect. 3.2) which are controlled via aileron, elevator, and rudder deflection (Malz et al., 2019).

Figure 12 Ampyx Power AP2 reference kite aerodynamic lift c_L (a) and drag c_D,total coefficients (b) (Malz et al., 2019; Ampyx Power BV, 2020), including tether drag according to Eq. (2), for a projected wing surface area A=20 m² and tether diameter of d=7.8 mm. Tether length varies between 250 and 1000 m. Panel (c) shows the pitch moment coefficient c_m as a function of angle of attack α. The bottom figures display lift over drag (d), glide ratio over angle of attack (e), and $c_{\mathrm{R}}^{\mathrm{3}}/c_{\mathrm{D},\mathrm{total}}^{\mathrm{2}}$ over angle of attack (f).

The longitudinal dynamics of the tether is controlled via the tether jerk $\stackrel{\dots }{l}$ from which the tether acceleration $\ddot{l}$, reeling speed $\dot{l}$, and length l are determined. The tether is modeled as a single solid rod which supports neither compressive forces nor bending moments (De Schutter et al., 2019). The rod is divided into 10 segments. Tether drag is calculated individually for each segment, using the local apparent wind speed (Bronnenmeyer, 2018). The tether drag of every segment is equally distributed between the two endpoints. This leads to an underestimation of total tether drag at the kite (Leuthold et al., 2018). The ground station dynamics are not modeled explicitly but are implemented using a set of constraints. These constraints serve as an example of a system rather than representing a fully optimized design. A reel-in speed of ${\dot{l}}_{\mathrm{in}}=\mathrm{15}$ m s⁻¹ and reel-out speed of ${\dot{l}}_{\mathrm{out}}=\mathrm{10}$ m s⁻¹ are assumed to be realistic winch motor constraints based on information provided by a ground station manufacturer and literature review. This results in a reel-out to reel-in ratio of $\frac{\mathrm{2}}{\mathrm{3}}$. A maximum tether acceleration of $\ddot{l}=\mathrm{20}$ m s⁻² is imposed to comply with generator torque limits. The tether diameter is selected to be able to withstand 3 times (safety factor of SF = 3) the maximum tether tension of $F_{\mathrm{tether}}^{\max }=\mathrm{50}$ kN. This results in a rated average cycle power of about P_rated≈260–300 kW, according to AWEBox simulations.

3.2 Aerodynamic model

The presented model utilizes the Ampyx Power AP2 aerodynamic coefficients from De Schutter et al. (2020), Malz et al. (2019), and Ampyx Power BV (2020). The AP2 reference is scaled from a projected wing surface area of A_AP2=3 m² to A=20 m² to generate results for more realistic and probable devices, while the aspect ratio is kept constant at AR=10. The total combined drag coefficient of the aircraft and tether c_D,total,

$$\begin{array}{}\text{(2)}& c_{\mathrm{D},\mathrm{total}}=c_{\mathrm{D},\mathrm{kite}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}{\displaystyle \frac{ld}{A}}{c}_{\mathrm{D},\mathrm{tether}},\end{array}$$

depends on the diameter d and length l of the tether, as well as the projected surface area A and the aerodynamic drag coefficient c_D,kite of the kite. To illustrate the effect of a longer tether, we utilize a simple tether drag estimation for a cylindrical tether with constant diameter and an aerodynamic tether drag coefficient c_D,tether of 1.0. This value would be even higher for braided tethers. Assuming a constant and uniform wind speed, the line integral along the tether results in a total effective drag coefficient of $ld/\mathrm{4}A$ and accounts for the different reference areas for c_D,kite and c_D,tether. See Houska and Diehl (2007), Argatov and Silvennoinen (2013), and van der Vlugt et al. (2019) for details. Figure 12 visualized the tether drag impact on the aerodynamic coefficients for tether lengths up to l=1000 m.

Figure 13 Published actual (circle) and anticipated (square) aircraft mass scaling provided by Makani (red color scheme) (Echeverri et al., 2020) and Ampyx Power (blue color scheme) (Ampyx Power BV, 2020; Kruijff and Ruiterkamp, 2018). Diamond-shaped data (green color scheme) depict the mass of scaled-up versions of Ampyx Power prototypes used in research papers (Haas et al., 2019; Eijkelhof et al., 2020; van Hagen et al., 2023). For most data, mass scales with a scaling exponent between κ=2.2–2.6 (gray area). The chosen mass-scaling exponent of κ=2.4 is represented by a dashed line, and the investigated scaled AP2 design is highlighted by a red X.

The lift coefficient c_L (Fig. 12a) is approximated as a second-order polynomial function of angle of attack α to simulate stall effects. A polynomial description is necessary for the entire range of angle of attack, as the optimization algorithm requires a two-times differentiable function. For the sake of simplicity, a piecewise, continuous, and differentiable function has not been implemented. As a result, the implemented c_L (blue) slightly exceeds the linear (orange) lift coefficient $c_{\mathrm{L}}^{\mathrm{ref}}$ of the AP2 reference (Malz et al., 2019) between $-\mathrm{5}\le \mathit{\alpha }\le \mathrm{10}$^∘. The side slip angle β is included in the model, but variations in aerodynamic coefficients due to β are neglected. The pitch moment (Fig. 12c) is assumed to behave linearly. Changes in the drag coefficient (Fig. 12b) are approximated by a second-order polynomial. Tether drag is independent of α and therefore added to the zero-lift drag coefficient. The resultant aerodynamic force coefficient c_R is represented as

$$\begin{array}{}\text{(3)}& c_{\mathrm{R}}=\sqrt{c_{{\mathrm{L}}^{\mathrm{2}}}+c_{\mathrm{D},{\mathrm{total}}^{\mathrm{2}}}}.\end{array}$$

The drag polar in Fig. 12d depicts the relationship between the kite's lift coefficient c_L and total drag coefficient c_D,total for the tethered aircraft. The maximum values of the glide ratio $c_{\mathrm{L}}/c_{\mathrm{D},\mathrm{total}}$ (Fig. 12e) and the ratio $c_{\mathrm{R}}^{\mathrm{3}}/c_{\mathrm{D},\mathrm{total}}^{\mathrm{2}}$ (Fig. 12f), which is one of the main determining factors of AWES power (Loyd, 1980; Schmehl et al., 2013), decrease significantly with tether length and shift towards higher angles of attack. The impact of tether drag on the total drag coefficient is less significant for larger kites because its impact decreases with the size of the aircraft.

3.3 Aircraft mass model

The aircraft dynamics are described by a single rigid body of mass m_kite and moment of inertia J, subject to aerodynamic forces and moments. The inertial properties m_kite and J are determined by upscaling the AP2 reference kite from A_AP2=3 m² to A=20 m². The mass m_scaled and moment of inertia J_scaled of a fixed-wing aircraft scale as functions of the wing span b and aspect ratio AR, which is kept constant and whose impact on scaling is neglected here, with a mass-scaling exponent κ (Noth and Siegwart, 2006):

$$\begin{array}{}\text{(4)}& {\displaystyle }& {\displaystyle }{m}_{\mathrm{scaled}}=m_{\mathrm{ref}}{\left({\displaystyle \frac{b}{b_{\mathrm{ref}}}}\right)}^{\mathit{\kappa }},\text{(5)}& {\displaystyle }& {\displaystyle }{\mathbf{J}}_{\mathrm{scaled}}={\mathbf{J}}_{\mathrm{ref}}{\left({\displaystyle \frac{b}{b_{\mathrm{ref}}}}\right)}^{\mathit{\kappa }+\mathrm{2}}.\end{array}$$

Pure geometric scaling of solid bodies, in contrast to aircraft structures that use a lightweight structural frame, corresponds to Galileo's square–cube law with κ=3. In reality, as has been seen for the development of conventional WTs, design and material improvements occur over time. An appropriate mass-scaling factor was determined based on a review of the available literature. Figure 13 depicts actual (circle) and anticipated (square) aircraft mass scaling provided by Makani (red color scheme) and Ampyx Power (blue color scheme). The diamond-shaped data points (green color scheme) are scaled-up versions of Ampyx Power prototypes used in several research papers (Haas et al., 2019; Eijkelhof et al., 2020; van Hagen et al., 2023). The gray area encompasses most of the data points with κ=2.2–2.6. We chose κ=2.4 based on a curve fit of the available published sizing study data. This appears to be an ambitious goal for rigid kites but attainable for flexible ones. The mass of these hollow tensile structures filled with air mostly scales with the wing surface area, leading to significantly lower mass-scaling exponents and more beneficial mass scaling. Sommerfeld et al. (2022) examined the impact of various size, mass, and aerodynamic scaling factors on performance.

4 Optimal control model

AWESs need to dynamically adapt to changing wind conditions to maximize power generation and ensure safe operation. Section 4.1 introduces the dynamic trajectory optimization toolbox AWEBox (De Schutter et al., 2020). We describe the most important boundary conditions in Sect. 4.2 and initial conditions in Sect. 4.3. Section 4.4 explains the implementation of the previously described wind profiles (Sect. 2).

4.1 AWES model overview

Maximizing the average cycle power can be formulated as a trajectory optimization problem, which takes into account the interaction between the tether, kite, and ground station. This study analyzes the mechanical power produced by a single aircraft tethered with a straight line throughout one production cycle, including reeling in and out, while disregarding take-off and landing. Power production is intrinsically linked to the aircraft's flight dynamics, as the AWES never reaches a steady state over the course of a power cycle. Generating dynamically feasible and power-optimal flight trajectories is nontrivial, given the nonlinear and unstable system dynamics and the presence of various flight envelope constraints. Optimal control methods are a natural candidate to tackle such problems, given their inherent ability to deal with nonlinear, constrained multiple-input–multiple-output systems (De Schutter et al., 2019; Leuthold et al., 2018). This trajectory optimization is a highly nonlinear and non-convex problem which can have multiple local optima. The initial and final states of each trajectory must be equal to ensure periodic operation but are freely chosen by the optimizer. In periodic optimal control, an optimization problem is solved by computing periodic system states and control inputs that maximize a performance index (here average power output $\stackrel{\mathrm{‾}}{P}$) while satisfying the system's dynamic equations and constraints. We use this approach to generate a variety of realistic trajectories from WRF-simulated wind velocity profiles. Any wind data sets, such as wind atlas data, lidar, or met mast measurements, can be implemented into the optimization model via a two-times differentiable function, depending on the scope and purpose of the investigation.

4.2 Constraints

Several important constraints define the operational envelope. The most important constraints such as tether length, tether reeling speed, and tether force are summarized in Table 3. The following constraints define a representative and not design-optimized AWES. The power of a ground-generation AWES is limited by the tether force, which is defined by the tensile strength ${\mathit{\sigma }}_{\max }^{\mathrm{tether}}$, tether diameter d, and the tether reeling speed $\dot{l}$. The tether diameter is chosen such that the maximum tether tension is approximately $F_{\mathrm{tether}}^{\max }=\mathrm{50}$ kN with an additional safety factor of SF = 3. This produces a peak power of P_peak≈500 kW, with a maximum reel-out speed of $\dot{l}=\mathrm{10}$ m s⁻¹. This corresponds to a rated average cycle power of approximately P_rated≈260–300 kW. We assume a reel-out to reel-in tether reeling speed ratio of $\frac{\mathrm{2}}{\mathrm{3}}$ to be within winch design limitations. The tether length constraint is relatively lenient to allow the optimizer to investigate a wide range of possible operating heights z_oper. The flight envelope is constrained by limitations on the aircraft's acceleration and roll and pitch angles (to prevent collision with the tether), as well as the angle of attack α and side slip angle β. A minimal operating height of $z_{\min }=\mathrm{50}\mathrm{m}+\frac{b}{\mathrm{2}}$ is imposed for safety reasons.

Table 3 Selected AWES design parameters for the original AP2 reference system (Malz et al., 2019) and the scaled-up A=20 m² design analyzed in this study. Values in square brackets represent the upper and lower bounds, which are implemented as inequality constraints.

4.3 Initialization

The trajectory optimization process is highly nonlinear and non-convex, resulting in multiple local optima. These solutions depend on the chosen initial conditions. Some of the locally optimal solutions may be feasible and within the constraints but may have undesirable characteristics such as looping maneuvers during reel-in or excessively high operating altitudes. As a result, it is necessary to evaluate the quality of all solutions. To solve the complex optimization problem, initial guesses are generated using a homotopy technique similar to Gros et al. (2013). This technique initially fully relaxes the dynamic constraints using fictitious forces and moments to reduce model nonlinearity and coupling, which improves the convergence of Newton-type optimization techniques. The constraints are then gradually re-introduced until the relaxed problem matches the original problem. The trajectory optimization is initialized with a circular trajectory in downwind direction (positive x direction) with a fixed number of n_loop=5 loop maneuvers at a 30^∘ elevation angle, an initial tether length l_init=500 m, and an estimated aircraft speed of v_init=10 m s⁻¹ along the entire initial trajectory. This initialization is kept constant for all vertical wind velocity profiles. The number of loop maneuvers is not part of the objective function and remains unchanged during all optimization runs. Further investigation is needed to determine the impact of the number of loops. However, previous analyses have shown that the average cycle power estimated by AWEBox is relatively unaffected by the number of loops. It is likely beneficial to reduce the number of loops with wind speeds because the system can reel out faster at higher wind speeds and reach maximum tether length faster.

Figure 14 Representative WRF-simulated vertical onshore wind speed profiles U (a) and hodograph (top view) up to 1000 m (c). The highlighted sections indicate Lagrangian polynomial fit of the wind velocity at operating height. Panel (b) and panel (d) show the side and top view of the corresponding AWEBox-optimized trajectories. The reference wind speed in the legend is ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}=\stackrel{\mathrm{‾}}{U}(\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m})$. The results correspond to the time series shown in Fig. 16.

4.4 Wind profile implementation

To reduce the computational cost while maintaining an adequate representation, we only implement three wind velocity profiles from each cluster into the trajectory optimization toolbox. More profiles could be chosen for an in-depth analysis. The power for a total number of 60 wind profiles, three profiles for each of the k=20 clusters (Sect. 2.3), for each location is maximized. The three selected profiles correspond to the 5th, 50th, and 95th percentiles of average wind speed $\stackrel{\mathrm{‾}}{U}(z_{\mathrm{ref}}=\mathrm{100}$–400 m) within each cluster. We assume that these profiles represent the cluster's spectrum of wind conditions at operating height z_oper. The AWEBox includes a simplified atmospheric model based on international standard atmosphere to account for the variation in air density.

The vertical wind velocity profiles U^′ are rotated such that the main wind direction, which is defined as the average direction between 100 and 400 m, points in the positive x direction (Fig. 4). As a result the wind velocity components at every height consist of a main component u in the x direction and transverse component v in the y direction. The results can be seen in the hodographs of Figs. 14c and 15c.

The AWEBox toolbox uses the gradient-based MA57 solver (HSL, 2020) in IPOPT (Waechter and Laird, 2016) to solve the nonlinear control problem. Therefore, it is necessary to interpolate the vertical wind velocity profiles with a two-times continuously differentiable function. We chose to use Lagrangian polynomials (Abramowitz and Stegun, 1965) because the resulting polynomials pass through the input data points. To avoid over-fitting a limited number of data points are implemented. These data points are chosen based on the anticipated z_oper to best represent the wind conditions at relevant heights.

For comparison, logarithmic wind speed profiles,

$$\begin{array}{}\text{(6)}& U_{\log }=U_{\mathrm{10}}\left({\displaystyle \frac{{\log }_{\mathrm{10}}\left(z/z_{\mathrm{0}}\right)}{{\log }_{\mathrm{10}}\left(z_{\mathrm{10}}/z_{\mathrm{0}}\right)}}\right),\end{array}$$

with a roughness length of $z_{\mathrm{0}}^{\mathrm{onshore}}=\mathrm{0.1}$ and $z_{\mathrm{0}}^{\mathrm{offshore}}=\mathrm{0.001}$ are implemented into the trajectory optimization toolbox.

The reference wind speed U₁₀ at the reference height z₁₀=10 m varies from 3 to 20 m s⁻¹ with a step size of ΔU₁₀=1 m s⁻¹.

5 Reference models

This section introduces a simplified quasi-steady-state AWES reference model (QSM) (Sect. 5.1) and a steady-state WT model (Sect. 5.2). The so generated analytical solutions will serve as a reference for the optimization results.

5.1 AWES reference model

The QSM estimates the mechanical power of ground-generation AWES based on the assumption that the trajectory of the tethered aircraft can be approximated by a progression through steady equilibrium states where tether tension and total aerodynamic force are aligned. We simplify the QSM by approximating the reel-out and reel-in trajectory with a single state and neglecting the effects of gravity. The QSM, based on Argatov et al. (2009) and generalized by Schmehl et al. (2013) and van der Vlugt et al. (2019), approximates the aircraft as a point mass. Its position is described in terms of spherical coordinates, i.e., the radial distance from the ground station, the elevation angle ε, and azimuth angle ϕ relative to the direction of the mean wind velocity vector U. For lightweight soft-wing kites, this is a reasonably good approximation because the low mass of the kite leads to very short acceleration times. The model includes losses caused by the misalignment of the tether and the wind velocity vector. The same design parameters, system constraints, and wind conditions (Sect. 2) apply to the optimization model (Sect. 4.2), as well as the QSM reference model. We maximize the cycle average power P_QSM by varying l, $\dot{l}$ and z and assuming an optimal ratio $c_{\mathrm{R}}^{\mathrm{3}}/c_{\mathrm{D},\mathrm{total}}^{\mathrm{2}}$.

The average cycle power P_QSM,

$$\begin{array}{}\text{(7)}& P_{\mathrm{QSM}}={\displaystyle \frac{P_{\mathrm{out}}{t}_{\mathrm{out}}-P_{\mathrm{in}}{t}_{\mathrm{in}}}{t_{\mathrm{total}}}}=P_{\mathrm{out}}{\displaystyle \frac{{\dot{l}}_{\mathrm{in}}}{{\dot{l}}_{\mathrm{out}}+{\dot{l}}_{\mathrm{in}}}}-P_{\mathrm{in}}{\displaystyle \frac{{\dot{l}}_{\mathrm{out}}}{{\dot{l}}_{\mathrm{out}}+{\dot{l}}_{\mathrm{in}}}},\end{array}$$

can be estimated from the reel-out power P_out, the power losses during reel-in P_in, the reel-in time t_in, and reel-out time t_out. We assume reel-in power losses P_in to be zero because optimal reel-in tether tension is negligible. This reduces the average cycle power by up to 30 %, depending on wind speed. Due to the cyclic nature of the trajectory, we can determine the ratio of the reel-in time t_in and reel-out time t_out to the total cycle time t_total from the reel-in speed ${\dot{l}}_{\mathrm{in}}$ and reel-out speed ${\dot{l}}_{\mathrm{out}}$. ${\dot{l}}_{\mathrm{out}}$ depends on the wind speed, while ${\dot{l}}_{\mathrm{in}}=-\mathrm{15}$ m s⁻¹ is assumed to be the maximum reel-in speed.

During the reel-in and reel-out phases, we assume that the tether force F_tether and reeling speed remain constant. The time it takes to transition between these two phases is not taken into account. P_out is calculated from the product of tether reeling speed $\dot{l}$ and tether tension F_tether:

$$\begin{array}{}\text{(8)}& P_{\mathrm{out}}=F_{\mathrm{tether}}{\dot{l}}_{\mathrm{out}}={\displaystyle \frac{\mathit{\rho }}{\mathrm{2}}}A{U}_{\mathrm{app}}^{\mathrm{2}}{c}_{\mathrm{R}}{\left({\displaystyle \frac{c_{\mathrm{R}}}{c_{\mathrm{D},\mathrm{total}}}}\right)}^{\mathrm{2}}{\dot{l}}_{\mathrm{out}}.\end{array}$$

Tether tension is a function of the apparent wind speed U_app, air density ρ, and the resultant aerodynamic force coefficient c_R (Eq. 3). The apparent wind speed can be nondimensionalized by

$$\begin{array}{}\text{(9)}& {\displaystyle \frac{U_{\mathrm{app}}}{U\left(z\right)}}=(\cos \mathit{\epsilon }\cos \mathit{\varphi }-f)\sqrt{\mathrm{1}+{\left({\displaystyle \frac{L}{D}}\right)}^{\mathrm{2}}}.\end{array}$$

The tether reeling speed $\dot{l}$ is nondimensionalized by defining the reeling factor:

$$\begin{array}{}\text{(10)}& f={\displaystyle \frac{\dot{l}}{U\left(z\right)}}.\end{array}$$

The elevation ε and azimuth angle ϕ constrain f≤cos εcos ϕ because the magnitude of the apparent wind speed cannot be negative. Combining Eqs. (8) and (10) results in

$$\begin{array}{}\text{(11)}& P_{\mathrm{out}}={\displaystyle \frac{\mathit{\rho }}{\mathrm{2}}}AU\left(z{)}^{\mathrm{3}}{c}_{\mathrm{R}}{\left({\displaystyle \frac{c_{\mathrm{R}}}{c_{\mathrm{D},\mathrm{total}}}}\right)}^{\mathrm{2}}f\right(\cos \mathit{\epsilon }\cos \mathit{\varphi }-f{)}^{\mathrm{2}}.\end{array}$$

The optimal reeling factor $f_{\mathrm{opt}}=\frac{\mathrm{1}}{\mathrm{3}}\cos \mathit{\epsilon }\cos \mathit{\varphi }$ can be obtained from Eq. (11) through an extreme value analysis. We assume an average reel-out trajectory represented by a single crosswind state instead of tracking the actual trajectory. The trajectory center is aligned with the main wind direction (ϕ=0^∘). The elevation angle ε is determined using the tether length l and operating height z_oper. F_tether is constrained by the tether diameter d, the tensile strength ${\mathit{\sigma }}_{\max }^{\mathrm{tether}}$, and the safety factor SF:

$$\begin{array}{}\text{(12)}& F_{\mathrm{tether}}\le {\displaystyle \frac{d^{\mathrm{2}}}{\mathrm{4}}}\mathit{\pi }{\mathit{\sigma }}_{\max }^{\mathrm{tether}}.\end{array}$$

The power-harvesting factor ζ (Diehl, 2013) is an AWES performance metric:

$$\begin{array}{}\text{(13)}& \mathit{\zeta }={\displaystyle \frac{P}{P_{\mathrm{area}}}}={\displaystyle \frac{P}{\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }A{\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}^{\mathrm{3}}}}.\end{array}$$

The harvested power P is expressed relative to the kinetic wind energy flow rate P_area, P_area being a mathematical concept rather than a physical power flux, through an area equivalent to the wing surface area A to nondimensionalize the power. ζ can be derived from Eq. (8) by setting the elevation angle ε and the azimuth angle ϕ to zero.

5.2 WT reference model

This section introduces a simplified steady-state reference WT model that calculates power as

$$\begin{array}{}\text{(14)}& P_{\mathrm{WT}}={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{c}_{\mathrm{p}}^{\mathrm{WT}}{A}_{\mathrm{WT}}{U}^{\mathrm{3}}\left(z_{\mathrm{WT}}\right),\end{array}$$

with a hub height of z_WT=100 m for both onshore and offshore conditions. The rotor diameter D_WT≈26.9 m is sized such that an equivalent rated power of P_rated=260 kW is reached at a rated wind speed of $U_{\mathrm{rated}}(z_{\mathrm{WT}}=\mathrm{100}\mathrm{m})=\mathrm{12}$ m s⁻¹, assuming a constant power coefficient of $c_{\mathrm{p}}^{\mathrm{WT}}=\mathrm{0.45}$. The power is kept constant above the rated wind speed. The performance of the WT model, dynamic optimization toolbox, and QSM is estimated using the same sampled WRF-simulated wind conditions (Sect. 2).

6 Results and discussion

This section analyses the optimization results and compares them to the reference models. Section 6.1 investigates power-optimal trajectories and the time series of operational parameters. Section 6.2 examines operating height statistics, tether length, and elevation angle trends. Section 6.3 visualizes the impact of wind speeds at different reference heights on the power curve. We compare three different wind speeds: the wind speed at a reference height of 100 m ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}=U(z_{\mathrm{ref}}=\mathrm{100}\mathrm{m})$, the average wind speed at pattern trajectory height ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}=\stackrel{\mathrm{‾}}{U}(z_{\mathrm{ref}}=z_{\mathrm{PTH}})$, and the average wind speed at an a priori guess of pattern trajectory height ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}=\stackrel{\mathrm{‾}}{U}(\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m})$. The investigated power curves do not represent design-optimal performance. Section 6.4 examines the variation in average cycle power caused by realistic wind profiles simulated by WRF and compares them to reference power estimates based on logarithmic profiles. All results are subject to the same constraints and design parameters introduced in Sects. 3 and 4.

Figure 15 Representative WRF-simulated vertical offshore wind speed profiles U (a) and hodograph (top view) up to 1000 m (c). The highlighted sections indicate Lagrangian polynomial fit of the wind velocity at operating height. Panel (b) and panel (d) show the side and top view of the corresponding AWEBox-optimized trajectories. The reference wind speed in the legend is ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}=\stackrel{\mathrm{‾}}{U}(\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m})$. The results correspond to the time series shown in Fig. 17.

6.1 Flight trajectory and time series results

Figure 14 compares representative onshore and offshore power-optimal flight trajectories. These results have been chosen to visualize typical performance-optimized trajectories for realistic wind conditions determined with the AWEBox. The reference wind speed ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}$ in the legend is the average wind speed at the a priori guess of the pattern trajectory height $z_{\mathrm{ref}}(\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m})$.

Figure 14a shows the magnitude of the vertical wind velocity profile U. Figure 14c shows the corresponding top view of the wind velocity profile, rotated such that the main wind component (average wind direction between 100 and 400 m) u points in the positive x direction. The WRF-simulated wind profiles are shown in gray. The highlighted segments depict the Lagrangian polynomial fit (Abramowitz and Stegun, 1965) at operating heights. These polynomials that have been incorporated into the optimization toolbox provide a sufficient fit for the wind data. Figure 14b and d show a side (x–z plane) and top view (x–y plane) of the optimized trajectories. The optimization predicts an increase in tether length, operating height, and stroke length with wind speed. Figure 15 shows similar results for the offshore location. Figure 16 illustrates the corresponding temporal development of important operational parameters.

Figure 16 Time series of instantaneous tether tension (a), apparent wind speed (b), tether-reeling speed (c), angle of attack (d), power output (e), and tether length (f) over one pumping based on representative onshore WRF-simulated wind data. The results correspond to trajectories shown in Fig. 14.

Figure 17 Time series of instantaneous tether tension (a), apparent wind speed (b), tether-reeling speed (c), angle of attack (d), power output (e), and tether length (f) over one pumping based on representative offshore WRF-simulated wind data. The results correspond to trajectories shown in Fig. 15.

The optimizer maximizes tether tension by adjusting the reel-out speed and angle of attack (Fig. 16a) during reel-out even for lower wind speeds and adjusts the reel-out speed (Fig. 16c) to maximize average cycle power. This causes the reeling factor to exceed its optimal value of $f_{\mathrm{opt}}=\frac{\mathrm{1}}{\mathrm{3}}\cos \mathit{\epsilon }\cos \mathit{\varphi }$ at high wind speeds, resulting in an increase in power (Fig. 16e) even when the maximum tether force has been reached. The low-wind-speed example ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}=\mathrm{5.4}$ m s⁻¹ (blue) seems to be just above cut-in wind speeds. The tether reeling speed decreases to zero for a prolonged period during the reel-out phase in order to generate enough lift to keep the aircraft aloft. The production period remains almost constant (t≈60 s) for the moderate- and high-wind-speed trajectories (orange, green, and red), while the reel-in period increases with wind speed due to the increased reel-out length caused by a higher average reel-out speed. There are significant power losses during the transition between the production and retraction phases when the tether is reeled in, and the tension remains high because the aircraft is unable to de-power quickly enough. During the reel-in phase, the tether reeling speed reaches its limit and the tether tension decreases to zero. This happens as the aircraft reduces its angle of attack and lift, as shown in Fig. 16d. At higher wind speeds the optimizer increases the elevation angle and reduces angle of attack to stay within the constraints. This can result in odd or unexpected trajectories, even though these local minima are feasible solutions within the system constraints. Tether length (Fig. 16f) generally increases with wind speed as the system reels out faster, increases its elevation angle, and operates at higher altitude. Similar results for the offshore location are shown in Fig. 17.

6.2 Tether length, elevation angle, and operating altitude

Figure 18a illustrates the range of tether lengths l for each of the 60 onshore (blue) and offshore (orange) wind velocity profiles. The maxima and minima are highlighted by circles and plotted over the reference wind speed $\stackrel{\mathrm{‾}}{U}(z_{\mathrm{ref}}=\mathrm{100}$–400 m).

Figure 18 Tether length range (a) and frequency distribution of operating height z_oper (b) over reference wind speed U(z_ref=100–400 m) based on AWEBox trajectory optimizations of k=20 onshore (blue) and offshore (orange) clusters.

Figure 19 Onshore (a) and offshore (b) AWES power curve approximations with wind speeds at z_ref=100 m (blue), $\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m}$ (orange), and z_ref=z_PTH (green) reference heights. The dashed lines represent least-square spline interpolations that have been added to aid in visualization.

None of the optimizations reach the maximum tether length constraint of l_max=2000 m because a longer tether would not be advantageous due to the added drag and weight, which would decrease performance. Both locations show a trend towards longer tethers up to rated wind speed, where the reel-out speed and tension are almost constant and close to their respective constraint (Fig. 16). The maximum tether length remains almost constant above rated wind speed, while the minimum tether length increases slightly, reducing the total stroke length. The elevation angle (Fig. 18b) decreases as the tether length increases. The optimizer tries to keep the elevation angle low in order to reduce misalignment (cosine) losses between the tether and the wind velocity vector. The onshore elevation angle is slightly higher because of the increased wind shear, which makes higher operating altitudes more justifiable. This can also be seen in Fig. 18c which shows the frequency distribution of operating altitude. Of the optimal operating heights, 78.6 % onshore and 74.7 % offshore are below 400 m, confirming the findings in Sommerfeld et al. (2019a, b). Larger or multikite AWESs could benefit from higher operating altitudes due to their higher lift to tether drag ratio and lift to tether weight ratio, but more detailed analyses are required.

6.3 Impact of reference height on power curve

The power curve of wind energy converters quantifies the power that can be harvested at a given reference wind speed. For conventional WTs the wind speed at hub-height is commonly used as the reference wind speed. Whether this is appropriate for ever-growing towers and longer rotor blades is debatable (Van Sark et al., 2019; Wharton and Lundquist, 2012; European Wind Energy Association, 2012). Defining a reference wind speed for AWES is not trivial, as the operating height depends on the shape and magnitude of the vertical wind speed profile. The choice of reference wind speed impacts the power curve representation. The AWE glossary (Airborne Wind Europe, 2021) recommends using the wind speed at pattern trajectory height z_PTH, which is the expected or logged time-averaged height during the power production phase, as the reference wind speed. We estimate $\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m}$ as an a priori guess of the pattern trajectory height. Figure 19 compares onshore (Fig. 19a) and offshore (Fig. 19b) average cycle power over U(z_ref=100 m) (blue), $\stackrel{\mathrm{‾}}{U}(z_{\mathrm{ref}}=z_{\mathrm{PTH}})$ (green), and an a priori guess of the wind speed at pattern trajectory height $\stackrel{\mathrm{‾}}{U}(\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m})$ (orange).

Each data point corresponds to one of the sampled WRF-simulated wind velocity profiles U. The dashed lines, which are only added as visual aid, are a least-square spline interpolation of the approximately 60 data points. Based on these results, we can conclude that the selection of the reference height is more important for onshore conditions. The onshore wind conditions with their higher number of non-monotonic wind profiles and higher wind shear lead to larger deviations from the typical power curve shape described in Licitra et al. (2019), Airborne Wind Europe (2021), and others. The higher wind shear onshore leads to a shift towards lower wind speeds for a reference height of z_ref=100 m. The a priori pattern trajectory height guess of $\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m}$ is relatively close to the actual z_PTH, especially for lower wind speeds. At very high wind speeds above ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}\ge \mathrm{20}$ m s⁻¹ the z_PTH power shifts towards higher wind speeds indicating an increased operating height.

Figure 20 Average cycle power $\stackrel{\mathrm{‾}}{P}$ and power-harvesting factor ζ for the onshore (a, c) and offshore (b, d) locations as a function of average wind speed ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}$ between 100 and 400 m. The data points obtained from AWEBox (blue), QSM (orange), and the WT model (green) are based on WRF wind data and are compared to AWEBox data derived from standard logarithmic wind speed profiles (red).

The more homogeneous offshore wind conditions result in less power variation. The three different reference heights have almost no impact on the offshore power curve up to the rated wind speed. Above ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}\ge \mathrm{20}$ m s⁻¹ the power curves diverge and the average cycle power decreases. This seems to be a result of the AWEBox optimization and its initialization with a fixed number of loop maneuvers. As the wind speed and reel-out speed increase, the aircraft cannot complete all the loop maneuvers before reaching the maximum tether length and transitioning into reel-in. Therefore, one of the loop maneuvers is performed when already reeling in, leading to an increase in tether tension (Fig. A1a) and additional losses during the reel-in period (Fig. A1e). The corresponding trajectories are shown in Fig. A2.

6.4 Reference model power comparison

Figure 20 compares the variation in the power curve for a reference wind speed of ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}(\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m})$ based on sampled WRF-simulated wind data (blue) and power estimates based on standard logarithmic wind speed profiles (red). These results are verified against the QSM (Sect. 5.1, orange) and WT reference models (Sect. 5.2, green).

The QSM and WT reference model use the same sampled WRF-simulated wind data. No cut-out wind speed is defined. The cut-in wind speed of ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}\approx \mathrm{5}$ m s⁻¹ is the result of unconverged optimizations. The optimization algorithm was not able to find a feasible trajectory for these low wind speeds, indicating that the wind is insufficient to keep the AWES aloft and produce power. The QSM and WT model estimate power for these wind speeds. Rated power is achieved around U_rated≈12–15 m s⁻¹, depending on the wind speed profile shape. At these wind speeds the reel-out speed is almost constant, while a constant reel-out tension is already achieved at lower wind speeds (Fig. 16).

The logarithmic wind speed profiles (Eq. 6) use roughness lengths of $z_{\mathrm{0}}^{\mathrm{onshore}}=\mathrm{0.1}$ and $z_{\mathrm{0}}^{\mathrm{offshore}}=\mathrm{0.001}$. Onshore, the power predicted based on WRF wind data is often higher than the power predicted using logarithmic profiles (Fig. 20a). This is likely due to higher than predicted wind shear and the presence of LLJs that are not represented by logarithmic profiles. Offshore, the logarithmic and WRF data are in close agreement with the logarithmic results because most of the simulated wind profiles are more monotonic.

The ζ trends for both onshore (Fig. 20c) and offshore (Fig. 20d) conditions show a decrease with increasing wind speed and are consistent with the QSM. WT power fluctuates due to the choice of reference height. AWESs outperform WTs up to rated wind speed, particularly onshore where AWESs can take advantage of higher wind speeds aloft. Lower wind shear offshore reduces the need to operate at higher altitudes, reducing the benefit of AWESs. The QSM predicts the highest power, as anticipated, due to its simplified assumptions such as constant power during reel-out and reel-in and neglected mass.

7 Conclusions and outlook

This research compares the optimal performance of a single-aircraft, ground-generation AWES using mesoscale wind data simulated by the WRF model to the performance calculated using standard logarithmic wind profiles. It also describes trajectories, instantaneous performance, and trends in tether length and operating height. These analyses use 1 year of onshore wind data at Pritzwalk in northern Germany and 1 year of offshore wind data at the FINO3 research platform in the North Sea to drive the AWEBox optimization, which determines dynamically feasible, power-optimal trajectories. The annual wind data are categorized into k=20 clusters of vertical wind velocity profiles U using a k-means clustering algorithm. To decrease the computational expense, three profiles based on the 5th, 50th, and 95th percentiles of wind speed between ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}(\mathrm{100}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m})$ for each cluster are incorporated into the performance optimization model. The performance model uses a scaled Ampyx Power AP2 aircraft with a wing surface area of A=20 m² and is subject to realistic tether and operational constraints. Our investigation into the impact of wind speed at reference height found that the a priori guess of $\mathrm{100}\le z\le \mathrm{400}$ m is a good guess for the investigated AWES design and size. Optimal average cycle power is compared to a quasi-steady-state AWES model and a steady-state WT model.

The optimization model is able to determine power-optimal trajectories for complex, non-monotonic wind velocity profiles. The optimized results are only marginally lower than those obtained using the simplified QSM, which neglects the effects of gravity and only simulates a single reel-in and reel-out state instead of the entire trajectory. The predicted onshore AWES power exceeds the WT reference model. This is because AWES can adapt their operating altitude to benefit from higher wind shear or LLJs. Offshore wind velocity profiles are generally more monotonic and exhibit higher wind speeds with less turbulence and wind shear. As a result, offshore winds produce average power that is similar to their logarithmic approximation. Due to the initialization of the AWEBox with a fixed number of loop maneuvers, which is not a variable in the objective function, high-wind-speed trajectories show loops during the reel-in period which reduce the average cycle power. This can lead to a deterioration of the trajectory at high wind speeds, as the optimizer struggles to stay within the tether tension and tether reeling speed constraints.

An investigation of the time series data shows that the optimizer first maximizes tether tension by adjusting reel-out speed and angle of attack. As the wind speed increases, the tether reel-out speed approaches its maximum limit and becomes more constant, whereas at lower wind speeds, the reel-out speed varies more. Up to rated wind speed, when average tether tension and tether reeling speed are maximized, the optimizer increases the deployed tether length and reduces the elevation angle to operate at optimal height. At higher wind speeds, the elevation angle increases to de-power the system and stay within design constraints. As a result, approximately 75 % of the optimal onshore and offshore operating heights are below 400 m. The offshore power curve appears to be independent of the reference height due to the lower number of non-monotonic wind speed profiles and lower wind shear. In contrast, the choice of reference height is more important for the onshore power curve.

The mesoscale wind simulations, which include a year’s worth of wind data with a temporal resolution of 10 min, are analyzed and categorized for both onshore and offshore locations. The annual wind roses for heights of 100 and 500 m confirm the expected wind speed increase and clockwise rotation at both locations. Offshore shows a lower wind shear and veer than onshore. Annual wind speed statistics reveal that low wind speeds still occur at a fairly high probability up to 1000 m at both locations. The k-means clustering algorithm is able to categorize the wind regime and identify LLJs, as well as various non-logarithmic and non-monotonic wind profiles. The primary factor in assigning a profile to a cluster appears to be wind speed, while the shape of the profile seems to have a lesser impact. Individual clusters produce coherent groups of similar wind profiles whose probability correlates with seasonal, diurnal, and atmospheric stability variation. The k-means clustering method provides good insight into the wind regime, especially for higher altitudes where classification by Obukhov length is inadequate.

As a continuation of this study, the power curves and realistic wind conditions described here could be utilized to calculate AEP estimations. Further research is required into AWES power curves and their reference wind speed, which could be accomplished by deriving shape-specific power curves from normalized wind speed profiles or by considering the correlation between wind speeds at different reference heights. Future work should include a variable number of loop maneuvers as a variable in the optimization objective function. Using the same data and model, it is possible to investigate the annual and diurnal AWES power variation in comparison to WT performance. A parallel sizing study (Sommerfeld et al., 2022) using the same wind clustered wind data investigated the impact of mass and aerodynamic efficiency on AWES performance. Adding a design optimization to the AWEBox model could enable location-specific aircraft and tether investigation.

Appendix A: Figures

Figure A1 Time series of instantaneous tether tension (a), apparent wind speed (b), tether-reeling speed (c), angle of attack (d), power output (e), and tether length (f) over one pumping based on high-wind-speed WRF-simulated offshore wind data. The results correspond to the trajectories shown in Fig. A2.

Figure A2 High-wind-speed WRF-simulated vertical offshore wind speed profiles (a) and hodograph (top view) up to 1000 m (c). The highlighted sections indicate Lagrangian polynomial fit of the wind velocity at operating height. Panel (b) and panel (d) show the side and top view of the corresponding AWEBox-optimized trajectories. The reference wind speed in the legend is ${\stackrel{\mathrm{‾}}{U}}_{\mathrm{ref}}=\stackrel{\mathrm{‾}}{U}(\mathrm{100}\mathrm{m}\le z_{\mathrm{ref}}\le \mathrm{400}\mathrm{m})$. The results correspond to the time series shown in Fig. A1.