1 Introduction

The Arctic is believed to be the region with the most abundant undiscovered hydrocarbon deposit in the world. Even though the Arctic encompasses only 6% of the Earth's surface, 30% of unexplored gas and 13% of unexplored oil can be found there, particularly in the Russian Arctic (Gautier et al., 2009). However, this cold area has extreme oceanic environmental conditions, such as low temperature, long polar nights, polar lows, dense fog, and extensive ice and snow cover. In particular, the formation of sea ice can not only influence the marine environment, atmospheric circulation, and climate change (e.g., Budikova, 2009; Mori et al., 2014; Zhao et al., 2019), but also threaten shipping navigation, drilling platform operation, and other marine activities (e.g., Timco and Weeks, 2010; Marchenko, 2014; Gudmestad, 2018).

The southern Kara Sea (SKS) is a semi-closed shelf sea and contains 39% potential oil-gas resources, ranking first in the Russian Arctic. As shown in Figs.1 and 2, the SKS is surrounded by different lands and seas, with the west separated from the Barents Sea by Novaya Zemlya, the south and east bounded by the Eurasian continent, and the north connected to the northern Kara Sea by the boundary (denoted by the dotted line at the center) between Cape Zhelaniya and Dikson Island. The bathymetry is also favorable to oil-gas exploitation in terms of the average depth of 110 m, with 40% less than 50 m. In a Russian plan, 60% of oil-gas production in 2035 is expected from the SKS fields (Bekker et al., 2015). In general, sea ice is the predominant natural threat to oil-gas exploitation and marine structures in the Arctic. Thus, spatiotemporal analysis and intensity assessment of sea ice conditions in the SKS have become important research topics.

In general, sea ice is an inevitable challenge to oil-gas exploitation in the SKS (Efimov et al., 2020). Different from the western Barents Sea, which is toward an ice-free zone even in the winter months (Onarheim and Årthun, 2017; Duan et al., 2018), the fixed offshore platform can operate all year round without ice threats. Subjected to the blocking effect of Novaya Zemlya, large amounts of North Atlantic warm waters cannot freely enter the SKS (Duan et al., 2019). Thus, the SKS maintains longer cold and frozen conditions. Ice-free periods are limited to 2–3 months on average and have prominent interannual variability (Duan et al., 2019; Efimov et al., 2020). In this seasonally ice-covered area, the mobile floating platform is a feasible choice (see Fig.2). In open waters, the mobile floating platform can operate safely; in light ice conditions, it can still work with the protection of ice-resistant structures, such as shockproof walls; and in the case of hazardous ice loads, it can be promptly disconnected and evacuated from the drilling site (Kutvitskaya and Ryazanov, 2013).

At present, the regulations for polar sea ice engineering in China are still incomplete. Therefore, the determination of the intensity of sea ice for the design and operation of ice-resistant floating platforms in the SKS requires comprehensive knowledge. Level ice thickness needs to be determined and carefully considered. The mechanical properties of ice and their associated theoretical formulas are closely dependent on their physical properties, such as compressive and flexural strength (Timco and Weeks, 2010; Kovalev et al., 2019; Chai et al., 2021). However, most relevant investigations concentrated primarily on the spatiotemporal sea ice extent trends and their associated relationship with metocean factors in the SKS (e.g., Belchansky et al., 1995; Divine et al., 2004; Rodrigues, 2008; Zubakin et al., 2008; Cavalieri and Parkinson, 2012; Ahn et al., 2014; Matishov et al., 2014; Bushuk and Giannakis, 2017; Duan et al., 2019, 2020; Efimov et al., 2020). Few studies have emphasized the sea ice thickness distributions across the SKS because of the scarcity of long-term successive level ice thickness observations.

Generally, several methods can be possibly used to determine ice thickness in the Arctic. 1) Drilling holes is the most accurate approach, as demonstrated in the fieldwork conducted by Romanov (2004) and Kovalev et al. (2019). 2) Airborne electromagnetic observation is an efficient, fast, and high-precision method that has been successfully used to measure the pack sea ice areas in March and April (Haas et al., 2010). 3) Subsea upward-looking sonar can obtain measurements over relatively long distances (Birch et al., 2000). However, these three methods require considerable labor and material costs and are often limited to strict geopolitics. Furthermore, the observations usually vary by region at different times of the year, making it impossible to convert the scattered data into an area-averaged spatial ice thickness picture. 4) Satellite monitoring: Since 2002, with the launch of ENVISAT, CryoSat-2, and SMOS satellites, large grid-scale winter ice thickness can be determined using remote sensing technology. However, the retrievable data show possible uncertainties and are easily contaminated by the poor atmospheric environment (Huntemann et al., 2014; Ricker, 2014). 5) Stefan's empirical formula: In this hypothesis, the level ice thickness is assumed to be proportional to the square root of the sum of freezing degree days (Ashton, 1986). This equation only considers the air temperature over the sea surface, not the sea ice dynamics.

To overcome the aforementioned defects, numerical simulation is another effective and powerful technology used to reproduce the historical sea ice freezing–melting processes and fill the gap of long-term successive level ice thickness in the SKS. A short-term attempt during the winter of 2008–2009 was conducted in several coastal sites using the thermodynamic sea ice model HIGHTSI (Cheng et al., 2013; Similä et al., 2013), but it did not examine the ice thickness spatial pattern across the SKS. Some reanalysis data, such as TOPAZ4, PIOMAS, ORAP5, ORAS5, CFSR, and CFSv2, have been available. The TOPAZ4 product was well assessed (Xie et al., 2017) and even recommended for the medium-range predictability of early summer sea ice thickness distribution in the East Siberian Sea (Nakanowatari et al., 2018). However, the TOPAZ4 product did not assimilate sea ice thickness from the CryoSat-2 and SMOS satellite until 2014. No reliable conclusions on sea ice thickness in the SKS were obtained because no direct comparison with in situ observations, particularly in periods without assimilation, was conducted. Moreover, the Los Alamos Sea Ice Model (CICE) (Hunke et al., 2015), which integrates the complete sea ice physical theories and parameterizations, has been widely applied to conduct scientific studies of Arctic sea ice (Wu et al., 2015). The CICE provides three types of thermodynamics, i.e., zero-layer module (Semtner, 1976), Bitz-Lipscomb module (Maykut and Untersteiner, 1971; Lipscomb, 1998; Bitz and Lipscomb, 1999), and mushy module (Feltham et al., 2006), and two types of dynamics, i.e., elastic-viscous-plastic (EVP) module (Hunke, 2001; Hunke and Dukowicz, 1997, 2002, 2003; Bouillon et al., 2013) and elastic-anisotropic-plastic (EAP) module (Wilchinsky and Feltham, 2006; Tsamados et al., 2013). In general, most studies show that the simulated sea ice extent interannual variability coincides well with the satellite observations in the Arctic. However, the simulation performance is primarily evaluated on a large global scale (i.e., regarding the central Arctic and all of the marginal seas as a unified whole) and focuses mainly on the ice extent rather than the ice thickness (e.g., Hunke and Bitz, 2009; Flocco et al., 2012; Wang and Su, 2015; Wu et al., 2015; Urrego-Blanco et al., 2016; Chu et al., 2019). Few studies have explored the simulation capability on the regional and local scales (i.e., separate marginal seas, such as the Kara Sea and the Laptev Sea). Different thermodynamic coupling selections with different physical parameter sets can result in significant differences in the simulation (Hunke, 2010; Wu et al., 2015; Urrego-Blanco et al., 2016). Therefore, in marginal and regional seas, such as the SKS, the local model configuration should be thoroughly tested, given that we particularly emphasize the intensity of level ice thickness for ice-resistant offshore structures.

This study constitutes the first quantification of the level ice thickness and associated concept on ice-resistant floating platforms across the SKS. Three critical topics are studied. First, we quantitatively explore configurations in the CICE model to achieve better ice simulation at the local SKS scale. Second, we supplement a 30-year level ice thickness database from 1980 to 2009 and determine its spatiotemporal characteristics. Third, for floating platform design and operation, we propose a reasonable ice regime grade classification criterion and estimate hostile ice conditions at longer return periods based on extreme value statistics. We hope that the present work could be a useful reference for a feasibility study of the potential oil-gas exploitation and ice-resistant structural design in the SKS, particularly for researchers with considerable interests in the Arctic but are far from it geographically.

2 CICE Model Setup and Validation 2.1 Ice Thickness Distribution (ITD) Function in the CICE

The CICE defines an ITD function g = g(x, t, h) to simulate the sea ice pack evolution in time and space. In this function, x = x(i, j) denotes the ocean surface grid, i denotes the longitude, j denotes the latitude, t denotes the time, and h denotes the ice thickness. At a given time and location, g(x, h, t)dh is the fractional area covered by ice in the thickness range [h, h + dh]. Thus, the ITD function is numerically regarded as the approximation of the basic thickness evolution (Thorndike et al., 1975), which can be expressed as follows:

$ \frac{{\partial g}}{{\partial t}} = - \nabla \cdot \left({g\vec u} \right) - \frac{\partial }{{\partial h}}\left({fg} \right) + \varphi, $

(1)

where, on the right side of Eq. (1), the first term represents horizontal transport, with ∇ = (∂/∂x, ∂/∂y) and $ \vec u = $$ \vec ui + \vec uj $, updated by the EVP or EAP theory; the middle term represents the ice thermodynamic transport, with the growth rate f; and the final term φ represents the ice mechanical ridging redistribution.

To solve Eq. (1), the CICE discretizes the ice pack into different thickness categories in each grid cell. The number of categories N can vary from 1 to 15, with the upper boundary H_n for each category n determined using the following recursive equation (Lipscomb, 2001):

$ {H_n} = {H_{n - 1}} + {c_1} + {c_2}\left\{ {1 + \tanh \left[ {{c_3}\left({\frac{{n - 1}}{N} - 1} \right)} \right]} \right\}, $

(2)

with H₀ = 0, c₁ = 3/N, c₂ = 15c₁, and c₃ = 3.

2.2 Model Configuration

Sea ice evolution processes are closely dependent on the metocean conditions. Based on previous research (Wang and Su, 2015; Wu et al., 2015; Wang et al., 2020), the following meteorological and oceanographic forcing fields are selected to achieve better simulation performance in the SKS. For atmospheric conditions, the sea air temperature, wind, specific humidity, and precipitation are obtained from the Common Ocean Reference Experiments version 2 dataset (Griffies et al., 2009), whereas the cloud fraction is obtained from the Arctic Ocean Model Intercomparison Project dataset (Hunke and Holland, 2007). For oceanic conditions, the sea surface tilt, current, and thermal flux are obtained from the Community Climate System Model climate dataset (Collins et al., 2006), whereas the sea surface salinity and temperature are obtained from the Polar Science Center Hydrographic Climatology version 3.0 dataset (Steele et al., 2001; Wang and Su, 2015).

The present simulation is conducted based on the inherent gx1 displaced pole grids in the CICE. Five categories of ice thickness, namely, 0– 0.64, 0.64 – 1.39, 1.39 – 2.47, 2.47 – 4.57, and > 4.57 m, are selected (Hunke and Bitz, 2009). The EVP dynamics and mushy thermodynamics are coupled. Moreover, the bubbly brine algorithm for thermal conductivity, linear remapping scheme for advection, and Delta-Eddington method for shortwave are selected. Table 1 lists the primary inputted run parameters. The sea air temperature, specific humidity, precipitation, cloud fraction, thermal flux, sea surface salinity, and temperature are inputted into the mushy thermodynamics module, whereas the wind, current, and sea surface tilt are inputted into the EVP dynamics module. The initial condition is set as a cold start with no ice, given the lack of data. To achieve better model stability, the total run duration is from January 1969 to December 2009 with a 1-h time step. For the ice-resistant design of floating platforms, assessments of the modeled level ice thickness during 1980 – 2009 in the SKS are emphasized (Fig.1).

2.2 Model Validation

Verification analysis is performed between simulated level ice thickness and scattered observations in the SKS. The in situ measurement data are obtained from the Soviet Union's historical Sever airborne and North Pole drifting station program (Borodachev and Shilnikov, 2002; Romanov, 2004). During this program, the level ice thickness was discretely measured by drilling ice holes on the aircraft runway in the Russian Arctic, where the sea ice cover was completely frozen in March, April, and May from 1928 to 1989. First, the age and partial concentration of ice were determined using low-flying aircraft. After landing, three to five measurements of thickness at 150–200 m intervals were made on the runway. Then, the values were averaged and regarded as the level ice thickness in this location. Thus, one data record exists for each location. Herein, data during 1980–1986 in the SKS were selected. Fig.3 shows these scattered in situ locations in March and April. Each location denotes one valid measurement. For comparison, the simulated data were gridded into the observed scattered locations in time and space by applying a bilinear interpolation.

Fig.4 illustrates the comparison between drilled and simulated level ice thicknesses in March and April. The results indicate the good consistency between both data. The bias, mean absolute error, and root-mean-square error (RMSE) are relatively small, whereas Pearson's correlation coefficient (R) can be up to 0.80 and 0.76. Therefore, the simulated level ice thicknesses in these locations are acceptable, and we deduce that the spatiotemporal variations of level ice thickness in the SKS could be approximately reproduced well based on the present CICE configuration. However, the simulation results for April are somewhat higher than the observations. We infer that the simulated thicker sea ice is largely caused by the inputted overestimated low sea air temperatures or underestimated inflowing heat in the thermal flux. However, because of the lack of measured hydrometeorological data in the Arctic region, the reanalysis data for forcing fields in the present CICE run were not fully assimilated to reduce these deviations. Thus, to decrease the inconsistencies between simulations and observations, both the more accurate inputted metocean forcing data and reasonable sea ice parameterizations of the model should be thoroughly explored.

3 Spatiotemporal Variation of Level Ice 3.1 Annual Distribution of Level Ice Thickness

The 30-year hindcast data during 1980 – 2009 in the CICE model are selected to determine the mean annual distribution of level ice thickness in the SKS. This spatial pattern is illustrated in Fig.5.

As shown in Fig.5, the total domain is covered by level ice with varying thickness. In general, the mean annual level ice thickness primarily ranges from 0.5 m to 1.0 m across the SKS. The thinnest level ice (< 0.6 m) is primarily distributed in the southern zone, the median (0.7–0.8 m) is primarily distributed in the central zone, and the thickest level ice (> 0.9 m) is primarily distributed in the northwest zone. Moreover, a relatively thinner zone (< 0.7 m) is observed in the eastern estuary adjacent zone. Overall, the level ice gradually becomes thicker from southwest to northeast with the increasing latitude.

3.2 Seasonal Distribution of Level Ice Thickness

The seasonal mean spatial distributions of the level ice thickness are also analyzed herein. Considering the spatiotemporal characteristics of the ice cover in the Kara Sea (Duan et al., 2019), in this study, January-February-March is defined as winter, April-May-June is defined as spring, July-August-September is defined as summer, and OctoberNovember-December is defined as autumn. These seasonal patterns are shown in Fig.6.

As shown in Fig.6, the four seasons exhibit different level ice thicknesses in the SKS. In spring, the level ice is the thickest, i.e., primarily between 1.1 m and 1.5 m from south to north. In winter, the level ice is the second thickest, i.e., mainly between 0.7 m and 1.4 m. In autumn, the level ice thickness is no more than 0.4 m. In summer, the SKS has the relatively thinnest level ice thickness of less than 0.3 m in most regions. Thus, the level ice thickness has an obvious seasonal variation in the SKS.

3.3 Monthly Distribution of Level Ice Thickness

The monthly distribution analysis is informative and indispensable to floating platform design suitable for oil-gas exploitation and can directly determine the operation duration and ice-resistant requirements. The simulation results show that level ice thickness has a typical increasing-decreasing cycle and a noticeable monthly variation in the SKS, indicating that monthly sea ice tends to be 1-year level ice.

In detail, Fig.7 plots the contours of the monthly spatial patterns of the mean level ice thickness in the SKS from January to December. In terms of the temporal variability, the entire domain is covered with thicker level ice in the winter and spring months, whereas large ice-free zones occur in the SKS in the summer and autumn months. The level ice approaches the thickest, i.e., primarily between 1.2 m and 1.5 m, in May. In April, June, and March, the level ice exceeds 1 m in most areas. February and January are the two other months that have relatively thicker level ice, i.e., between 0.6 m and 1.0 m. Meanwhile, July with level ice thickness of 0.6 m is the transitional month when ice begins to break up. In August, the thickness even decreases to 0.1 m. Things are quite different in September and October. During both months, the SKS is characterized by large open waters. Then, in November, seawater refreezes again but cannot form 0.4 m level ice thickness. In December, the level ice thickness mainly varies from 0.3 m to 0.7 m. The monthly variations of the simulated level ice thickness in the present CICE configuration also exhibit good consistency with the calculated values from an earlier thermohaline and ice model established by the Eco-System Company during 1990–1995 (Arhipov et al., 1996). In this model, the maximum level ice thickness was equal to the observed 1.4 m in May at Baidara Bay, which is located on the southeastern part of the SKS. Therefore, we can conclude that September and October are the two optimal months for oil-gas exploitation without severe ice threats. To improve economic efficiency, floating platform operations may be extended in several ice-covered months.

4 Grade Divisions and Extreme Ice Conditions 4.1 Grade Divisions of Level Ice

Based on the monthly distribution analysis, 1-year level ice can be up to 1.5 m in the SKS. The World Meteorological Organization (WMO, 1985) categorizes the climatologic ice types according to the ice formation and development stages as follows:

1) Nilas ice, with an ice thickness of less than 0.1 m;

2) Young ice, with an ice thickness of 0.1–0.3 m;

3) Thin ice, with an ice thickness of 0.3–0.7 m;

4) Medium ice, with an ice thickness of 0.7–1.2 m;

5) Thick ice, with an ice thickness of more than 1.2 m.

The ISO 19906 code (International Organization for Standardization, 2010) also employs this classification criterion (hereafter WMO criterion) for the Arctic 1-year level ice. Thus, based on this grade division method, level ice types in each month across the SKS can be summarized as shown in Table 2. These results are acceptable for climatologic analysis. However, this grade division method does not completely consider the engineering characteristics of the floating platform, such as the installation, operation, movability, and associated oil-gas transportation and storage. Therefore, by combining the WMO criterion with the spatiotemporal ice concentrations (Duan et al., 2019), we propose a more reasonable grade division method of level ice for floating platforms in the SKS. As shown in Table 3, ice condition grades are categorized into five groups, namely, excellent, good, moderate, severe, and catastrophic. For floating platforms in the SKS, the first two grades are best, the moderate grade is acceptable if with good ice-resistant protections, and the last two grades are deemed as threatening and disastrous cases.

4.2 Extremal Analysis Method

Floating platforms and other marine structures may suffer from extreme sea ice loads but should survive. Thus, extreme ice conditions at longer return periods must be predicted. In engineering practice, one can fit the extreme value distribution function to a set of maxima derived from subsets of the sample data (e.g., daily, monthly, seasonal, annual maxima) (Coles, 2001). Several mathematical modeling functions, such as Gumbel (Type Ⅰ), Fréchet (Type Ⅱ), and Weibull (Type Ⅲ), have been considered. Theoretically, these three distributions are special cases of the generalized extreme value (GEV) distribution (e.g., Coles, 2001).

The cumulative GEV function can be expressed as follows:

$ {F_{{\text{GEV}}}}\left(x \right) = \exp \left\{ { - {{\left[ {1 + \xi \left({\frac{{x - a}}{b}} \right)} \right]}^{ - \frac{1}{\xi }}}} \right\}, $

(3)

where a, b, and ξ denote the location, scale, and shape parameters, respectively. Herein, the maximum likelihood is utilized to estimate these three parameters. The corresponding equations can be expressed as follows:

$ \left\{ \begin{array}{l} - \frac{{\partial {\text{ln}}L}}{{\partial a}} = \frac{Q}{b} = 0 \hfill \\ - \frac{{\partial {\text{ln}}L}}{{\partial b}} = \frac{1}{b}\frac{{P + Q}}{\xi } = 0 \hfill \\ - \frac{{\partial {\text{ln}}L}}{{\partial \xi }} = \frac{1}{\xi }\left({R - \frac{{P + Q}}{\xi }} \right) = 0 \hfill \end{array} \right., $

(4)

and

$ \left\{ \begin{array}{l} P{\text{ = }}n - \sum\limits_{i = 1}^n {\exp \left({ - {y_i}} \right)} \hfill \\ Q = \sum\limits_{i = 1}^n {\exp \left({ - {y_i} + \xi {y_i}} \right)} - \left({1 - \xi } \right)\sum\limits_{i = 1}^n {\exp \left({\xi {y_i}} \right)} \hfill \\ R = n - \sum\limits_{i = 1}^n {{y_i}} + \sum\limits_{i = 1}^n {{y_i}\exp \left({ - {y_i}} \right)} \hfill \\ {y_i} = \frac{1}{\xi }{\text{1n}}\left[ {1 - \xi \left({\frac{{x - a}}{b}} \right)} \right] \hfill \end{array} \right., $

(5)

where 1nL denotes the likelihood function and (x₁, y₁), (x₂, y₂), ···, (x_n, y_n) denote the level ice thickness samples. These equations can be solved by numerical iterations. With the obtained values of a, b, and ξ, the extreme level ice thickness higher than the observed maxima can be estimated.

4.3 Hostile Ice Condition Estimation

Using the aforementioned grade division method, the extreme level ice thicknesses in November and December across the SKS are estimated based on the monthly maxima of the 30-year simulated ice database.

Fig.8 plots the contours of the spatial patterns of the extreme conditions (i.e., 25, 50, and 100 yr) for level ice thickness in the SKS in November and December. The contours seem to have similar patterns but have different values. Regarding the 25 yr case, level ice thickness primarily ranges from 0.2 m to 0.6 m in November and from 0.5 m to 0.8 m in December. Regarding the 50 yr case, level ice thickness mainly ranges from 0.3 m to 0.6 m in November and from 0.5 m to 0.9 m in December. Regarding the 100 yr case, level ice thickness predominantly ranges from 0.4 m to 0.7 m in November and from 0.6 m to 1.0 m in December. Specifically, the level ice thickness at a return period of 100 yr can be classified as severe. This hostile ice regime analysis could be a reference for the ice-resistant design of floating platforms and affiliated marine structures in the SKS.

5 Conclusions and Future Work

Floating platforms are suitable for oil-gas exploitation in the seasonally ice-covered SKS. However, long-term successive level ice thickness measurements across the entire SKS are still lacking. To fill this gap, the thermodynamic CICE model coupling the EVP and mushy modules is established at the local SKS scale for the first time. Compared with the available measurements in March and April in 1980 – 1986, the RMSE of simulated level ice thickness is no more than 0.2 m. Based on the hindcast data, the spatiotemporal variations of level ice thickness are analyzed annually, seasonally, and monthly. In general, the ice regime belongs to the 1-year level ice and shows an obvious monthly cycle. Thicker level ice appears during the spring months (thickest in May, i.e., between 1.2 m and 1.5 m), whereas large ice-free zones occur in September and October. For floating platforms, a new ice grade criterion with five classifications, namely, excellent, good, moderate, severe, and catastrophic, is proposed. The first two grades are optimal for offshore activities, particularly from August to October, and the moderate grade is acceptable if with ice-resistant protections. Furthermore, hostile ice conditions are estimated based on the GEV distribution. At a return period of 100 yr, level ice is primarily from 0.6 m to 1.0 m in December.

The present investigation could be a useful reference for a feasibility study of the oil-gas exploitation and ice-resistant structural design in the SKS. Furthermore, the sea ice engineering characteristics, such as ridges or bergy bits, should be considered. These ice features threaten the safety of Arctic offshore structures. The CICE model can simulate the sea ice ridging process. Regrettably, the currently available field data cannot adequately and reliably validate the simulated ridged ice thickness across the entire SKS. More fieldwork should be conducted in the future.

Acknowledgements

The study is supported by the National Key Research and Development Program of China (No. 2016YFC0303401), the National Natural Science Foundation of China (No. 51779236), and the National Natural Science Foundation of China–Shandong Joint Fund (No. U1706226).