Study on the Wind and Wave Environmental Conditions of the Xisha Islands in the South China Sea

: Wind and waves are the main factors of environmental loading on ships and offshore structures. Thus, detailed understanding of wind and wave conditions can improve the design and maintenance of these structures. This paper developed a validated long-term wind and wave hindcast database covering the recent 32 years from 1988 to 2019. The spatial distribution of wind and wave characteristics for the whole Xisha Islands’ domain were analyzed. Frequency and directional distributions of wind speeds and significant wave heights were investigated at several locations around typical islands. Extreme value models were used to estimate the wind speed for 100-year return levels, whereas environmental contour approaches were utilized to establish the extreme sea-state parameters for 50- and 100-year return periods. It was found that the Weibull distribution was better fitted to the significant wave heights of the Xuande Atoll’s sites in the open sea, while the exponential Weibull distribution provided a better fit at the Yongle Atoll’s sites where waves are sheltered.


Introduction
Remote islands and reefs are critical for human society, allowing us to venture into the deep sea, deliver supplies from the sea, conduct maritime search and rescue, provide humanitarian relief, and enact disaster prevention and mitigation measures. In the South China Sea, the Xisha Islands serve as a major transportation hub, connecting the Chinese mainland with other countries. More and more offshore infrastructures have been planned on the Xisha Islands in order to fulfill various functional requirements. For both human activities and construction operations in this area, the environmental conditions are critical and serve as the most important guiding criteria [1].
Wind and wave data covering a sufficiently long time period can provide a statistical description of seasonal variations, as well as multi-year extreme conditions, for the appropriate design and operation of offshore structures. The long-term wind and wave characteristics have been widely investigated in the South China Sea, especially at some locations near the mainland [2][3][4][5][6][7]; however, due to a lack of observation data, the wind and wave environmental conditions of the Xisha Islands have rarely been investigated. Fortunately, numerical hindcast techniques for wind and ocean waves have been well-developed in recent years and validated in different engineering areas, including the Red Sea [8], the Baltic Sea [9][10][11], the Mediterranean Sea [12,13], the Adriatic Sea [14], the Black Sea [15], the south Atlantic Ocean [16], the Persian Gulf [17], and the Italian coast [18].
Although the Xisha Islands comprise the biggest archipelago in the Northern South China Sea, they only occupy a very small part of the map. Observation data and numerical hindcast data are rarely reported in this region [19,20]. The wave conditions at the Xisha Islands are mainly affected by swells from distant storms in the South China Sea and young waves generated by winds passing through the islands [21]. Similar conditions are encountered in the Madeira Archipelago [22], Hawaii [23], and Australian and Pacific Islands [24]. Considering the complex bathymetry surrounding the different islands and reefs, the spatial resolution used in previous studies of the South China Sea is too low to resolve the features around the Xisha Islands. As a result, standards, guidelines, and research reports on environmental conditions related to the development and utilization of the Xisha Islands are still severely lacking [25].
Utilizing a modern atmospheric model and a spectral wave model, the purpose of this work was to develop a long-term wind and wave data set for the Xisha Islands. Systematic methods were applied to analyze the marine engineering environmental conditions. Accordingly, the remainder of this paper is divided as follows: Section 2 introduces the data and methods used in this paper. Section 3 introduces the wind and wave climate analysis results. The mean and seasonal distributions of wind and wave data in the whole domain and specific locations are presented. In Section 4, an extreme value analysis of the wind speed and significant wave height for specific return periods is implemented using the annual maxima GEV method and the Peaks-Over-Threshold GPD method [26]. Extremes of wave environmental variables for engineering design are also determined using the full sea-state environmental contour method [27]. Finally, a summary and our conclusions are given in Section 5.

Wind Data
For driving ocean wave models, it is critical to have reliable and precise high-resolution surface wind data both temporally and geographically. The initial wind data source for this study comes from the Global Data Assimilation System. The typhoon initial field was formed using an upgraded National Center for Atmospheric Research-Air Force Weather Agency typhoon bogus technique paired with Weather Research Forecast Variational Data Assimilation System. The forcing field of wind was defined covering the period from 1988 to 2019 (wind speed of 10 m above sea level), with a spatial resolution of 0.1° × 0.1° and a temporal resolution of 3 h. The coverage space ranged from 105° E to 125° E and 9° N to 24° N. Validation of the wind field data can be found in [28].

Unstructured Grid for Wave Modelling in the South China Sea
As shown in Figure 1, the Xisha Islands form the biggest archipelago in the Northern South China Sea. Coral reefs and atolls are widely distributed in the Xisha Islands; in particular, the Xuande Atoll and the Yongle Atoll are two main atolls where human activities are frequent. The Xuande Atoll lies in the northeast of the Xisha Islands, at 16°53′ N 112°17′ E. The group consists of low, narrow islands with sand cays and enclosed shallow lagoons connected by reefs. The Yongle Atoll lies about 70 km southwest of the Xuande Atoll, at 16.5° N 111.7° E. It consists of islands and reefs that form a crescent-like structure from west to east, enclosing a deep central lagoon.
Similar to the issues encountered in the Hawaiian Islands [29], wave transition modelling from the South China Sea to the reefs and islands requires high geometric resolution in this area. As the Xisha Islands only occupy a very small area of the South China Sea, grid modelling becomes a challenge, due to the large span of spatial bathymetry. To tackle this issue, a triangular unstructured grid model was constructed in the third-generation spectral wave model WAVEWATCH-III (WW3). As shown in Figure 2, the resolution of the unstructured numerical grid varied from 3 km to 300 m, in order to distinguish the effects of the reefs and islands. The Xisha Islands are distributed over a maritime area of around 15,000 square kilometers, with a land area of approximately 7.75 square kilometers. A gird resolution of 300 m was enough to resolve the islands and reefs. The boundary condition of the WW3 model was set as an open boundary at the right side and bottom side of the domain. Bulk wind and wave parameters were outputted and archived for all the grid points spanning the Xisha Islands, comprising the data source used in this study.

Validation of Hindcast Wave Characteristics
Since 2014, an observation system has been installed in a typical Crescent Group reef off the coast of the Xisha Islands. At the three locations marked with black stars in Figure  3, instrumented wave observations were carried out to capture wind-generated wave states. The measurement devices consist of two Datawell directional waverider devices (Buoy1 and Buoy2), as well as an ADCP (acoustic Doppler current profiler; Buoy3). Detailed information of our observation system can be found in [21].
Buoy recordings were used to ensure that our hindcast data were accurate. Figure 4 shows a comparison of the time-series Hs data at three Buoy locations in different measuring periods, from which it can be seen that a high level of consistency was reached. Table 1 provides the statistical comparisons between wave WW3 hindcast data and buoy data. In general, Pearson's correlation coefficients indicated that the model and buoy Hs were in good agreement. The Pearson's correlation coefficients were 0.881, 0.904, 0.922. The Bias were 0.021 m, 0.076 m, and 0.025 m of Buoy1, Buoy2, and Buoy3, respectively.
These results indicate that the unstructured wave forecast model can be used to accurately simulate wave propagation near islands and reefs. The computational efficiency of the unstructured numerical grid was also tested in [21]. Figure 4 shows an overall good agreement between the simulation and observation. It is also important to know the agreement in severe conditions. In Figure 5, a detailed comparison of Hs during typhoon Kalmeage waves is shown. The simulated peak Hs matches well with the observed peak Hs during the typhoon passage at both Buoy1 and Buoy2. In shallow water, the simulated peak Hs is almost the same at Buoy1. In deep water, the simulated peak Hs is 0.18 m larger than the observed peak Hs at Buoy2.

Extreme Value Analysis Methods
Extreme value theory has been widely used to estimate the regression values of meteorological and ocean data, such as wind speed, significant wave height, tide level, precipitation, and so on [26]. The extreme value analysis of time-series data can be divided into two main methods: the first involves pre-processing the data by taking the block maximum value for a sufficient period of time, such as the annual maxima value (AM). The block maximum is theoretically calculated with the Generalized Extreme Value distribution function (GEV) [30]: where μ, σ, and ξ are location, scale, and shape parameters, respectively. The second method is to perform Peak Over Threshold processing on time-series data, using the Generalized Pareto distribution for probability-distribution fitting. The Peaks-Over-Threshold method is a sampling technique that uses exceedances above a given threshold to create a sample. The Generalized Pareto Distribution (GPD) can be used to approximate a sample of exceedances beyond a high threshold: where σ and ξ are scale and shape parameters, respectively, and u is the threshold. When designing marine structures, such as wave energy devices and offshore wind turbines, it is necessary to consider their survivability in extreme sea conditions. The main purpose of extreme value analysis is to estimate the probability distribution of extreme value variables based on empirical sample data, in order to calculate the extreme sea state that represents a small probability event once in n years. The environmental contour approach has been widely used in the process of designing offshore structures, as well as in standards and guidelines such as the DNVGL-RP-C205 [31]. The idea is to define contours in the sea-state space (usually Hs, Tz), along which extreme responses with a given return period should occur [32]. This method finds the most likely design point of structural failure on the environmental contour line of a specific return period, which is input as a design parameter to obtain the same return period structural response, allowing us to avoid long-term response forecasting [27]. The process of developing an environmental contour typically involves statistical modelling of the joint probability density function of the variables of interest, then constructing a contour based on that joint distribution [33,34]. In the first step, the marginal distribution of Hs is typically modeled with the two-parameter Weibull distribution: The exponentiated Weibull distribution is a novel distribution choice for sea-state models. It extends the two-parameter Weibull distribution with a second shape parameter, δ , used as an exponent of the cumulative distribution function (CDF): for Hs > 0; here, β > 0 is the first shape parameter, δ > 0 is the second shape parameter, and α > 0 is the scale parameter of the distribution. In the case δ = 1, the exponentiated Weibull distribution becomes the two-parameter Weibull distribution. The zero up-crossing period follows a conditional log-normal distribution. For the second step, the inverse first-order reliability method (IFORM) [32] was used for the construction of the contour.

Spatio-Temporal Variations of Wind Characteristics
The wind speed and wind direction change with time and height above sea level. The standard reference height is 10 m above the sea surface. The common practice in ocean engineering is to give the average duration of wind speed and reference elevation, combining them with local wind-related parameters such as the wind-profile function, gust coefficient, and autocorrelation wind spectrum, in order to describe wind conditions at different temporal and spatial scales reasonably [29]. Figure 6 shows the spatial distribution characteristics of the hourly averaged wind speed (WSD) covering the Xisha Islands. Note that, direction in this paper is defined as those directions where the wind and waves come from. In general, the wind field did not change significantly in such a small spatial scale over the Xisha Islands.  Figure 7 shows the spatial distribution of the seasonal mean WSD at the Xisha Islands. The mean WSD in winter was the highest, exceeding 8.5 m/s in most sea areas of the Xisha Islands; meanwhile, the mean WSD in spring was the lowest, ranging from 5-6.5 m/s. In summer, the mean WSD was higher in the southeast of the Xisha Islands. In autumn and winter, the mean WSD was higher in the northeast of the Xisha Islands.

Spatio-Temporal Variations of Wave Characteristics
Wave conditions are critical for the design and deployment of marine structures such as wave energy converters and floating wind turbines. Through the bulk wave parameters output by WW3, the mean values of the significant wave height (Hs) at each grid point were calculated. Figure 8 shows the distribution of the mean Hs over the Xisha Islands. Due to the variations in reef topography, the mean Hs in the Xisha Islands was unevenly distributed. In the deep-water area of the Xisha Islands, the mean Hs ranged between 1.5-1.8 m, while it decreased rapidly when approaching the islands and reefs: the mean Hs inside the atoll was only 0.2-0.4 m. The overall wave conditions of the Xuande Atoll were worse than those of the Yongle Atoll.  From the components of the wave frequency spectrum, waves can be divided into wind waves, swell waves, and mixed waves. In the calculation process of WAVEWATCH-III, the wind wave and swell parts can be separated by their dynamic frequency, based on the wave age. This numerical method of spectral separation further divides the wave spectrum into wind wave components and swell components. Then, the wind wave ratio coefficient, representing the contribution of the wind wave component to the total wave energy, can be derived. Figure 10 shows the distribution of the average annual wind wave ratio coefficient in the Xisha sea area. In the open sea, wind wave components were dominant, with most of the sea area ratio coefficients being above 85%; in the sheltered sea areas affected by islands and reefs, the wind wave components were reduced, and the impact of swells on these areas was relatively large. An exemplary area is the southwest side of the Xuande Atoll, where the wind wave component in this area only accounts for 50-60%.

Long-Term Wind Characteristics
The 10-m wind speed and direction data at a random location were extracted from the hindcast wind field data to represent the wind characteristics in the Xisha Islands. The observation period was from 1988 to 2019, updated once every three hours. Table 2 shows the wind speed classification results according to the Beaufort wind force scale. From the table, it can be seen that the main wind speed frequency was concentrated from Class 3 to Class 6, accounting for 64.11%. The frequency of strong winds above Class 6 was 4.08%, mainly in November and December. The maximum wind speed at the site was 40.57 m/s, which occurred in September. The highest mean wind speed was in December, reaching 10.49 m/s. The lowest wind speed standard deviation was observed in March, while the largest standard deviation was in November.  Figure 11 shows a rose plot of the wind speed and direction at the Xisha Islands throughout the entire hindcast period. The dominant wind direction in the rose plot appeared in the N-E and S directions. The wind speed distribution in the N-E direction was wide, while the wind speed distribution in the S direction was relatively concentrated and small.  Figure 12 shows the rose plots of wind speed and wind direction divided by month. From March to April, the wind direction was relatively scattered and distributed in N-E to S directions; meanwhile, from May to August, the wind direction was concentrated in the S direction. September was a transitional month. From October to February of the following year, the wind direction was stable and concentrated in the N-E direction.

Long-Term Wave Characteristics
As the Xisha Islands consist of a large number of islands and reefs, the wave conditions may significantly differ from island to island. The high-resolution wave data calculated by the unstructured WW3 model allowed us to calculate long-term wave statistics and conduct analysis on specific islands and reefs. As shown in Figure 13, six sites were selected near important islands with busy human activities. These sites include Island A, Island B, and Island C in the Xuande Atoll, as well as Island D, Island E, and Island F in the Yongle Atoll. Among these characteristic sites, Island E is located inside the lagoon, with a water depth of 30 m. The other sites are distributed in the open sea, near islands and reefs. The water depth at these sites ranges between 40 m and 60 m. The summary information of these sites is provided in Table 3.  The long-term distribution of wave conditions can be represented by a scatter diagram of sea-state parameters. The scatter diagram reflects the historical occurrence frequency of a given sea state. In this paper, Hs and Tz (zero up-crossing period) are used for the description of the sea state. Figures 14-19 show the respective wave scatter diagrams at the reference sites. The wave conditions at Islands A, B, and C of the Xuande Atoll were similar, while the wave conditions at Islands D, E, and F of the Yongle Atoll differed. These results can be used as a reference condition for the fatigue analysis of offshore engineering structures.    .2 1,223.5 1,522.1 1,600.6 1,511.7 1,164 Figure 20 shows rose plots of the Hs and mean wave direction for the characteristic sites. It can be seen that the main wave directions at the three characteristic sites of the Xuande Atoll were all in the N-E direction, while Island A and Island B were also under the influence of waves from the SSE direction. The directional distributions of waves near the different sites of the Yongle Atoll were quite different. The main wave direction at the Island D site was the N-E direction. Meanwhile, as Island E is inside the Yongle Atoll lagoon, there was no obvious main wave direction. The wave directions with a higher frequency at Island E were NNW and S. Island F was dominated by the waves from the E direction.

Extreme Value Analysis on Wind Speed
The AM-GEV method was used for the wind speed extreme-value analysis. The first step was to select the maximum wind speeds appearing in the 32-year wind speed time history data as a sample, and then perform empirical function fitting. Figure 21 shows the relationship between the calculated wind speed regression period and the regression value. The 100-year wind speed calculated by the AM-GEV method was 50.22 m/s. A 95% confidence interval for the wind speed was obtained as [42.94 m/s, 55.68 m/s]. To perform the peak over threshold method, the storm process affecting the Xisha Islands was first determined. The storm period was extracted from the best track data of the Joint Typhoon Warning Center. The total number of storms was found to be 220. Figure 22 shows the paths of all tropical storms affecting the Xisha Islands from 1988 to 2019, classified by months, while Figure 23 shows the monthly distribution of tropical storm numbers. The number of tropical storms had obvious seasonal characteristics. The typhoon season started from April, and the number of storms gradually increased. The frequency and intensity of typhoons reached a peak value in September. From December to March of the next year, almost no tropical storms occurred. From the general perspective of the storm path, most tropical storms entered the South China Sea from the Northwest Pacific Ocean. Summer tropical storms tend to land on the coast of China, while winter storms tend to land on the coast of Vietnam. Taking the position of the Xisha Islands as a reference position, most of the typhoon tracks were spread from the northwest to the southeast of the Xisha Islands.  Once the start and end times of the storms were determined, the maximum wind speed during the 220 storms could be extracted from the time history data, which can be used as an independent sample of storm extreme values. Of the 220 storm maxima samples, some of them had a limited impact on the Xisha Islands. Therefore, a threshold was used to remove these samples. According to the recommendation in [27], we used 1.5 times the average wind speed as the threshold. Figure 24 shows the regression curve of the GPD using the Peaks-Over-Storm threshold. The calculated wind speed of the 100year return value was 45.19 m/s. A 95% confidence interval for the wind speed was obtained as [39.77 m/s, 49.36 m/s].

Extreme Value Analysis on Significant Wave Height
Considering that the length of our simulation wave data is 32 years long, the block period is chosen as a year. The estimation of the GEV and GPD model parameters are based on the maximum likelihood estimator. The return level plots using the AM-GEV method are shown in Figure 25. Of all the sites, Island A had the largest storm Hs larger than 10 m during the last five years. The GEV distribution was reasonably well-fitted except for Site F. To solve a critical problem in the POT-GPD method, an appropriate threshold needs to be selected. Different thresholds will lead to a large variety of return level estimates. Following [21], threshold selection was done on the basis of the evaluation of the mean residual life (MRL) trend in the function of a defined threshold range. More precisely, the threshold should be chosen as the point in whose proximity the trend changes from approximately linear to unstable. The mean residual life plots at six sites are shown in Figure  26, with thresholds at 4.3 m, 4.5 m, 4.4 m, 2.6 m, 1.0 m, and 2.1 m, respectively. The return level estimation results using the POT-GPD method are shown in Figure 27. After proper threshold selection, the POT data can be well fitted by the GPD model.
The summary and comparison of the 100-year return levels estimated by the GEV approach and POT-GPD approach are summarized and listed in Table 4. The confidence intervals legitimize the differences among the different approaches. The best fitted results of the two approaches were basically similar at Site D and Site E. The POT-GPD method tended to underestimate the return levels at Site A, Site B, and Site C, while the GEV method tended to underestimate the return levels at Site F.    Figures 28 and 29 show the fitting results of significant wave heights using the twoparameter Weibull distribution and the exponentiated Weibull distribution, respectively, at different sites. In Figure 28, the QQ-plots show that the Weibull distributions provide good model fit at Islands A, B, and C. The fitting parameters for these sites were similar. However, the Weibull distributions predicted too low wave heights at high quantiles for Islands D, E, and F. In Figure 29, the results are the complete opposite. The QQ-plots show that the exponentiated Weibull distributions provided a fairly good fit at Islands D, E, and F, while predicting too low wave heights at high quantiles for Island A, B, and C. The fitting parameters are listed in Table 5. It should be noted that the sea states at Site E and Site F are much calmer than the other sites at normal times due to the sheltering effect of the surrounding reefs. The significant wave heights at Site E and Site F are relatively smaller and compact. As a result, the value of the scale parameter α estimated using the exponentiated Weibull distribution is relatively small compared with the other sites.

Extreme Sea-State Estimation at Selected Sites
Generally, it can be concluded that the Hs data of the three sites in the open sea at the Xuande Atoll fit well with the two-parameter Weibull distribution, while the exponentiated Weibull distribution fitting results were better for the sites at the Yongle Atoll, which is sheltered by reefs and islands. After determining the well-fitted distributions, the corresponding contour plots of the sea-state parameters were calculated and are shown in Figure 30. After the environmental contour method was estimated, the maximum Hs on the contour line was selected. The single-value extreme value analysis results were compared with the environmental contour results in Table 6. Generally, the 100-year return level Hs estimated using the Environmental Contour method was close to the GEV results. An advantage of the Environmental Contour method was that it could avoid the unrealistically large values that could appear in the extreme-value analysis results.

Conclusions
The wind and wave environmental conditions in the Xisha Islands of the South China Sea were described qualitatively and quantitatively in this research. First and foremost, the numerically modeled long-term wind and wave data employed in this study were thoroughly explained and validated. The wind data had a spatial resolution of 0.1° and covered the whole South China Sea. The resolution of unstructured wave hindcast data varied from 3 km to 300 m, enabling the examination of near-island sites in the South China Sea at the Xisha Islands.
First, the spatial wind and the wave environment were analyzed. In general, the mean wind speed ranged between 7-8 m/s, with little spatial variation. The spatio-temporal variations in the wave characteristics for the Xisha Islands were determined in terms of annual mean Hs, seasonal mean Hs, and annual mean wind wave ratio coefficient. The impact of swell components was more evident in the sheltered areas, affected by islands and reefs.
Then, the frequency and directional distribution of wind and wave characteristics at various sites were analyzed. The main wind speed frequency throughout the year was concentrated from Class 3 to Class 6, accounting for 64.11%. The dominant wind directions in the rose plot appeared in the N-E and S directions. Six sites near the main islands were analyzed in order to provide detailed and structured insights into the wave conditions of the Xisha Islands. The scatter diagrams and rose plots at these sites facilitated the quantitative description of the sea-state conditions in normal times.
In the extreme wind and wave analysis section, the extreme wind speed was first evaluated, with respect to a 100-year return period, using the AM-GEV method and the POT-GPD method. The 100-year wind speed calculated using the AM-GEV method was 50.22 m/s, while the wind speed calculated with a 100-year return value using the POT-GPD method was 45.19 m/s. The best track data from the Joint Typhoon Warning Center were also used to study tropical storm events. There were 220 tropical storms over the last 32 years in the study area, where most tropical storms made their way into the South China Sea from the Pacific Northwest. Summer tropical storms tend to hit the Chinese coast, whereas winter storms are more likely to hit the Vietnamese coast.
Finally, the extreme wave environment was analyzed. The extreme Hs was evaluated at six sites using the traditional AM-GEV method and the POT-GPD method. It was found that these two approaches have different performances at different sites. Despite the similar results at Site D and Site E, the POT-GPD method tended to underestimate the return levels at Site A, Site B and Site C, while the GEV method tended to underestimate the return levels at Site F.
Extreme wave environmental contours were then established for multi-variable extreme analysis at the six considered sites. Joint probability distributions for significant wave heights and mean wave periods were developed. It was found that the Hs data of the three sites in the open sea fit well with the two-parameter Weibull distribution, while the exponentiated Weibull distribution fitting results were better for sites sheltered by reefs and islands. Accordingly, generally, the 100-year return level Hs estimated using the Environmental Contour method was close to the GEV results. An advantage of the environmental contour method was that it could avoid the unrealistically large values that could appear in extreme value analysis results.
This research can aid in gaining a better understanding of the wind and wave climate in the Xisha Islands, which can then be used for the design, operation planning, and maintenance of marine-related engineering objects.