2. Experiments

The experiments are performed in a rectangular cuboid RB cell (height $H=0.12$ m, aspect ratio $\varGamma = L/H = 2$, width–height ratio $W/H = 0.5$, i.e. a so-called quasi-two-dimensional system). An expansion vessel is attached to compensate for the volume change during the freezing process and so to monitor the mushy phase growth; see figure 1(a). The top and bottom plate temperatures are measured with thermistors and controlled with circulating baths. Sidewalls are isolated with foam material to approach ideal adiabatic boundary conditions. The cell is initially filled with degassed salty water of initial salinity $S_i$. Here, we use sodium chloride (NaCl), the principal solute of sea water (Lyman & Fleming Reference Lyman and Fleming1940), as the salt, and $S_i$ is varied from 0 % (fresh water) to 3.5 % (average ocean salinity), with the corresponding freezing point $T_{0i}$ (Hall, Sterner & Bodnar Reference Hall, Sterner and Bodnar1988) varying from $0^{\circ }$C to $-2.1\,^{\circ }$C. The (hot) bottom plate temperature $T_b$, which is equal to the initial liquid temperature, ranges from $0.4^{\circ }$C to 12.9 $^{\circ }$C, and the (cold) top plate temperature $T_t$ is set 10 K below the freezing temperature, hence $T_t = T_{0i}-10$ K (from $-10^{\circ }$C to $-12.1\,^{\circ }$C, depending on $S_i$). A dimensionless bottom superheat parameter can then be defined as $\varTheta _i = (T_b-T_{0i})/(T_b-T_t)$. This parameter characterizes the initial ratio of the admissible temperature difference in the liquid phase over the total gap across the cell. It can also evolve slightly during the experiment. This is because the melting temperature depends on the actual salinity, and at the equilibrium state, $S_e$ and so $T_{0e}$ might differ from the initial values (actually $S_e$ can be as high as $7\,\%$, with $T_{0e}\approx -4.2\,^{\circ }$C; see figure 5b). Since the total quantity of solute is constant, one expects an enhancement of salinity in the liquid due to salt ejection from the mush (and so an increase of $\varTheta$).

Figure 1. Sketch of the experimental system and temperature–salinity parameter space. (a) Cartoon of the experimental set-up. The height of the cell is $H=0.12$ m, the aspect ratio is $\varGamma = L/H = 2$, and the width–height ratio is $W/H = 0.5$. An expansion vessel is attached to compensate for the volume change and to monitor the ice/mushy phase growth. The bottom plate temperature is kept at $T_b$, and the top plate temperature $T_t$ is maintained at 10 K below the freezing point $T_{0i}$ at the initial salinity $S_i$. A dimensionless bottom superheat parameter is defined as $\varTheta _i=(T_b-T_{0i})/(T_b-T_t)$. (b) Temperature–salinity parameter space. In this work, we test $\varTheta _i$ values $1/5$ (purple), $1/3$ (red), $3/8$ (green), $3/7$ (pink), $1/2$ (blue) and $3/5$ (brown). The initial salinity $S_i$ varies from 0 % to 3.5 %; $T_b$ (dots) is shown for each case. We also plot $T_{0i}$ (black squares) and $T_t$ (grey triangles) as functions of $S_i$.

In this work, we first test the states with $\varTheta _i$ values $1/3$, $3/8$ and $1/2$ with salinities $S_i \in [0, 3.5]\%$. Cases with $\varTheta _i$ values $1/5$, $3/7$ and $3/5$, and $S_i = 3.5$%, are also performed to test the results in a wider temperature range. A visualization of the temperature–salinity parameter space explored in the experiments is provided in figure 1(b). It reports the temperature values $T_b$, $T_{0i}$ and $T_t$ for different $\varTheta _i$ and $S_i$. In this experimental configuration, the ice/mushy phase grows from the top plate until the system reaches an equilibrium. We assume that this state is attained when the mean mushy phase thickness $h(t)$ (see figure 1a) changes less than 1 % in eight hours, with the maximum experiment duration being three days. The ice or mushy phase thickness and its porosity $\phi (t)$ (i.e. the volume fraction of residual liquid in the ice or mushy phase) can be obtained, respectively, by identifying the ice–liquid interface in photos of the system, and by using the expansion vessel readings together with the assumption of total mass conservation in the system. The evaporation from the expansion vessel is negligible as it causes only a mass loss less than 0.1 % of the total mass of the system. We refer to Appendix A for more details on the measurement protocols.

3. Mushy phase morphology and growth dynamics

Figures 2(a,c) show photos of the experiment at its equilibrium state in two typical cases ($S_i=0$, 3.5 % and $\varTheta _i=3/8$). One can readily appreciate the different thicknesses of the iced domains (thinner in the salty water case) and the different interface shapes between the two cases. However, their dissimilarity is not only volumetric and morphological, but also structural. To illustrate this, we hollowed out a tiny (5 mm) hemispherical cavity on the top horizontal surface of the ice matrix on which we deposited a drop of red dye. Figures 2(b,d) show two lateral view photos of the stain after 5 min to form its deposition. We observe that when $S_i=0$ (figure 2b), the iced phase is transparent and the dye does not spread into it, indicating the dense structure of the matrix. On the contrary, at $S_i=3.5\,\%$ (figure 2d), the iced phase looks opaque and the dye spreads anisotropically, which highlights a complex porous medium filled with liquid (i.e. a mushy phase structure) with internal vertical channels.

Figure 2. Images of the iced layer together with its internal structure, and evolution of the freezing front. (a–d) Photos of ice/mushy phases without (a,b) and with (c,d) porosity. The experimental conditions are (a,b) $S_i=0\,\%$, and (c,d) $S_i=3.5\,\%$, both with $\varTheta _i = 3/8$. In (b,d), the red dye of the same amount and concentration is deposited into a small hemispherical pit carved at the same location on the upper surface of the corresponding solid matrix in (a,c). The ice in (a,b) is transparent, and the dye does not spread. The mushy phase with porosity in (c,d) is relatively opaque, and the dye spreads anisotropically, possibly indicating the presence of channels. (e) Dynamics of the freezing front. The mushy phase thicknesses $h(t)$, non-dimensionalized by $H$, are plotted as functions of time non-dimensionalized by the diffusive time scale $t^* =H^2/\kappa _{ice}=1.22\times 10^4$ s (where $\kappa _{ice}$ is ice thermal diffusivity) for three experimental runs. (f) The same data as in (e) are compensated by the diffusive behaviour $h(t)/\sqrt {t}$ (normalized by the first data to compare different cases).

It is noteworthy that the different iced phase structures correspond to different freezing dynamics. This is well illustrated in figures 2(e,f), where the temporal growth of the freezing front, $h(t)$, is shown for three sets of experimental conditions, corresponding to pure water, intermediate salinity and sea water salinity, all at $\varTheta _i=3/8$. Though all three cases display an initial period of diffusive ice phase growth, i.e. $h \propto \sqrt {t}$, later the growth rate decreases for the case $S_i=0$, 2 %; while when $S_i=3.5\,\%$, the instantaneous growth rate increases and even exceeds the diffusive growth rate before reaching the equilibrium. It is remarkable also to note that the equilibrium height is not simply proportional to the salinity. The above highlighted temporal evolution of the ice growth in a salty solution depends on the internal structure of the mushy phase that in turn might also change in time (e.g. due to the formation of brine channels Wettlaufer et al. Reference Wettlaufer, Worster and Huppert1997; Worster & Wettlaufer Reference Worster and Wettlaufer1997). A quantitative description of this process is challenging because it requires detailed knowledge of the mushy microstructure. The same is true for any attempt to characterize the morphology of the solid–fluid interface. For this reason, in the following we focus on a simpler question: we try to tackle the globally averaged equilibrium state that is reached asymptotically in time. In other words, we aim at understanding the factors that control the average thickness, denoted $h_e$, of the iced mushy phase.

4. Global properties of the equilibrium state

We display in figure 3 two important spatially averaged properties of the mushy phase at the equilibrium (from now on denoted with a subscript $e$): they are the dimensionless thickness $h_e/H$ and the porosity $\phi _e$. Besides this, we report the saturation time (or time to equilibrium) denoted as $t_{e}^{90\,\%}$, here for practical reasons defined as the time it takes to reach 90 % of the equilibrium thickness non-dimensionalized by the diffusive time scale $t^*=H^2/\kappa _{ice}=1.22\times 10^4$ s. The dependence of the thickness $h_e$ as a function of the initial salinity $S_i$ is very different for distinct values of the control parameter $\varTheta _i$ (see figure 3a): $h_e$ exhibits a similar non-monotonic (decrease–increase–decrease) pattern for $\varTheta _i = 1/3, 3/8$, and is monotonically decreasing for $\varTheta _i = 1/2$. These trends are in part reflected on the $\phi _e(S_i)$ dependence, which presents considerable fluctuations for intermediate salinity values, although its overall trend is increasing for all $\varTheta _i$ cases (see figure 3b).

Figure 3. Global properties of the equilibrium state from experiments. The initial superheats are respectively $\varTheta _i$ values $1/3$ (red), $3/8$ (green) and $1/2$ (blue), and the initial salinity $S_i$ varies from 0 % to $3.5\,\%$. The error bars are measured by repeated experiments in two cases with $S_i = 3.5\,\%$ and $\varTheta _i = 1/3$, $1/2$, respectively. (a) The mean mushy phase thickness $h_e$ initially decreases with $S_i$ for all three superheats, and is followed by a sharp increase for the $\varTheta _i = 1/3$ and $3/8$ cases. The shaded areas indicate the minimal and maximal spatial variations of mushy phase thicknesses. The thicknesses are non-dimensionalized by $H$. (b) The mean mushy phase porosity $\phi _e$ exhibits a non-monotonic yet overall increasing trend with $S_i$. (c) Thickness $h_e$ does not exhibit a statistically strong correlation with $\phi _e$. (d) The saturation time $t_e^{90\,\%}$ is defined as the time to reach 90 % of the equilibrium mushy phase thickness $h_e$, and is non-dimensionalized by the diffusive time scale $t^*=1.22\times 10^4\,{\rm s}$. The variations of $t_{e}^{90\,\%}$ with $S_i$ are consistent with those of $h_e$. (e) Thickness $h_e$ shows a clear positive correlation with $t_{e}^{90\,\%}$. A linear fitting gives $h_e/H = 0.0912 \times t_{e}^{90\,\%}/t^*$ (dash-dot line), with correlation coefficient $R^2=0.9846$.

The fact that a more porous ice matrix (high $\phi _e$) can release more salt into the liquid and as a consequence lower the final freezing point of the solution, $T_{0e}$, hence reducing the size of $h_e$, might give a rationale for the globally observed trends. However, this alone is not enough to explain any of the observed fluctuations or non-monotonic dependencies in $h_e$ and $\phi _e$. Moreover, it is important here to notice that the existence of a correlation between $\phi _e$ and $h_e$ is not evident from the present measurements (see figure 3c). If, for instance, mass conservation and geometry were simply ruling the relation between the equilibrium mushy layer thickness and the mushy layer porosity, then one would have $h_e(1-\phi _e) = const.$ (where the constant represents the height of the ice layer in the case of zero salinity and correspondingly no porosity). This implies an increasing trend of $h_e$ with $\phi _e$, while figure 3(c) suggests just the opposite (one can notice a decreasing $h_e(\phi _e)$ for fixed values of $\varTheta _i$). We conclude that the $h_e(\phi _e$) relation is a non-trivial function of the system control parameters $S_i$ and $\varTheta _i$.

Finally, we observe that the trends for $t_{e}^{90\,\%}$ are instead consistent with those of $h_e$ (figure 3d), and a clear linear correlation can be seen between the two quantities (figure 3e). This correspondence can provide a way to estimate the time to equilibrium just by focusing on the magnitude of $h_e$. Therefore, in the following, we will focus on modelling the physical mechanisms responsible for the variations of $h_e$. The question that we wish to address is: can the equilibrium global properties be understood in terms of a model that considers just (i) the role of porosity, and (ii) the role of thermal stratification, on the steady state of the heat transfer process in the system, without taking into account its complex transient dynamics?

5. Multi-layer heat flux model

To answer the above question, we need to understand how $h_e$ and $\phi _e$ are related and how the heat transport across the system is influenced by the convection in the mushy phase and by the combination of stable stratification and unstably stratified turbulence. We will show that a one-dimensional (or integral) multi-layer heat flux model that relies on scaling laws known for natural convection in fluid and porous media can explain $h_e$ with good accuracy once $\phi _e$ is provided as input from the experiment. On the other hand we observe that a quantitative description of $\phi _e$ requires a model for the evolution of the microscopic structure of the porous ice matrix. This is a challenging task that currently goes beyond our capabilities and that we leave for future research. In summary, the multi-layer heat flux model presented in the following allows us to explain the complex relation linking $h_e$ and $\phi _e$, when the latter is considered as an independent variable.

We assume that the position of the solid–fluid interface in the system at the equilibrium $h_e$ is the result of the balance between the mean heat flux across the mushy phase, $\mathcal {F}_m$, and the one in the liquid, $\mathcal {F}_l$ (see figure 4). Both these quantities require a parametrization, which we detail in the following. For conciseness and better readability, we refer the reader to table 1 in Appendix D for the denomination of the many physical quantities and material properties involved in the model.

Figure 4. Synoptic representation of the averaged multi-layer heat flux model. Row 1 provides the adopted parametrization for the global heat flux expression in the mushy phase, $\mathcal {F}_m$, together with the definitions of $Nu_m$ and $Ra_m$. Rows 2–4 report the expression for the global heat flux in the liquid phase, $\mathcal {F}_l$, together with the $Nu_l$ and $Ra_l$ definitions (column 4). Three different combinations of stratifications may exist: (a) unstable, (b) mixed, (c) unstable. The discriminating criterion for their occurrence is the monotonicity of liquid density within the temperature range $[T_{0e}, T_b]$, i.e. the relative magnitude of $T_{0e}$, $T_b$ and the maximum-density temperature $T_{max}$ (column 2), which depends on the salinity $S_e$. Column 3 gives the sketches of these combinations, together with the layer thicknesses and the interface temperatures. The colour gradient indicates the density gradient in the liquid phase, with red representing lighter liquid, and blue representing heavier liquid. The model is supplemented with five input parameters, $T_t$, $T_b$, $S_i$, $\phi _e$ and $H$, and leads to an estimation of $h_e$ (row 5) through the assumption of a statistically steady heat flux balance among layers.

Table 1. Symbols used in this paper.

We express the heat flux across the mushy phase in the statistically steady condition in the form

(5.1)\begin{equation} \mathcal{F}_m = Nu_m \, \varLambda _m(T_m)\,\frac{T_t-T_{0e}}{h_e}. \end{equation}

Equation (5.1) involves the mushy phase Nusselt number $Nu_m$, i.e. the ratio between the total and the conductive heat transfer across the mushy phase. Its dependence as a function of the Rayleigh number $Ra_m$ is known from previous studies on convection in porous media. Here, we use the parametrization by Gupta & Joseph (Reference Gupta and Joseph1973):

(5.2)\begin{equation} Nu_m=\left\{ \begin{array}{@{}ll} 1 , & Ra_m < Ra_{mcr}, \\[2pt] 1.3338+0.0099\,Ra_m, & Ra_m \ge Ra_{mcr}, \end{array} \right . \end{equation}

where the mushy-critical Rayleigh number is $Ra_{mcr}=4{\rm \pi} ^2$. The Rayleigh number for the flow in a porous medium is defined as

(5.3)\begin{equation} Ra_m=\frac{ \varPi \it}{\phi_e}\,\frac{g( \rho(S_t, T_t)- \rho(S_e, T_{0e}))h_e}{\nu(T_m)\,\kappa(T_m)\, \rho(S_e, T_m)}. \end{equation}

This definition involves the evaluation of the permeability $\varPi$, a quantity describing the conductivity of a porous medium to fluid flow (Pal, Joyce & Fleming Reference Pal, Joyce and Fleming2006). The permeability is a function of the structure of the porous ice matrix and is notoriously difficult to measure, in particular for ice, because it is highly sensitive to the environmental conditions and to the history of the medium. Several previous studies, including Kawamura et al. (Reference Kawamura, Ishikawa, Takatsuka, Kojima and Shirasawa2006), Petrich, Langhorne & Sun (Reference Petrich, Langhorne and Sun2006) and Polashenski et al. (Reference Polashenski, Golden, Perovich, Skyllingstad, Arnsten, Stwertka and Wright2017), revealed that there exists a threshold porosity ($\phi _c$) below which the ice behaves as impermeable. Above such a threshold, the permeability is well approximated by a power-law function of the porosity. We adopt here a similar functional parametrization to relate $\varPi$ to $\phi _e$:

(5.4)\begin{equation} \varPi=\left\{ \begin{array}{@{}ll} 0 , & \phi_e < \phi_{cr}, \\[2pt] C(\phi_e - \phi_{cr})^\alpha, & \phi_e \ge \phi_{cr}, \end{array} \right. \end{equation}

where the three adjustable parameters are chosen as $C=7 \times 10^{-8}$ m$^2$, $\alpha =3$ (as at the leading order of the Carmen–Kozeny relation; Wells et al. Reference Wells, Hitchen and Parkinson2019) and $\phi _{cr}=0.054$ to best fit $h_e$ measured experimentally. Although these values are ad hoc for our study, the modelled $\varPi$ is of the same order of magnitude as the ones reported in sea ice laboratory and field researches (Kawamura et al. Reference Kawamura, Ishikawa, Takatsuka, Kojima and Shirasawa2006; Petrich et al. Reference Petrich, Langhorne and Sun2006; Polashenski et al. Reference Polashenski, Golden, Perovich, Skyllingstad, Arnsten, Stwertka and Wright2017).

Similarly, the mean heat flux in the liquid phase can be expressed in the approximate form

(5.5)\begin{equation} \mathcal{F}_l = \left\{ \begin{array}{@{}ll} Nu_l \, \varLambda_l(T_l)\,\dfrac{T_{0e}-T_b}{H-h_e}, & T_{max} \le T_{0e}, \\ \varLambda_l(T_s)\,\dfrac{T_{0e}-T_{max}}{h_{max}} = Nu_l \, \varLambda_l(T_u)\,\dfrac{T_{max}-T_b}{H-h_e-h_{max}}, & T_{max}\in (T_{0e},T_b),\\ \varLambda_l(T_l)\,\dfrac{T_{0e}-T_b}{H-h_e}, & T_{max} \ge T_b. \\ \end{array} \right . \end{equation}

Note that the above expression is composed of three distinct cases according to the value, $T_{max}$, of the temperature corresponding to the maximal density of the liquid solution at the equilibrium with respect to $[T_{0e}, T_b]$. This reflects in three distinct cases: (a) unstable stratification, (b) a stable layer above an unstable one, and (c) a single stable layer; see figure 4. The scale $h_{max}$, appearing in case (b), is the height of the stably stratified layer in the intermediate case. It can be determined by equating the two fluxes across the stable and unstable layers of liquid. The dimensionless heat flux in the liquid is accounted for by the Nusselt number $Nu_l$, dependent on the Rayleigh number $Ra_l$, a relation that is relatively well known. Wang et al. (Reference Wang, Calzavarini, Sun and Toschi2021c) proposed the parametrization below by fitting the data from their simulations on the freezing of fresh water in an RB convection system, which leads to theoretical predictions agreeing well with their experiment results. The same parametrization is adopted in the present study:

(5.6)\begin{equation} Nu_l =\left\{ \begin{array}{@{}ll} 1 , & Ra_l < Ra_{lcr}, \\ 0.12+0.88\,Ra_l / Ra_{lcr}, & Ra_l \in [Ra_{lcr}, 1.23\,Ra_{lcr}], \\ 0.27(Ra_l -Ra_{lcr})^{0.27}, & Ra_l > 1.23\,Ra_{lcr}, \end{array} \right . \end{equation}

where the critical Rayleigh number for the onset of convection is $Ra_{lcr}=1708$. We note that the definition of the Rayleigh number in the liquid is needed only for the two cases that involve a density unstable stratification (cases (a) and (b) in figure 4). We use

(5.7)\begin{equation} Ra_{l} = \left\{ \begin{array}{@{}ll} \dfrac{g( \rho(S_e, T_{0e})- \rho(S_e, T_b))( H-h_e )^3 }{\nu(T_l)\,\kappa(T_l)\,\rho(S_e, T_l)}, & T_{max}\le T_{0e},\\ \dfrac{g( \rho(S_e, T_{max})- \rho(S_e, T_b))( H-h_e-h_{max}) ^3 }{\nu(T_u)\,\kappa(T_u)\,\rho(S_e, T_u)}, & T_{max}\in (T_{0e},T_b). \end{array} \right . \end{equation}

Finally, the model tuning requires the numerical values of material properties (see Appendix C for the empirical formulas adopted), which are evaluated at the equilibrium salinity $S_e$. Given the facts that in the present system no salinity source exists, that liquid in the mushy phase contains only a small salt fraction, and that salt is mixed in the liquid either vigorously by convection or slowly by diffusion, we can assume that $S_e$ is uniform in the liquid layer. Hence $S_e$ can be related implicitly to $h_e$ and the inputs of the model by the mass conservation of the salt:

(5.8)\begin{equation} S_e\,\rho(S_e, T_{mean})\,(H-h_e(1-\phi_e))=S_i\,\rho(S_i, T_{b})\,H. \end{equation}

It is worth noticing that (5.8) is suitable only for the present system (when applied to other systems, e.g. open systems, a measurement of $S_e$ must be supplied). The evaluation of $Ra_m$ requires also the value of $\rho (S_t, T_t)$. While $T_t$ is known as it is imposed externally, $S_t$ is not. We assume here that the latter takes the value such that $T_t$ is a local freezing point, $T_0(S_t)=T_t$; in other words, we assume local thermodynamical equilibrium within the mushy phase (Wells et al. Reference Wells, Hitchen and Parkinson2019). Finally, by solving simultaneously the algebraic equations of heat flux balance $\mathcal {F}_m = \mathcal {F}_l$ and the mass conservation of the salt (5.8), using the measured value $\phi _e$, the model leads to the equilibrium mushy phase thickness $h_e$. We refer to Appendix B for a more detailed explanation of how the model is implemented numerically.

The model results are all illustrated in figure 5. Figure 5(a) shows the mean mushy phase thickness $h_e/H$ as a function of $S_i$ together with its comparison with the experiment results. At $\varTheta _i=1/2$ (blue lines and symbols), where ice thickness decreases for increasing salinity, the agreement with the model is fairly good. The model suggests that in this condition, the heat transfer through the mushy phase remains diffusive (mushy diffusion, MD state). Increasing the salinity, the freezing point is lowered, and this slightly reduces the frozen layer thickness. The role of porosity is marginal as here it stays approximately constant with $S_i$ (see figure 3b). For $\varTheta _i=1/3$ and $3/8$ (red and green lines and symbols), the non-monotonic behaviour is also captured by the model.

Figure 5. Model results for ice/mushy phase thickness compared with the experiments, and identification of heat transport modes. MD/MC stand for diffusive/convective heat transfer in the mushy phase; similarly, LD/LC stand for diffusive/convective heat transport in the liquid phase, while LPC stands for diffusive and convective heat transport due to penetrative convection in the liquid phase, i.e. a combination of diffusive and convective heat transport, as a stable stratification lies above an unstable one. (a) The initial superheat values are respectively $\varTheta _i$ values $1/3$ (red), $3/8$ (green) and $1/2$ (blue), and the initial salinity $S_i$ varies from 0 % to $3.5\,\%$. The calculated $h_e$ (lines with open circles) is close to the experimental results (dots). The spatial variations of mushy phase thicknesses from the experiments (shaded areas) are also included, together with results of reduced models ignoring mushy phase convection ($Nu_m \equiv 1$, dash-dot lines) and ignoring liquid stratification (linear density dependence on temperature $\rho (S_e, T)=\rho (S_e, T_b)\,(1-\beta (S_e, T_b)\,(T-T_b))$, dash-double-dot lines). (b) Phase space diagram based on the two effective control parameters $Ra_m$ and $Ra_l$. The third effective control parameter $S_e$ is not presented as an axis, for readability. Instead, different circle sizes are used to distinguish the different $S_e$. There are in total six heat transport modes of the system. Our experiments fall within five cases. Different circle colours are used to distinguish different $\varTheta _i$. The highest salinity ($S_e \gtrsim 2.7\,\%$) is expected to suppress the LPC regime as it removes the density anomaly of the solution. (c–g) Images from the experiment for typical cases corresponding to the five observed heat transport modes.

It is relevant here to also consider two reduced versions of the model, one ignoring the occurrence of convection in the mushy phase (by assuming $Nu_m \equiv 1$, dash-dot lines in figure 5a), the other neglecting the non-monotonic density stratification in the liquid phase (by assuming linear density dependence on temperature, i.e. $\rho (S_e, T)=\rho (S_e, T_b)\,(1-\beta (S_e, T_b)\,(T-T_b))$, where $\beta (S_e, T_b)$ is the liquid thermal expansion coefficient at the bottom plate, dash-double-dot lines). The comparison of the calculated $h_e$ stemming from the two reduced models allows us to interpret the ‘jump’ or sudden increase in the ice thickness as being associated with the occurrence of convection in the mushy phase (mushy convection, MC state) with increased $\phi _e$ (see figure 3b). Convection strengthens the heat transport, promoting a larger mushy phase growth rate, because $\text {d}h/\text {d}t \propto \mathcal {F}_m - \mathcal {F}_l$. The equilibrium mushy layer thickness, as a result of the competitive balance between $\mathcal {F}_m$ and $\mathcal {F}_l$, also increases because of the enhanced heat transfer efficiency in the mushy phase.

However, the mushy phase convention alone is not sufficient to explain the tendency of $h_e$ (see the dash-double-dot lines in figure 5a). Especially large deviations from experiment results can be seen when $S_i$ is low as they fall into LD (liquid diffusion) and LPC (liquid penetrative convection) regimes (see figure 5b). In the LD state, only diffusive heat transport occurs in the liquid phase, either because the whole liquid phase is stably stratified or because the liquid phase is too thin to support convection. In the LPC state, there is always a stable stratification in contact with the mushy phase above the unstably stratified turbulence, also known as penetrative convection. The presence of this stable layer limits the below convective motion and prevents the mushy phase from a direct penetration of the convection. As an effect, the heat transport across the liquid phase is relatively low, enabling the mushy phase to grow thick. When $S_i$ is increased, the temperature differences decrease in the stably stratified region (from $T_{0e}$ to $T_{max}$) and increase in the unstably stratified turbulence (from $T_{max}$ to $T_b$). As a consequence, the convective motion intensifies, becoming even more turbulent, leading to strong heat transport across the liquid phase and so decreasing $h_e$. Similarly, the temperature difference in the unstable stratification is smaller when $\varTheta _i$ is lower. The unstably stratified turbulence is weaker, and the thermal resistance in the stable stratification above plays a more important role in determining the heat transport of the whole liquid phase. This accounts for the relatively larger errors between the dash-double-dot lines (assuming the entire liquid phase to be unstably stratified) and the experiment results in figure 5(a) when $\varTheta _i= 1/3$ and $3/8$, compared to cases with $\varTheta _i= 1/2$ and the same $S_i$. By further increasing $S_i$, the stable stratification layer disappears, and the heat transport in the liquid phase is purely convective (liquid convection, LC state; see figure 5b).

The remaining discrepancies between the model and the experiments are likely to be related to the modelling of the mushy phase, in particular from neglecting its non-homogeneous internal structure and the uncertainties in the parametrization of its permeability. It is also worth mentioning that the heat transport in the LPC regime can also be clearly described by means of two alternative parameters introduced in Wang et al. (Reference Wang, Reiter, Lohse and Shishkina2021a). They are the so-called density inversion parameter (which here reads $(T_{max}-T_{0e})/(T_b-T_{0e})$) and the Rayleigh number relative to the whole fluid layer. It is an alternative but equivalent approach with respect to the one adopted here.

In summary, there are six possible heat transport modes in the system, depending on whether the heat transport is purely diffusive (MD) or convective (MC) in the mushy phase, and whether the heat transport in the liquid phase is purely diffusive (liquid diffusion, LD), diffusive and convective in distinct superposed layers (LPC), or convective (LC). The occurrence of these states is determined by the magnitude of the two control parameters $Ra_m$ and $Ra_l$. The boundaries between purely diffusive (MD and LD) and convective (MC, LPC and LC) heat transport are the critical Rayleigh numbers, respectively $Ra_{mcr}=4{\rm \pi} ^2$ in the mushy phase, and $Ra_{lcr}=1708$ in the liquid. On the other hand, the transition from LPC to LC is due to the disappearance of temperature anomaly in the equation of state of the salty water solution. This occurs when $T_{max}< T_{0e}$, or equivalently, roughly when $S_e \gtrsim 2.7$%. The experimental measurements can be mapped into this posterior three-dimensional phase space, identified by the triplet $(Ra_m, Ra_l, S_e)$; see figure 5(b). It is worth noticing that while $S_i$ varies in $[0,3.5]\,\%$, $S_e$ can vary in a range from 0 % to approximately 7 %. We observe that five out of six possible equilibrium state freezing modes (i.e. the $2 \times 3$ combinations of MD, MC states with LD, LPC, LC ones) can be realized in our experiments and identified by our model. Figures 5(c–g) show photos of the interface profile in these typical cases. We may observe that non-flat interfaces are always associated with convective motions in either the liquid or the mushy phase, or both.