Snow effects on brash ice and level ice growth

Abstract

Brash ice formation and accumulation occur at a faster rate in ship channels, harbours and turning areas compared to the surrounding level ice. Accurate prediction of brash ice thickness plays an important role in addressing operational challenges and optimisation of ice management strategies. This study enhances existing brash ice growth models by considering the effects of snow and accounting for brash ice expulsion towards the sides of ship channels at each passage. To validate the influence of these critical factors on brash ice thickness, three distinct ship channels located in the Bay of Bothnia, Luleå, Sweden, were investigated. For two test channels formed for study purposes, the slower growth rate of brash ice caused by snow insulation was more prominent than the brash ice growth acceleration caused by the snow–slush–snow ice transformation. In the third channel characterised by frequent navigation, the transformation of slush into snow ice played a more substantial role than snow insulation. In both test channels, the brash ice growth model performed optimally, assuming a 10% expulsion of brash ice sideways at each vessel passage. In the third, wider and more frequently navigated channel, a 1.2% brash ice expelling coefficient predicted well the measured brash ice thicknesses.

1 Introduction

Brash ice forms because of regular navigation in ice-infested waters. This phenomenon is particularly common in ship lanes that are frequently used for navigation within land-fast ice. Two types of ship channels can be differentiated by their formation phase, as discussed in studies by Riska and others (2019) and Zhaka and others (2022a). These two distinct formation processes depend on the navigation frequency, and when the first ship passage occurred. The first formation process occurs when ships consistently follow the same route throughout the year, navigating in thin ice at the beginning of winter. This results in the formation of slush, which, through subsequent freezing and breaking cycles, evolves into more rounded ice pieces. The second formation process is initiated by ship navigation through level ice that has already grown in thickness. The first ship passage will form ice floes. Subsequent recurrent freezing and breaking processes then contribute to the creation of rounded brash ice pieces. Following the initial breaking events, the ship channel becomes partially filled with ice pieces and partially with water (Ashton, 1974). After several passages, the brash ice channels will approach a near-constant water fraction at the surface. Huang (1988) suggested that the time required to attain this steady value primarily depends on the navigation frequency and freezing air temperatures. It can be inferred that when the surface becomes covered with brash ice pieces, the formation phase ends, and the main brash ice consolidation phase starts. The final phase involves the melting of brash ice, which commences in spring or summer due to an increase in solar radiation (Riska and others, 2019).

After each ship passage, part of the brash ice remains in the channel, and part is moved sideways and piled under the level ice to form the side ridges. Previous studies have explored the impact of vessel geometry and speed on the process of brash ice expulsion along the channel's sides (Kitazawa and Ettema, 1985; Ettema and others, 1998). Between two ship passages, the brash ice consolidates within the water-filled and slush-filled voids, resulting in the development of columnar and granular ice (Zhaka and others, 2023a). Several key factors influence the consolidation process of brash ice. These factors include the navigation frequency, cumulative freezing air temperatures experienced between ship passages, macroporosity, lateral ice movement, incoming snowfall and solar radiation (Sandkvist, 1980, 1986; Huang, 1988; Ettema and Huang, 1990; Riska and others, 2014, 2019). Earlier brash ice growth models were primarily founded on heat conduction principles and the ice growth assumptions established for level ice by Stefan (1889). An illustration of this is the simple model devised by Sandkvist (1980; 1986), which assumes that, at each breaking event, all brash ice remains within a ship channel, with its width equivalent to the vessel's beam. This model solely considers the cumulative freezing air temperatures (θ) between ship passages and utilises an empirically derived growth coefficient (α), which was validated using two full-scale test channels located in the Bay of Bothnia, Luleå. Ashton (1974) accounted in his model a constant open water fraction. Another approach was taken by Ettema and Huang (1990), who considered a constant macroporosity and a fixed expelling coefficient. The expelling coefficient described the fraction of ice that was expelled sideways at each ship passage. It is worth noting, however, that this latter model was not validated using full-scale ship channel results. The most recent brash ice growth model (Riska and others, 2019), introduced a more comprehensive approach. Before each ship's passage, the model incorporated three distinct layers: the freeboard, the consolidated brash ice layer, and the unconsolidated brash ice layer. An initial macroporosity was also factored into the model, which was subsequently assumed to decrease after each breaking event due to thermal inertia effects. It can be noted that this model did not include the snow effect, the lateral movement of brash ice, or radiation effects.

The incoming snow has two effects on the brash ice growth as highlighted in Zhaka and others (2022a, 2022b, 2023a). Between two ship passages, the snow accumulates on the brash ice surface insulating the ice and decreasing the brash ice growth rate. Thus, one effect of snow is the brash ice growth retardation between passages. Following a ship passage, this accumulated snow is submerged in water, creating slush within the brash ice macropores. This slush partially undergoes melting if the water temperature is above the freezing point, and the remaining slush will freeze into snow ice (Zhaka and others, 2023a; 2023b). Laboratory investigations on the freezing of slush and water have shown that snow ice develops more rapidly than congelation ice, mainly due to the porous nature of slush (Zhaka and others, 2023b). This process indicates that the slush content within the brash ice may increase the consolidation rate. Therefore, the second effect of snow on brash ice is the acceleration of growth, induced by the variation in freezing rates between freezing slush and freezing water within the brash ice macropores. A study employed image analysis of thin ice sections to quantify the snow–ice content in 200 mm cores sampled along the cross-section of ship channels (Zhaka and others, 2023a). The results indicated that in a channel navigated nine times in total, the snow–ice content constituted 21% of the total ice volume sampled in the middle of the channel and 58% in the side ridges. In another frequently navigated channel, the microstructure of which was studied during winter 2023, the corresponding figures were 41 and 43% for the middle and side ridges, respectively. These observations highlight the significant influence of snow on various components of the brash ice channel.

The main objective of this study was to investigate the two effects of snow on the brash ice growth in land-fast ice ship channels, expanding the scope of the model proposed by Riska and others (2019) by introducing a fourth layer of snow. The new model takes into consideration the insulation caused by snow between two ship passages, as well as the transformation of snow to slush and finally to snow ice following each ship passage. To assess the relative significance of the effects of snow – the insulative effect and the growth acceleration effect – on brash ice growth, a comparison is made between model results that incorporate both effects, one considering only the snow insulative effect, and the original model (Riska and others, 2019) which disregards the effects of snow. These comparative analyses offer insights into the effect of snowfall on brash ice growth. The model was validated using field data from three full-scale channels studied during the winters of 2020–21 and 2021–22 in the Bay of Bothnia, Luleå, Sweden.

The second objective of this study was to examine the lateral movement of brash ice and to estimate the fraction of brash ice that is expelled from the channel at each ship passage. An expelling coefficient is determined for all three channels, and an envelope function is then fitted to establish the correlation between this coefficient and the number of breaking events. This envelope function is subsequently integrated into the four-layer model and its validity is tested against the full-scale measurements. Furthermore, to estimate the expelling coefficient at each vessel passage for all three channels, various constant values for the expelling coefficient are considered within the four-layer model, and the model results are then validated using full-scale thickness measurements. The study also investigates the growth of level ice adjacent to the ship channel. The examination includes an analysis of how two different calculation approaches for the snow–slush transformation impact the estimation of level ice growth.

An overview of the field study is given in the following section. Sections 3 and 4 describe the theoretical aspects of the level ice and brash ice growth. Section 5 presents and discusses the outcomes of this study.

2 Field study

The field investigation was conducted in the land-fast ice of the Bay of Bothnia, Luleå, Sweden, during the winters of 2020–21 and 2021–22. The evolution of brash ice was investigated in three ship channels within distances ranging from 160 to 400 m from the shoreline. Two channels were intentionally created by the icebreaker Tug Viscaria (Luleå port), for research purposes and are hereafter referred to as test channels (TCh01 and TCh02). The third channel is established annually and serves as a primary route for navigation by tugs, merchant vessels, and icebreakers, and is hereafter referred to as the main channel (MCh). The study's geographic location is shown in Figure 1.

Figure 1. Modified optical (Sentinel) satellite images of the Bay of Bothnia and Luleå archipelago. (a) The research site in Luleå is indicated by an orange mark and the SMHI meteorological station is shown with a black mark. (b) Ice conditions in the Bay of Bothnia, Swedish coast, on 2022-03-22. (c) The test and main channels are denoted by orange and purple location markers, while the level ice measurement location is marked in blue.

A comprehensive description of the study site, navigation schedule, measurements and meteorological parameters can be found in Zhaka and others (2023a, 2023c). The air temperature (T _A) and snow thickness (H _S) were measured at the Swedish meteorological station at the Luleå airport and used in this study as input parameters in the level ice and brash ice growth predictions. These are illustrated in Figure 2. In addition, this figure illustrates the calculated cumulative freezing air temperatures (θ).

Figure 2. (a) Air temperature (T _A), and (b) snow thickness (H _S) measured during winters 2020–21 and 2021–22 at the SMHI meteorological station in Luleå. (c) Cumulative freezing air temperatures (θ).

2.1 Thickness measurements

The test channels (TCh01 and TCh02) were investigated in two winters and were positioned at the same location, ~200 m from the shore. These channels were formed by breaking an initial layer of level ice with an approximate thickness of 15 cm. The main channel, on the other hand, was only studied during the second winter and was positioned parallel to the test channel, ~400 m from the shore, as indicated in Figure 1c. This navigation route is used all year round, and the brash ice channel was formed in early winter from navigation through thin ice. The development of level ice at the shoreside of the test channels was examined in both years. The test channels were drilled every metre along a selected cross-section. At each location, measurements were taken for the overall brash ice thickness, side ridges, and the level ice in the channel's vicinity. Additionally, the dimensions of water-filled or slush-filled macro pores, the freeboard, as well as the thicknesses of snow and slush were measured. Cross-section measurements were carried out after the brash ice was partly consolidated.

The first channel was navigated 9 times, and 18 cross-sections were investigated in total. The second test channel was navigated 10 times, and 9 cross-sections were measured. The main channel was navigated 116 times until the final measurement, and the channel's cross-section was measured three times. The specific dates and times of the breaking events and corresponding measurements are detailed in an earlier study (Zhaka and others, 2023a, 2023c). For each cross-section that was measured, several parameters were calculated. These include the average thickness of level ice at the vicinity of the channel, the total equivalent thickness of brash ice which assumes that all the brash ice remains in a ship channel, the average brash ice thickness remaining in the channel, and the equivalent thickness of the side ridges. Illustrative examples of cross-section profiles for TCh02 after the second and tenth breaking events are presented in Figure 3. These cross-sections demonstrate the progression in the geometry of TCh02. For example, the thickness, width and macroporosity increased from the second to the tenth breaking event.

Figure 3. Cross-section profiles of the second test channel (TCh02): (a) after the 2nd breaking event (BE) and (b) after the 10th breaking event. The thickness of brash ice blocks above and below the WL is represented by dark grey bars, while the thickness of snow is illustrated with lighter grey bars. Macropores' dimensions are depicted using blue bars. The snow top is delineated by a solid cyan line, while the freeboard and the bottom of the brash ice are outlined with solid black lines.

On the shoreside of both test channels, measurements were taken of various level ice parameters, including level ice, snow ice, congelation ice, snow and slush thicknesses, as well as the freeboard. During the first winter, level ice measurements were conducted six times, within a surface region specified by a grid pattern with 30 m spacing along the channel and 20 m perpendicular to the shore. In the second winter, level ice measurements were carried out 12 times, covering a distance of 60 m along the channel and 40 m towards the shore.

3 Level ice growth model

The growth of level ice adjacent to the ship channels is usually estimated based on Stefan's (1889) analytical ice growth formulation. In this study, the primary objective was to investigate the impact of two different approaches for calculating snow–slush transformation resulting from flooding, on the growth of snow ice and level ice. To achieve this, a simple model was used to consider the snow–slush transformation while omitting factors such as thermal inertia, water heat flux, incoming radiation, and superimposed ice formation. The first approach used to calculate the snow–slush transformation is based on studies by Leppäranta (1983; 1993), Saloranta (2000), Cheng and others (2014) and Zhaka and others (2022a). The second approach to the snow–slush transformation relies on an empirical correlation established between the water flooding the ice/snow surface and the resulting slush thickness (Zhaka and others, 2023a).

The models check initially if snow is present on level ice and thereafter compare the snow load with the buoyancy of the ice. If there is no snow or the snow mass does not exceed the buoyancy of the ice then the ice growth is estimated according to Ashton's (1986, 1989, 2011) formulation, which considers the snow–air or ice–air coupling as:

(1)$$\displaystyle{{{\rm d}H_{\rm I}} \over {{\rm d}t}} = \displaystyle{1 \over {\;\rho _{\rm i}L_{\rm i}}}\cdot \displaystyle{{T_{\rm F}-T_{\rm A}} \over {H_{\rm I}/k_{\rm i} + H_{\rm S}/k_{\rm s} + 1/h_{\rm a}}},$$

where T _F is the freezing temperature of fresh water, H _S is the snow thickness and H _I is the thickness of the level ice. The snow and ice layer conductivities k _i and k _s, densities ρ _i and ρ _s, latent heat of freezing L _i and heat convective transfer coefficient h _a used in the current study are summarised in Table 1 given in the Appendix.

If the snow mass is larger than the buoyancy of ice i.e.

(2)$$\rho _{\rm s}H_{\rm s}\;> \;( {\rho_{\rm w}-\rho_{\rm i}} ) H_{\rm I}, \;$$

then snow/ice interface flooding is possible and consequently slush is assumed to form instantly at the snow/ice interface. The thickness of the slush can be calculated using the first approach, which neglects the capillarity of snow and considers a mass balance between the ice–slush–snow layers. Assuming that the reduction in the snow's mass due to slush formation is equal to the increase in buoyancy, the mass conservation equation is:

(3)$$\rho _{\rm s}H_{\rm S}-\rho _{\rm s}H_{{\rm SL}} = ( {\rho_{\rm w}-\rho_{\rm i}} ) H_{\rm I} + ( {\rho_{\rm w}-\rho_{{\rm sl}}} ) H_{{\rm SL}}.$$

where ρ _w and ρ _sl are the water and slush densities and H _SL is the slush thickness that forms due to flooding and is equal to:

(4)$$H_{{\rm SL}} = \displaystyle{{H_{\rm S}\rho _{\rm s}-H_{\rm I}( {\rho_{\rm w}-\rho_{\rm i}} ) } \over {\rho _{\rm s} + \rho _{\rm w}-\rho _{{\rm sl}}}}.$$

In the second approach, the slush thickness is estimated through a linear regression correlation obtained by the distance from the ice top to the water level (WL) and the corresponding slush thickness formed in different flooding events primarily next to the second test channel (Zhaka and others, 2023a). The slush thickness increased linearly with the increase in the WL, and the linear regression that described best this correlation is as follows:

(5)$$H_{{\rm SL}} = 0.0386 + 1.0452{\rm WL}.$$

The slush layer was on average 4 cm thicker than the WL on the ice surface, indicating the presence of a slush layer formed from snow capillarity. In this approach, snow capillarity is taken into consideration, and similar to the first approach, it is assumed that the slush forms instantly when snow mass exceeds the buoyancy of ice.

To determine the slush thickness from Eqn (5), the WL or freeboard at each flooding is initially calculated. The freeboard (H _Fb) is calculated as the difference between the ice thickness and the water draft:

(6)$$H_{{\rm Fb}} = H_{\rm I}-\;\left[{\;\displaystyle{{H_{\rm S}\rho_{\rm s}} \over {\rho_{\rm w}}} + \displaystyle{{H_{\rm I}\rho_{\rm i}} \over {\rho_{\rm w}}}} \right].$$

If the freeboard results in a negative value, the WL is calculated as the absolute value of the freeboard and is substituted into Eqn 5 to determine H _SL. Then the snow–ice formation at the snow/ice interface is calculated as follows:

(7)$$\displaystyle{{{\rm d}H_{{\rm SI}}} \over {{\rm d}t}} = \displaystyle{1 \over {v_{\rm w}\rho _{{\rm si}}L_{\rm i}}}\cdot \displaystyle{{T_{\rm F}-T_{\rm A}} \over {H_{{\rm SI}}/k_{{\rm si}} + H_{\rm S}/k_{\rm s} + 1/h_{\rm a}}}, \;$$

where H _SI is the snow–ice thickness, ρ _si is the density of snow ice, k _si is the conductivity of snow ice, v _w is the water content in slush. The maximum water content in slush can be calculated assuming that all the pores in snow are replaced by water when submerged (Adolphs, 1998). However, field experiments reveal that not all the air pores in snow are replaced by water (Weeks and Lee, 1958), some air pores remain in slush. The water fraction in slush is influenced by the snow's initial density, water temperature and the amount of air entrapped during snow submergence (Zhaka and others, 2022b, 2023b). Assuming a snow density of 250 kg m⁻³, the maximum porosity will be equal to 0.73. This value can also be considered the maximum slush porosity. For these calculations, the water content in slush is assumed to be 0.5. This implies that in a slush layer with a maximum porosity of 0.73, 0.5 of these pores are filled with water, and 0.23 with air.

It is assumed that the bottom ice growth does not occur simultaneously with the growth of snow ice. The transformation of slush into snow ice is calculated first, and it terminates when the thickness of the snow ice layer (H _SI) matches the slush thickness. During this transformation, the slush melting, and compaction are not considered. Once all the slush has frozen, it is then assumed that ice continues to grow, with the bottom ice growth estimated as:

(8)$$\displaystyle{{{\rm d}H_{\rm I}} \over {{\rm d}t}} = \displaystyle{1 \over {\;\rho _{\rm i}L_{\rm i}}}\cdot \displaystyle{{T_{\rm F}-T_{\rm A}} \over {H_{\rm I}/k_{\rm i} + H_{{\rm SI}}/k_{{\rm si}} + H_{\rm S}/k_{\rm s} + 1/h_{\rm a}}}.$$

4 Brash ice growth model

This section gives an overview of the brash ice growth model which is further detailed in the following subsections. Before a breaking event, the brash ice is assumed to be divided into four layers (Fig. 4) which are the following. The top layer consists of snow accumulated between two ship passages. The second layer is dry brash ice above the WL formed due to buoyancy, and its pores are filled with air. This layer is assumed to have a constant porosity and an unchanged thickness between ship passages. The third and fourth layers are the consolidated and wet brash ice, respectively. The consolidated layer continuously increases between ship passages while the wet brash ice thickness reaches zero if ice gets fully consolidated. The thickness of each layer before the jth ship passage are denoted with H _S(j−1), H _D(j−1), H _C(j−1) and H _W(j−1), for snow, dry, consolidated and wet brash ice, respectively. Figure 4 provides a visual representation of the initial two calculation steps of the model.

Figure 4. (a) Brash ice thickness prior to a breaking event is divided into four layers, snow (H _S(j−1)) which has a porosity p _s; the dry layer (H _D(j−1)) above the WL with air-filled pores and a porosity p ₀; the consolidated brash ice (H _C(j−1)); and the water-filled wet brash (H _W(j−1)) with a porosity p_{j −1}. After each breaking incident, the brash ice is assumed to have a constant layer porosity p ₀, a new brash ice layer (H _B) and a uniform temperature (T _av). The snow submerges in the pores forming slush-filled water pores (SL) on the top layer. The bottom layer consists of water-filled pores (W). (b) The second step shows the redistribution of ice above and below the WL and the porosity change from p ₀ to p_j.

In the first step, a ship passage breaks all four layers and forms a mixture of ice water and slush, which has a new equivalent thickness H _B, a porosity p ₀ and an average temperature T _av. After the breaking event, snow transforms into slush within the macropores between the ice pieces. For modelling purposes, it is assumed that the slush migrates to the top of the wet brash ice layer as shown with blue background in Figure 4. The water temperature in the wet brash ice layer is assumed to always be at the freezing point. Consequently, no melting of slush is considered. Therefore, the thickness of the brash ice layer, where the macro voids are filled with slush instead of water, will depend on the initial snow thickness submerged in water and the brash ice porosity (p ₀). Here it is assumed that p ₀ instantly after each ship passage is equal to 0.25. This value was estimated from the measured macroporosity of the ship channels (see Zhaka and others, 2023c). Thus, the brash ice layer occupied by slush-filled voids will have a thickness equal to the ratio between the submerged snow thickness and the porosity (H _S/p ₀). However, if the slush layer exceeds the total brash layer, all the slush beneath the brash ice is assumed to melt instantly. A laboratory study conducted by Zhaka and others (2023b) observed that only 35% of the initial snow thickness that was submerged in water transformed into snow ice, highlighting the variability of the transformation process in different conditions. Thus, to account for different scenarios, the model is also validated for cases where the slush formed in brash ice is 30, 50 or 70% of the initial snow thickness.

In the second step, the model assumes the redistribution of brash ice below and above the WL. The fraction of brash ice pieces floating above the WL and the fraction of ice pieces that become submerged in the water are determined based on buoyancy principles (Riska and others, 2019). The porosity of brash ice above the WL is assumed to remain constant and similar to the initial porosity of brash ice after each breaking event (p ₀ = 25%) since there is no water to freeze within the pores. On the other hand, the porosity of brash ice below the WL is assumed to change instantly from p ₀ to p _j due to ice formation driven by the temperature difference between the cold ice pieces and the freezing water temperature. It is assumed that the slush completely submerges in the macropores of wet brash ice.

In the third step, the model estimates the brash ice consolidation, which occurs within the slush-filled or water-filled macropores. When accounting for slush freezing instead of water freezing, the water content (v _w) in the slush is considered. Assuming that the air pores in snow are completely substituted by water, then the maximum water content can be estimated as the ratio of snow density (ρ _s) to pure ice density (ρ _i) (Adolphs, 1998). For modelling purposes, ρ _s was assumed constant and equal to 250 kg m⁻³, based on density measurements during both winters. The effect that v _w has on the brash ice growth was tested by varying this parameter from a value of 0.3 to its maximum value. Additionally, the insulative effect of the incoming snow between two ship passages is considered. It is assumed that the snow remains as a uniform layer above the dry brash ice layer and will not submerge in water until the next breaking event. This assumption can be reasonable when the surface of the channel is covered completely with dry brash ice and the open water fraction on the surface is negligible. A schematic illustration of the final step is given in Figure 5.

Figure 5. Brash ice consolidation between two consecutive ship passages. The continuity equations for the heat flow at each interface are numbered 1–4.

4.1 Brash ice accumulation

The snow–slush transformation does not contribute to the overall thickness of the brash ice after a breaking event. Instead, the slush fills the wet macropores between the submerged ice pieces. As a result, the equivalent brash ice thickness along the ship channel after the breaking event is determined by the mass conservation of the ice before and after the ship passage (Riska and others, 2019). The equivalent brash ice thickness (H _B) after the breaking incident can be expressed as follows:

(9)$$H_{\rm B} = H_{{\rm D}( {\,j-1} ) } + \displaystyle{{H_{{\rm C}( {\,j-1} ) }} \over {( {1-p_0} ) }} + H_{{\rm W}( {\,j-1} ) }\displaystyle{{( {1-p_{( {\,j-1} ) }} ) } \over {( {1-p_0} ) }}.$$

Assuming that the average temperature of brash ice (T _av) is reached instantly, it can be calculated using the mass and heat conservation equations before and after the breaking event, taking into account the contributions from all four layers. This average temperature (T _av) can be expressed as:

(10)$$T_{{\rm av}} = \displaystyle{ { \eqalign{& { ( {T_{{\rm SA}} \! + \! T_{{\rm DS}}/2} ) ( {1 \!- \!p_{\rm s}} ) H_{{\rm S}( {\,j-1} ) }} { + ( {T_{{\rm SI}} \! + \! T_{{\rm DS}}/2} ) ( {1-p_0} ) H_{{\rm D}( {\,j-1} ) }} \cr& {+ ( {T_{\rm F} + T_{{\rm SI}}/2} ) H_{{\rm C}( {\,j-1} ) } + T_{{\rm F\;}}( {1-p_{\,j-1}} ) H_{{\rm W}( {\,j-1} ) }} }} \over {( {1-p_0} ) H_{{\rm D}( {\,j-1} ) } + H_{{\rm C}( {\,j-1} ) } + ( {1-p_{\,j-1}} ) H_{{\rm W}( {\,j-1} ) }}}, \;$$

where T _SA, T _DS and T _SI are the snow/air, dry brash ice/snow and dry brash ice/consolidated brash ice interface temperatures.

In the second step of the calculation, the initial porosity (p ₀) of wet brash ice is assumed to instantaneously change to a new porosity value p _i. The reduction in porosity occurs due to the freezing of the water fraction in the slush-filled voids or the water-filled voids, driven by the temperature gradient between the water and ice. If no snow is present, as described by Riska and others (2019), then the porosity can be expressed as:

(11)$$p_j = p_{0\;}-\displaystyle{{c_{\rm i}\Delta T} \over {L_{\rm i}}}( {1-p_{0}} ) ,$$

where ΔT is the temperature difference until equilibrium (T _F − T _av), and c _i is the specific heat capacity of ice, see Table 1 in Appendix.

If snow is present, and the slush formed after the ship passage completely fills the voids within the entire brash ice layer, then the energy and mass balance equations between the slush-filled pores and brash ice pieces, as well as the final porosity equation, can be expressed as:

(12)$$H_{\rm W}\rho _{\rm i}L_{\rm i}( {v_{\rm w}p_{0\;}-p_{\rm j}} ) = \rho _{\rm i}c_{\rm i}\Delta TH_{\rm W}( {1-p_{0}} ) ,$$

(13)$$p_j = v_{\rm w}p_{0\;}-\displaystyle{{c_{\rm i}\Delta T} \over {L_{\rm i}}}( {1-p_{0}} ).$$

If the voids within the brash ice layer are partially filled with slush and partially with water, then the thickness of brash ice with slush-filled voids is H _BSL = H _S/p ₀. The energy and mass conservation equation between the slush-filled and water-filled voids with the brash ice pieces, along with the derived porosity equation, are as follows:

(14)$$\begin{align}\displaystyle{{H_{\rm S}} \over {\,p_{0\;}}}&\rho _{\rm i}L_{\rm i}( {v_{\rm w}p_{0\;}-p_{\rm j}} ) + \left({H_{\rm W}-\displaystyle{{H_{\rm S}} \over {\,p_{0\;}}}} \right)\rho _{\rm i}L_{\rm i}( {\,p_{0\;}-p_{\rm j}} ) \cr & = \rho _{\rm i}c_{\rm i}\Delta TH_{\rm W}( {1-p_{0}} ) , \;\end{align}$$

(15)$$p_{\rm j} = \displaystyle{{H_{\rm S}( {v_{\rm w}-1} ) + H_{\rm W}( {\,p_{0\;}-c_{\rm i}\Delta T( 1-p_{0\;}) /L_{\rm i}} ) } \over {H_{\rm W}}}.$$

4.2 Brash ice consolidation

The third step of the model calculates the brash ice consolidation between two ship passages. Initially the model tests the presence of slush and estimates the thickness of brash ice occupied by slush (H _BSL = H _S/p ₀). If slush is present, the brash ice consolidation will occur only in the slush filled brash ice macropores and in particular in the water fraction of the slush.

If the entire wet brash ice column is filled with slush, then the freezing process will occur within the water fraction v _w of the slush-filled macropores until the wet brash ice fully consolidates, thereafter only water freezing is considered. In case that brash ice is filled partly with slush and partly with water, it is assumed that consolidation occurs initially in the H _BSL layer and then in the rest of wet brash ice.

The snow thickness that covers the ship channels between breaking incidents is calculated by the incoming snow thickness measured by the SMHI meterological station. The brash ice growth equation is derived by substituting the interface temperatures obtained from the first three continuous heat flux equations and integrating the fourth equation. These equations are illustrated in Figure 5. The brash ice growth equation is:

(16)$$H_{\rm C} = \sqrt {k_{{\rm bi}}^2 h_{\rm e}^2 -\displaystyle{{2k_{\rm bi}} \over {\rho _{{\rm bi}}p_jv_{\rm w}L_{\rm i}}}\theta } -k_{{\rm bi}}h_{\rm e}.$$

where p_jv _w accounts for the freezing of slush within the slush-filled macropores, ρ _bi is the brash ice density, v _w is the not considered for the part of wet brash ice that is solely filled with water. Additionally, p _j is assumed equal to one if the brash ice is fully consolidated and h _e is given as:

(17)$$h_{\rm e} = \displaystyle{{H_{\rm D}} \over {k_{\rm d}}}\;+\;\displaystyle{{H_{\rm S}} \over {k_{\rm S}}}\;+\;\displaystyle{1 \over {h_{\rm a}}}, \;$$

where k _d and k _s are the heat conductivities of the dry brash ice and snow, h _a is the surface convective heat transfer coefficient.

4.3 Expelling coefficient

The estimated brash ice equivalent thickness (H _B) can be compared with the brash ice equivalent thickness (H _eq) observed from all the measured cross-sections of the ship channel's (TCh01, TCh02, MCh). This equivalent brash ice thickness (H _eq) can be estimated as the total surface area of the ice accumulated on the channel (A _BI) and both side ridges (A _R) divided by the width of the vessel as:

(18)$$H_{{\rm eq}} = \displaystyle{{A_{{\rm BI}} + A_{\rm R}} \over {W_{{\rm BICh}}}}$$

A schematic illustration of the channels' cross-section and the surface areas of brash ice and side ridges, as well as a schematic illustration of the brash ice equivalent thickness are provided in Figures 6a, b. The channel's width (W _BICh) was considered to be 12 m for the TChs, which represents the beam of tug Viscaria. The main channel (MCh) was navigated by several vessels with beams varying between 12 and 31 m. Tug Viscaria would often navigate in the sides of this channel and not at the centre. Thus, the expelling coefficient for all three cross-section profiles of the MCh was calculated considering a constant width of 26 m representing the main navigation lane in the centre line. Also, the actual channel widths measured in each cross-section were used, these widths were 36, 30 and 42 m after 14, 73 and 116 BE, respectively.

Figure 6. (a) Typical cross-section of a ship channel assuming a width equal to the vessel's beam (W _BICh). The surface area covered by brash ice and the areas occupied by the side ridges are denoted as A _BI and A _R. (b) The equivalent brash ice thickness. A scheme of the first brash ice growth model (BIGM1) that does not consider the loss of ice sideways. (c) The average brash ice thickness that remains in the ship channel and the equivalent side ridge thickness. A scheme of the brash ice growth model that considers the ice loss sideways (BIGM2).

The trapezoid method was used to calculate the cross-sectional areas of brash ice and ridges (excluding the level ice thickness) for each cross-section measured in all three channels. If the measured thicknesses at each hole are H _i and the distance between these w _i then the surface area is equal to:

(19)$$A_{{\rm R/BI}} = \mathop \sum \nolimits^{} 0.5\;( {H_{\rm i} + H_{{\rm i\ + \ 1}}} ) w_{\rm i}.$$

H _i represents the thicknesses measured in each drilled hole along the channels' cross-section, while w _i is the distance between the drilled holes. In the case of the test channel, the thicknesses were measured every metre while in the main channel distances of 2, 2.5 and 5 m were used.

At each vessel passage, part of the brash ice is expelled sideways and piled under the level ice layer, forming side ridges on both sides of the channel, see Figure 6. The fraction of ice that has been moved sideways is calculated for each measured cross-section in all three ship channels assuming the following procedure. Initially, the average thickness of brash ice that remains in the channel is calculated as the surface area of the ice cross-section inside the channel divided by the channel's width, which is assumed to be equal to the vessel's beam. Thereafter, the equivalent thickness of the ice that is expelled sideways (H _Req) is estimated as the ratio between the surface area of ice cross-section in both side ridges with the channel's width. Then the cumulative expelling coefficient (χ_j) is calculated for each measured cross-section as:

(20)$$\chi _j = \displaystyle{{H_{{\rm Req}}} \over {H_{{\rm eq}}}}.$$

This cumulative expelling coefficient χ_j represents the total fraction of ice that is expelled sideways from the first ship passage until the jth ship passage when the cross-section measurements were conducted. The difference χ = χ_j − χ_j−₁ will give the fraction of ice that is expelled from the channel at each ship passage. It should be noted that the calculations of χ_j are based on the measured brash ice cross-sections, which were taken after the brash ice had undergone partial consolidation. As a result, the potential influence of the ice formation around the blocks on the value of χ_j at the time of measurement is neglected. χ_j can be included in the brash ice growth model discussed in the earlier section. The average thickness of the brash ice (H _BI) that remains in the channel and the side ridge equivalent thickness (H _R) can be estimated as follows:

(21)$$H_{{\rm BI}} = ( {1-\chi_j} ) H_{\rm B}, \;\;$$

(22)$$H_{\rm R} = \chi _jH_{\rm B}.$$

A schematic illustration of this concept is given in Figure 6c. For modelling purposes, the ridge porosities are assumed equal to the brash ice porosities. The analysis of the measured cross-sections showed that the maximum initial porosity, considering both brash ice and side ridges, was found to be about 25%, which is the value used in this study (see Zhaka and others, 2023c).

The brash ice and ridge thicknesses estimated by the model (H _BI and H _R) can be validated with the results obtained from all three ship channels. For the test channels, it is assumed that the equivalent ridge thickness represents the amount of ice that is expelled from the channel and accumulates in a width equal to the vessels beam (12 m). However, similarly to the equivalent brash ice thickness that represents the total volume of brash ice, the equivalent ridge thickness calculated from the model and the cross-sections of all three channels represent the volume of ice expelled from the channel.

In addition, to using the cumulative expelling coefficient χ_j to estimate the brash ice and side ridge thicknesses, the fraction of brash ice that is expelled at each ship passage (χ) is assumed constant and is incorporated in the model to continuously estimate the thickness of the brash ice that remain in the channel and the equivalent ridge thickness at each ship passage.

5 Model results and discussion

The estimated thicknesses of level ice (H _I) and snow ice (H _SI) using two different approaches for slush thickness calculation are validated against thicknesses measured on the shoreside of the test channels. The brash ice model's performance, considering the effects of snow, is evaluated with varying cumulative air freezing temperatures, snow–slush transformation rates and different water contents in the slush. The estimated equivalent brash ice thickness, considering the effects of snow, is compared with results from the original model (Riska and others, 2019), and the model that only considered the snow insulative effect. The model's performance when considering the snow effects and the expelling coefficient, is evaluated considering the cumulative expelling coefficient estimated from the measured ship channels, as well as various constant values of χ. Finally, the model results are validated using data obtained from all three ship channels.

5.1 Level ice adjacent to test channels

Figure 7a illustrates the measured and estimated thicknesses of snow, level ice and snow ice on the shoreside of TCh01. The snow–ice thickness adjacent to TCh01 was measured four times, while the thicknesses used to validate the model calculations were measured during each cross-section profile assessment. Figure 7b illustrates the measured and estimated thicknesses of snow, level ice and snow ice on the shoreside of TCh02. The level ice model best predicted the thicknesses for a surface heat convective coefficient of 10 W m⁻¹ °C⁻² in the first winter and 20 W m⁻¹ °C⁻² in the second winter. The air temperatures and snow thicknesses used as model input are given in Figure 2.

Figure 7. Measured and estimated thicknesses of level ice (H _I), snow ice (H _SI), and snow (H _S) adjacent to (a) the first test channel (2020–21) and (b) the second test channel (2021–22). The solid and dashed lines depict the model estimations where the snow–slush transformation was calculated using Eqns (4) and (5), respectively.

In the first winter, the predicted H _SI and H _S were in good agreement with the measured values for both approaches used to estimate the slush formation, the empirical linear regression and the mass conservation principles. A better agreement was achieved for the H _I when the slush thickness was calculated based on the mass conservation (Eqn 4), while using the empirical correlation between the H _SL and WL overestimated the level ice thickness from the 66th day, corresponding to 15 of February and thereafter.

During the second winter, the measured thicknesses of snow ice, and snow, similarly to the first winter were well predicted using both approaches, while the total level ice thickness was best predicted when H_SL was calculated from the empirical regression Eqn 5. On the other hand, when the slush thickness was estimated using Eqn 4, it led to an underestimation of the total ice thickness.

Using the empirical regression function (Eqn 5) resulted in slightly higher snow–ice thickness compared to the mass conservation approach. This difference started with the initial snowfall that occurred at the beginning of winter when the ice was thin. For example, in the second winter, the first 4 cm of incoming snow that covered the thin ice, based on the mass conservation approach, formed a slush layer of 1.6 mm. In contrast, the empirical regression function assumes that, on average, 4 cm of snow transforms into slush due to snow capillarity. Therefore, in this approach, all the initial snow thickness transforms into slush and snow ice. The difference in this initial snow–ice thickness between model approaches caused a difference in the total level ice thickness, which, in turn, affected the difference in snow–ice thickness between the two approaches, thereafter. With a higher level ice thickness and lower snow thickness (see the dashed lines in Fig. 7), there is a lower imbalance between the snow and ice layers, and therefore, the snow–ice thickness difference between the two approaches reduced over time.

5.2 Brash ice growth

The equivalent brash ice thickness was simulated using three different models: the brash ice growth model (BIGM) including two effects of snow, the model that considers only the snow insulative effect, and the original model which omits the snow effects or assumes a zero snow thickness. The brash ice thickness was simulated for a total of 150 days, using three levels of constant air temperatures, −5, −10 and −20°C, corresponding to a θ of 750, 1500 and 3000°C days. The total number of breaking events was 38 passages with a navigation frequency of one passage every 4 d. The snow thickness was assumed to increase linearly with time until it reached a total thickness of 0.75 m after 150 days. The effect that various ship frequencies have on the BIGM was assessed for a total of 298, 149, 75 and 38 ship passages. Finally, the effect that snow–slush transformation fraction and slush water content have on the brash ice equivalent thickness were assessed for a snow–slush transformation fraction varying from 0.3 to 1.0 and water contents in slush varying between 0.3 and 0.73.

5.2.1 Effects of snow

The equivalent brash ice thickness (H _B) obtained from the model that considers the snow effects for a total snow thickness of 0.75 m, is compared with the model results when only the insulative effect is considered for the same snow thickness input, and with the model that neglects the snow thickness, see Figure 8. The convective heat transfer coefficient used in all three models is equal to 10 W m⁻¹ °C⁻¹. The brash ice thickness was estimated for 38 ship passages and three different cumulative freezing air temperatures.

Figure 8. (a) Brash ice equivalent thickness (H _B) estimated from the BIGM which includes both effects of snow (H _SL), only the snow insulative effect (H _S), and the original model for a zero-snow thickness (H _Z). (b) The brash ice macroporosity (p _i) change simulated from the model that considers both effects of snow (P _SL), and the one considering only the snow insulative effect (P _S).

The slower growth rate caused by snow insulation is offset by slush freezing once the snow is submerged in water. The insulation effect of snow is more prominent than the increase in thickness due to the quicker growth of slush-filled pores compared to water-filled pores. The maximum difference in the total simulated brash ice thickness between bare ice and a scenario with a total snow thickness of 0.75 m is 10 cm for 750°C days, 24 cm for 1500°C days, and 45 cm for 3000°C days.

The brash ice thickness is lower when only the insulative effect is considered in the BIGM, compared to when both effects of snow are considered. The difference in H _B between the two models increases as the cumulative freezing air temperatures (θ) decreases. This indicates that in milder winters with relatively high snow thickness, neglecting the slush–snow–ice transformation but including the snow insulative effect leads to higher errors in estimating the brash ice thickness. The maximum differences in H _B between models were 20, 17 and 11 cm for 750, 1500 and 3000°C days, respectively.

A higher initial reduction in macroporosity is obtained when the snow effects are considered (see Fig. 8b), as the slush freezing till temperature equilibrium is quicker than freezing water for the same initial porosity (Zhaka and others, 2023b). The porosity reduction increases as the cumulative freezing air temperatures increases when only the insulative effect is considered but decreases when the slush–snow ice transformation is also included. When only the snow insulative effect is considered, the porosity is higher for lower θ due to a lower growth rate. However, when the slush–snow–ice transformation is considered, the same slush content occupied a larger fraction of macropores for lower H _B compared to H _B caused by higher θ. As a result, there is a higher porosity reduction with lower θ scenarios, when the column of brash ice has a greater proportion of slush.

5.2.2 Frequency of navigation, water fraction in slush, snow–slush transformation fraction

The snow effects BIGM is assessed for different frequencies of ship navigation, various slush water fractions and snow–slush transformation fractions. A θ of 1500°C days and a total snow thickness of 0.75 m is considered in all the simulations. Figure 9a shows the maximum brash ice thickness reached after 150 days when a channel is navigated twice a day (2BEs/day), once a day (1BE/day), once in two days (1BE/2 days) and once in four days (1BE/4 days). All the input parameters were consistent, with the only variation being the number of ship passages. The increase in the number of breaking events has resulted in a higher amount of ice accumulating in the channel, as brash ice consolidates more rapidly when it is frequently broken and mixed. Moreover, the frequent submergence of snow, forming slush and subsequently snow ice may have influenced the overall process. According to the model estimates, the snow accumulation on the ship channel between two passages ranges from 0.005 m to 0.02 m for a navigation frequency between 2BEs/day and 1BE/4 days. The brash ice thickness was 1.5 m higher when the ship channel was navigated twice a day compared to when it was navigated once every four days.

Figure 9. (a) Brash ice equivalent thickness (H _B) and (b) macroporosity variation for different frequencies of navigation. (c) The equivalent brash ice thickness (H _B) was estimated assuming various water content (v _w) in the slush (0.3–0.73) and different snow–slush transformation rates, e.g. H _SL = 0.7*H _S.

Figure 9b illustrates that the initial porosity reduction is the same for all the breaking frequencies, as the first passage occurs at the same time, and the submerged snow thickness is equal. However, with an increasing number of ship passages, the reduction in porosity until thermal equilibrium decreases. This is because the cumulative freezing air temperatures between ship passages decrease, and the temperature gradient between water and ice decreases as well. The blue solid line (2BE/day) in this figure consists of similar porosity reduction steps to the other lines given in the same figure. However, due to a high number of breaking events, these steps are closer, and thus, the line gives the impression of a filled surface area. Additionally, the snow accumulated between ship passages decreases with the increase in the number of breaking events (Figure 9b). This figure demonstrates that as the number of breaking events increases, less ice forms due to the temperature difference between the ice and water until thermal equilibrium is reached.

Figure 9c presents the maximum equivalent brash ice thickness after 150 days under varying water contents in slush (0.3, 0.4, 0.5, 0.6 and 0.73, shown as blue bars). In this scenario, a 100% snow–slush transformation is assumed, meaning that all submerged snow transforms into slush. The figure also includes model results using different snow–slush transformation fractions 0.3, 0.5 and 0.7, meaning that 30, 50 and 70%, of the submerged snow transformed to slush (represented by the orange bars). The maximum difference in equivalent brash ice thickness (H _B) between water fraction inputs of 0.3 and 0.73 was 31 cm. While comparing the scenario where H _SL is assumed equal to H _S and when it is assumed equal to 30% of the submerged H _S, the maximum difference in the brash ice equivalent thickness between the two cases is 15 cm.

5.3 Model validation

The brash ice equivalent thickness was estimated using three models: the one considering the two effects of snow, the one considering only the insulative effect of snow and the original model (Riska and others, 2019) that does not include the snow layer. These results were validated using data obtained from three ship channels located in the Bay of Bothnia, Luleå, Sweden. The same convective heat coefficient, equal to 10 W m⁻¹ °C⁻¹, was used in all models, except for the estimation of the equivalent brash ice thickness in the main channel, where the models performed better for a coefficient of 15 W m⁻¹ °C⁻¹.

5.3.1 Test channel 2020–2021

The equivalent brash ice thickness was calculated for each measured cross-section assuming a channel width equal to the vessel's beam, 12 m (Sandkvist, 1980; 1986). The predicted and observed equivalent thicknesses of brash ice in TCh01 are shown in Figure 10a.

Figure 10. (a) Observed and estimated equivalent brash ice thickness (H _B). (b) The brash ice macroporosity estimated with all three models and the snow thickness (H _S) accumulated on the channel between two ship passages.

The original model gave the highest brash ice thickness, while the snow's insulative effect model gave the lowest thickness. The value of H _B estimated with the model that considers the two effects of snow lies between the other two model results, as the quicker slush freezing compared to water freezing compensates for the reduction of the growth rate due to snow insulation. The maximum H _B estimated by the model considering the effects of snow was 26 cm lower than the original model and 10 cm greater than the snow insulative effect model. The initial amount of snow that was submerged at the first ship passage was 37 cm, and the corresponding reduction in macroporosity was 15%, as shown in Figure 10b. The macroporosity change was lower in the other two models, where snow–slush–snow ice transformation was neglected.

5.3.2 Test and main channels 2021–2022

The measured and estimated equivalent brash ice thicknesses for the second test channel (TCh02) and main channel (MCH) are illustrated in Figure 11. In TCh02, like TCh01, the results from the model which considered the two effects of the snow situated between the original model and snow insulative effect model predictions, with maximum differences of 30 and 14 cm. The predicted maximum H _B were 1.25 m for the two snow effects models, 1.55 m for the original model and 1.11 m for the snow insulative effect model.

Figure 11. Observed and predicted equivalent brash ice thicknesses (a) in the second test channel (TCh02) and (b) in the main channel (MCh).

However, this was not the case for the MCh, where H _B estimated from the two effects snow model was higher than those given by the other two models. This implies that in the main channel, the snow–slush–snow ice transformation had a greater impact on the equivalent brash ice thickness than the insulating effect of snow. This conclusion aligns with the observation that incoming snowfall in frequently navigated channels submerges quickly and does not insulate the ice growth for a long time, as seen in Figure 12. Omitting the snow's insulative effect only leads to a difference of 6 cm in the maximum thickness. The maximum difference in H _B between the snow effects model and the snow insulative effect model is 19 cm.

Figure 12. Brash ice macroporosity (p_j) and snow thickness (H _S) accumulated on the channels between two ship passages in the second test channel (a) and the main channel (b).

The equivalent brash ice thickness in the MCh was calculated considering a constant width, which was concluded to be equal to the centre track of the main channel (MCh cte). Also, the thickness was calculated considering the actual width of the channel (MCh ac) as measured. A wider channel results in a lower equivalent brash ice thickness. However, the H _B was better represented by the model for a constant width, equal to the central track of the channel.

The maximum porosity reduction in the test and main channels exceeds 8%, as shown in Figure 12. The porosity reduction estimated by the snow effects BIGM was higher than in the other two models for both channels. The maximum amount of snow submerged in TCh02 was 20 cm, while in MCh was 16 cm. This suggests that even in frequently navigated channels, the amount of snowfall that submerges at each passage can be significant.

5.4 BIGM including expelling coefficient

In this section, the cumulative expelling coefficient (χ_j) which represents the fraction of brash ice expelled sideways from the first ship passage until the jth ship passage is calculated for each measured cross-section of all three ship channels. An asymptotic envelope function was introduced to fit and represent the measurements. The asymptotic function was incorporated in the snow effects model. In addition, various levels of the ice fraction that is expelled from a channel at every ship passage (χ) were included in the BIGM model. Thus, the model's performance is assessed considering χ_j and various assumed constant values of χ. Finally, the estimated average thickness of brash ice (H _BI) that remains in the channel and the equivalent thickness of ridges (H _R) were validated against the results obtained from all three channels.

5.4.1 Expelling coefficient

The cumulative expelling coefficient (χ_j), calculated by Eqn 20 for all ship channels, and its progression with the number of breaking events is given in Figure 13. For TCh01, χ_j demonstrated a linear increase with the rising number of breaking events. However, for TCh02, the linear regression displayed a weak correlation. This can be attributed to the fact that the thickness of brash ice and ridges along the ship channel exhibit significant variability.

Figure 13. Calculated cumulative brash ice expelling coefficient (χ_j) for the test and main channels (TCh01, TCh02, MCh). The blue solid line illustrates the asymptotic envelope function (Eqn 23).

For the main channel, χ_j was calculated considering a constant width, representing the track at the centre of the channel, but was also calculated considering the actual channel's width. The results show that χ_j decreases with the increase of the channel width. This implies that in a wide channel, the fraction of ice that remains in the channel can be higher compared to a narrower channel. However, the value of χ_j will simultaneously depend on the position that the ship navigates along the channel and the ship's beam relative to the channel width. For example, if the ship navigates in the centre, the brash ice will be pushed to the sides but remain within the channel, but if the ship navigates on one side of the channel, the fraction of ice that is expelled out of the channel will increase.

Furthermore, an asymptotic envelope function was fitted to encompass all the calculated χ_j for each measured cross-section across all three ship channels. The data fit best with the asymptotic function, which can be expressed as:

(23)$$\chi _j = 0.58-0.6\;{\rm e}^{{-}0.3j\;}.$$

This function assumes that χ_j initially increases with the number of breaking events and subsequently levels off to a constant value of 0.58. This coefficient estimates the cumulative amount of ice expelled from the first breaking event (j = 1) to a certain number of BE (jth ship passage). From the cumulative expelling coefficient (χ_j), the expelling coefficient (χ) at each passage can be calculated as the difference in the cumulative expelling coefficient between two consecutive passages. The average expelling coefficient for the first 10 breaking events (BE) was ~0.09. Note that negative values of χ, which were obtained as a consequence of thickness variabilities along the channel, were neglected in these calculations.

5.4.2 Model estimations

The cumulative expelling coefficient (Eqn 23) for four navigation frequencies, and six various constant expelling coefficients (χ) from 5 to 50% were incorporated in the brash ice growth model that considers the two effects of snow (snow insulation and snow–slush–snow ice transformation), see Figure 14. According to the asymptotic function when χ_j approaches a constant value, the H _BI will account for 42% of the overall equivalent ice thickness, while the side ridges will constitute 58%. When the frequency of navigation increases from 1BE/4 days (38 BEs) to 2BEs/day (298 BEs), both the H _BI and H _R increase by a factor of 1.6. The model's results for χ ranging from 5% to 50% indicate that if the expelled ice after each ship passage is between 30 and 50%, H _BI after 5 BE remains nearly constant, ranging from 13 to 20 cm, while the ridge thicknesses continue to increase.

Figure 14. (a) Estimated average brash ice thickness (H _BI) and (b) ridge equivalent thickness (H _R) using the snow effects model that considered ice expelled sideways using Eqn (23), for four different navigation frequencies (i.e. H _ex38BE). Also using six different constant expelling coefficients (χ) varying from 5% to 50% (i.e. H ₁₀).

5.4.3 Model validation

The measured and estimated thicknesses of the brash ice that remains in the channel and the equivalent thicknesses of side ridges for all three channels are given in Figures 15 and 16. In the test channels, the corresponding thicknesses are estimated using the snow effects model for two constant expelling coefficients (χ) 0.1 and 0.2, and the cumulative expelling coefficient (χ_j), while for the main channel, the constant values of χ used are 0.01, 0.012 and 0.015. A convective heat transfer coefficient of 10 W m⁻¹ °C⁻² was employed for all three channels, as the model considering the snow effects performed best at this value.

Figure 15. Measured and estimated average thicknesses of brash ice (H _BI) and the equivalent ridge thicknesses (H _R) for the first test channel (TCh01). Two constant expelling coefficients of 10% and 20%, also the cumulative expelling coefficient (e.g. H _BIex) determined by Eqn (23) were used as input.

Figure 16. (a) Measured and estimated average thicknesses of brash ice (H _BI) and (b) the equivalent ridge thicknesses (H _R) for the second test channel (TCh02) and the main channel (MCh). The constant expelling coefficients (χ) used were 10 and 20% for the TCh01, also 1.0, 1.2 and 1.5% for the MCh. The cumulative expelling coefficient (χ _j) determined by Eqn 23 was also used as input (i.e. Tch_ex).

The BIGM predicts best the H _BI and H _R in both test channels for a constant expelling coefficient of 0.1, which is similar to χ calculated for the first 10BE of the test channels. In a previous study focused on brash ice formation in model-scale ship channels, the expelling coefficient calculated after each vessel passage was generally observed to be 0.1 or lower (Ettema and Huang, 1990). Employing the exponential function to account for the sideways movement of brash ice underestimated the average brash ice thickness and overestimated the amount of ice expelled sideways.

The brash ice thickness that remained in the main channel was not significantly affected by the channel width. This observation was made when calculating the brash ice thickness using two approaches: one considering a constant channel width (MCh_cte), representing the mid-lane of the channel, and the other considering the actual channels (MCh_ac). Refer to Figure 16a for details. This is because the increase in the amount of brash ice remaining in the channel was offset by the increase in the channel's width. The calculated equivalent ridge thickness was affected significantly when the channel width increased as the amount of ice expelled sideways decreased, see Figure 16b.

For the main channel, the BIGM performed well in predicting the average brash ice thickness that remains in the main channel using a expelling coefficient of 0.012. This constant coefficient implies that only 1.2% of the brash ice was expelled sideways at each ship passage. On the other hand, the exponential function underestimated the average brash ice thickness but was in good agreement with the ridge thickness, when assuming a constant channel width. However, it is important to note that this evaluation is limited to three measured cross-section profiles, which may not fully represent the actual development of the main channel. Also, the H_R calculated from the actual channel width are not well represented by the model.

The estimated brash ice thickness that was expelled from the test channels at each ship passage was about 10%, whereas in the main channel, ~1.2%. This difference between the two channel types is likely related to their width. While both test channels had a width comparable to the vessel's beam, the main channel had a width up to 42 m. As a result, in each ship passage in the test channel, the ice piled under the level ice, forming the side ridges. In contrast, in the main channel, ice was pushed to the sides but still remained within the channel, forming piles of ridged ice within the channel. Also, the main channel was navigated about 11 times more than the test channels till the last measurement.

Conclusions

This paper describes the effects of snow in the brash ice growth and the brash ice expulsion sideways in three ship channels. Two channels were established solely for research purposes during the winters of 2020–21 and 2021–22, and one regularly navigated channel occurring annually in Luleå was also studied. In addition, the level ice on the shoreside of both test channels was investigated.

The slush that formed because of flooding on the shoreside of the test channel was estimated using two different methods. The first method involved utilising the principle of mass conservation among the ice, slush and snow layers. The second method was based on an empirical linear correlation established between the WL and the corresponding slush thickness, which considered the capillarity of snow. Both approaches resulted in good predictions of the snow–ice thickness. However, employing the empirical linear regression, which considered the snow capillarity, yielded a slightly higher snow ice thickness. Also, it overestimated the level ice thickness next to the first test channel after the 15 of February but predicted well the level ice thickness next to the second test channel.

For both test channels, the estimated equivalent brash ice thickness (H _B) using the two snow effects model was lower than that obtained from the original model and higher than the estimation derived from the snow insulative effect model. In the case of the test channels, the snow insulation had a more significant influence compared to the accelerated consolidation resulting from the freezing of slush. For both channels, the decrease in brash ice equivalent thickness caused by the snow insulation was partially offset by the increase in brash ice thickness attributed to the freezing of slush.

In the main channel, the brash ice thickness estimated considering the effects of snow in the brash ice growth model was found to be 15 cm more than that estimated by the original model and 20 cm more than the snow insulative effect model. This suggests that the transformation of slush into snow ice played a significant role compared to the snow insulative effect, which was less pronounced as the snow remains on the brash ice for only a brief period.

Furthermore, it was observed that the difference in the estimated brash ice equivalent thickness between the snow insulative effect model and the snow effects model increased as the cumulative freezing air temperatures decreased. Consequently, in milder winters characterised by relatively high snow thickness, neglecting the formation of snow ice will lead to a higher underestimation of the brash ice thickness.

The analysis of brash ice sideways movement showed that the cumulative expelling coefficient, counting for the total ice expelled sideways from the first ship passage until the observed time, exhibited an increase in the expelling coefficient with the increase of breaking events. Thereafter, it reached a maximum value of 0.58, which remained unchanged despite the continued increase in breaking events. In both test channels, for a total of 10 breaking events, it was found that, on average, after each vessel passage, about 10% of brash ice was expelled sideways. In the main channel, calculations showed that cumulative amounts of brash ice expelled sideways decreased with the increase in the width of the channel.

The snow effects model provides the most accurate prediction for both the average brash ice thickness that remains within the channel after a ship passage and the equivalent side ridge thickness in both test channels, considering a consistent expelling coefficient of 10%. Additionally, this model performs well in the main channel when a constant expelling coefficient of 1.2% is applied. This difference between the two types of channels is mainly attributed to their width and frequency of navigation. In a narrow channel such as the test channels, the ice accumulates under the level ice while in a wide channel like the main channel the brash ice remains within the channel track after a ship passage. On the other hand, utilising the exponential function to accommodate the sideways movement of brash ice resulted in an underestimation of the average brash ice thickness in all three channels and an overestimation of the lateral ice expulsion in both test channels.