This dataset includes data and processing code from calculations of diffusive mushy layer growth and linear convective stability analysis of this system, and is made available under the GNU General Public License v3:
http://www.gnu.org/licenses/gpl.html
Note that thiscode has only been tested on original author's system, and we accept no responsibility for adverse consequences if used on other systems.
Use at your own risk.

In particular, be wary that MATLAB scripts supplied may include code to modify the system path/write data/clear data from a MATLAB workspace.  Use with caution.

Simulation results produced using code written by Joe Hitchen and edited by Andrew Wells for use with MATLAB. This code base was used to generate figures contained in the PhD thesis:
Hitchen (2017), “Understanding the growth and convective instability of mushy layers, with application to young sea ice”, PhD Thesis University of Oxford.

and two accompanying papers (in press at the Journal of Fluid Mechanics at the time of initial submission of this repository):

Joseph R. Hitchen & Andrew J. Wells, in press, J. Fluid Mech. doi:10.1017/jfm.2024.808 
Structure of mushy layers grown from perfectly and imperfectly conducting boundaries. Part 1. Diffusive solidification 
Joseph R. Hitchen & Andrew J. Wells, in press, J. Fluid Mech. doi:10.1017/jfm.2024.809
Structure of mushy layers  grown from perfectly and imperfectly conducting boundaries. Part 2. Onset of convection

Contact andrew.wells@physics.ox.ac.uk for further details.
Contents of this data repository follow below. 
Readme 

This zipped folder contains MATLAB code used to calculate the numerical results shown in the papers Hitchen \& Wells (2024a,b), Journal of Fluid Mechanics, hereafter Part 1 and Part 2. The code was originally written for use with MATLAB R2017a, but selected scripts have been edited for compatibility with MATLAB R2021b/R2022a following deprecation of certain plotting commands. 

The primary scripts to produce figures are as follows, which can be run from the MATLAB GUI command line by changing to the directory where this unzipped folder arrives, and typing the script name. All estimated run times are from a system with 3.5 GHz Dual-Core Intel Core i7 processor. 

Calculations are initially run using the notation from the thesis Hitchen (2017), “Understanding the growth and convective instability of mushy layers, with application to young sea ice”, PhD Thesis University of Oxford.

They are then converted into the revised notation in the papers Hitchen \& Wells (2024a,b) using the following transformations:
\[ \mathcal{C} = \frac{\Lambda}{\theta_{L\infty}}, \qquad \mathcal{S}_{new} = \frac{\mathcal{S}_{old} }{\theta_{L\infty} }, \qquad \theta_{\infty} = \frac{1}{\theta_{L\infty} }, \quad \RaylML_{,new} = \theta_{L\infty} \RaylML_{,old}  \]
where subscript “new” refers to the definitions in the papers Hitchen \& Wells (2024a,b) and subscript “old” refers to the thesis. 

Part 1:

CodeStore/DiffusiveSolidification/Part1Fig2a3abc.m % produces figures figures 2a, 3a, 3b, 3c of Part 1 by running relevant parameter sweeps that each calculate a solution of the self-similar growth problem, with varying concentration ratio and Stefan number. Estimated run time approximately 2 minutes. 

CodeStore/DiffusiveSolidification/Part1Fig2b.m % produces figures figure 2b of Part 1 by running relevant parameter sweeps that each calculate a solution of the self-similar growth problem, with varying concentration ratio and far field liquid temperature. Estimated run time approximately 10 minutes. 

CodeStore/DiffusiveSolidification/Part1Fig4.m % Produces Figure 4 of Hitchen and Wells (2024a) part 1 (plotting Liquid fraction versus scaled depth, and concentration ratio to illustrate collapse of Stefan limit for three different far-field temperatures). Loads data from ‘Fig3.10DataConcRat.mat
‘, which was computed by ConductingMush_InitialLqdus.m. Run time < 1 minute for Part1Fig4.m 

CodeStore/DiffusiveSolidification/ConductingMush_InitialLqdus.m % produces data used in figure 4 and 5 of Hitchen & Wells (2024a) by numerical integration of mushy layer equations with a perfectly conducting boundary.  Takes ~5 minutes to run, so data is saved in ‘Fig3.10DataConcRat.mat’ and ‘Part2NewNotation/Data/Fig2binvspace.mat’

CodeStore/DiffusiveSolidification/Part1Fig5.m % Produces Figure 5 of Hitchen and Wells (2024a) part 1 (plotting scaling collapses of growth rate versus \Omega in near eutectic limit, and versus Stefan number raised to -0.5 power in Stefan limit). Loads data from ‘Data/Fig2aData.mat' computed by ConductingMush_LqdusGrad.m and ‘Part2NewNotation/Data/Fig2binvspace.mat’, which was computed by ConductingMush_InitialLqdus.m. Run time < few minutes for Part1Fig5.m 

CodeStore/DiffusiveSolidification/ConductingMush_LqdusGrad.m % produces data used in figure 5 of Hitchen & Wells (2024a) by numerical integration of mushy layer equations with a perfectly conducting boundary and varying \Lambda from the thesis.  Takes nontrivial time to run, so data is saved in '/Part2NewNotation/Data/Fig2aData.mat' 

CodeStore/DiffusiveSolidification/Part1Fig6.m % produces figure 6 of Part 1 for an imperfectly conducting boundary, including the dashed line for the full numerical solution. For the numerical simulation, runs with different Biot numbers are spliced to find evolution with the time-dependent scaled Biot number. Estimated run time approximately 1 minutes.  Set to run with lower resolution by default than used in the paper to speed up initial run (but alternative values are provided on line 45 for a better resolved simulation taking a longer time to run). 

CodeStore/DiffusiveSolidification/Part1Fig678.m % produces data for figures 6, 7 and 8 of Part 1 for an imperfectly conducting boundary, by running relevant parameter sweeps that each calculate a solution of the growth problem with a conducting boundary, with varying concentration ratio and far field liquid temperature. Runs with different Biot numbers are spliced to find evolution with the time-dependent scaled Biot number. Estimated run time approximately 2 minutes. 

CodeStore/DiffusiveSolidification/Part1Fig9GraphicalAbstract.m % creates either Figure 9 (liquid fraction versus scaled Biot number and scaled depth) or graphical abstract (setting PSFRAG=true on line 82). Loads data from ‘CodeStore/Fig9Data/Fig9Data.mat’ that was pre-computed (slowly!) using CodeStore/ImperfectConduction_Mushy_Fig9.m. Estimated run time of Part1Fig9GraphicalAbstract.m is <1 minute. 

CodeStore/DiffusiveSolidification/ImperfectConduction_Mushy_Fig9.m % produces liquid fraction data versus depth and scaled Biot number needed for Figure 9. Note: at default ultra high resolution (Line 40) this takes many (4<N < 16) hours to run. Hence data has been saved in the file ‘CodeStore/Fig9Data/Fig9Data.mat’. 

CodeStore/DiffusiveSolidification/Part1Fig10.m % produces data for figure 10 of Part 1 for an example of sea ice growth with an imperfectly conducting boundary, and comparing to observational data of ice thickness from Notz & Worster (2008) cited in Hitchen & Wells (2024a). Estimated run time <1 minute. 

Note that the above scripts call a set of subsidiary code/functions stored in scripts in the folders:

CodeStore/CoreFunctions/

And hence are written to run within the supplied relative path tree from this directory. 

Part 2

CodeStore/MushyConvection/LinearMushyConvection/Part2Figs2678.m % Produces Figures 2, 6, 7 and 8 of Hitchen and Wells (2024b) part 2 
loads data from 'JFMData/Figs2_7Part2StabilityData.mat' that was originally computed by Thesis_LqdusGrad_Contours_JFM.m  which runs a sweep over a range of parameters, calculating diffusive mushy layer growth with a perfectly conducting isothermal boundary and then running the linear convective stability profile for each parameter value. (long computation time of ~30 mins, so data saved and loaded for plotting). Function called is available:

CodeStore/MushyConvection/Thesis_LqdusGrad_Contours_JFM.m

CodeStore/MushyConvection/LinearMushyConvection/Part2Fig3.m % Produces Figure 3 of Hitchen and Wells (2024b) part 2 
that investigates how the scaled Rayleigh number varies with concentration ratio for various parameter combinations. 
loads data from 'JFMData/' 
files: 'CritRayl-LqdusInf0.2.mat', 'CritRayl-LqdusInf0.4.mat', 'CritRayl-LqdusInf0.6.mat', 'CritRayl-LqdusInf0.8.mat'
These were originally computed by varying the parameter ‘Phys.LqdusInf’ on line 28 of
Thesis_LqdusGrad_Contours_JFM.m (each simulation ~30 minutes per value of Phys.LqdusInf

CodeStore/MushyConvection/LinearMushyConvection/Part2Fig4.m % Plots panels in figure 4 of Hitchen and Wells (2024b) part 2,
Examines how the background liquid fraction, background permeability, temperature perturbation and streamfunction perturbation vary with concentration  ratio for mushy layer growth from a perfectly conducting isothermal boundary. Run time <30s 

CodeStore/MushyConvection/LinearMushyConvection/Part2Fig5.m % Plots panels in figure 5 of Hitchen and Wells (2024b) part 2,
Examines how the background liquid fraction, background permeability, liquid fraction perturbation and streamfunction perturbation vary with concentration  ratio for mushy layer growth from a perfectly conducting isothermal boundary. Run time <30s. This code essentially uses the same code structure as Part2Fig4.m, but with temperature perturbation replaced by liquid fraction perturbation

CodeStore/MushyConvection/LinearMushyConvection/Part2Fig9.m  % Plots panels in figure 9 of Hitchen and Wells (2024b) part 2. 
Examines how the critical Rayleigh number, critical wavenumber and background temperature difference across the mush vary with scaled Biot number (or over time)   for mushy layer growth from an imperfectly conducting  boundary. Uses the reduced model for high concentration ratio. Run time <30s.

CodeStore/MushyConvection/LinearMushyConvection/Part2Fig10.m % Plots panels in figure 10 of Hitchen and Wells (2024b) part 2. Examines how the critical Rayleigh number, critical wavenumber and background liquid fraction profiles versus depth each vary with scaled Biot number (or over time)  for mushy layer growth from an imperfectly conducting  boundary. Uses the full mushy layer growth model for low concentration ratio.
Loads data 'JFMData/Fig10Part2-BiotStabilityData.mat' which has been precomputed using the script Thesis_Biot_Mushy_JFM.m at Ultra hi res (see line 33 of Thesis_Biot_Mushy_JFM.m to switch resolution) - this took many hours of run time so ouptut was saved to speed up plotting. Run time for this plotting script <30s.

CodeStore/MushyConvection/LinearMushyConvection/Thesis_Biot_Mushy_JFM.m - In thesis notation, calculates mushy layer growth as a function of scaled Biot number/time, then investigates the impact on the critical Rayleigh number and wavenumber. Generates data that was saved in 'JFMData/Fig10Part2-BiotStabilityData.mat’. Warning - extremely slow run time (between 4-24 hours) if run with high/ultra high resolution - see line 33-64  for options

CodeStore/MushyConvection/LinearMushyConvection/Part2Fig12.m % Plots panels in figure 12 of Hitchen and Wells (2024b) part 2 to examine stability in a field application (and a subsidiary liquid fraction versus depth/time plot).  Run time ~ 25 minutes at mid resolution ( see lines 68-81 for options).

Note that the above scripts call a set of subsidiary code/functions stored in scripts in the folders:

CodeStore/MushyConvection/StabilityFunctions/
CodeStore/MushyConvection/SolidificationFunctions/
CodeStore/MiscCode

And hence are written to run within the supplied relative path tree from this directory.