Crystal Nucleation Theory (CNT)

//One universal basic required here to get things going once loaded
window.onload = function () {
    //restoreDefaultValues(); //Un-comment this if you want to start with defaults
    Main();
};
//Any global variables go here


//Main is hard wired as THE place to start calculating when input changes
//It does no calculations itself, it merely sets them up, sends off variables, gets results and, if necessary, plots them.
function Main() {
    saveSettings();

    //Send all the inputs as a structured object
    //If you need to convert to, say, SI units, do it here!
    const inputs = {
        A: sliders.SlideA.value,
        gamma: sliders.Slidegamma.value/1e3, //mJ/m2 to J/m2
        v: sliders.Slidev.value*1e-30, //A3 to m3
        S: sliders.SlideS.value,
        Vol: sliders.SlideVol.value/1e6, //ml to m3
        TK: sliders.SlideT.value+273, //C to K
        tMax: sliders.SlidetMax.value,
        lS2Max: sliders.SlidelS2Max.value,
    }

    //Send inputs off to CalcIt where the names are instantly available
    //Get all the resonses as an object, result
    const result = CalcIt(inputs)

    //Set all the text box outputs
    document.getElementById('B').value = result.B
    document.getElementById('nstar').value = result.nstar
    //Do all relevant plots by calling plotIt - if there's no plot, nothing happens
    //plotIt is part of the app infrastructure in app.new.js
    if(result.plots) {
        for(let i = 0; i < result.plots.length; i++) {
            plotIt(result.plots[i], result.canvas[i]);
        }
    }

    //You might have some other stuff to do here, but for most apps that's it for CalcIt!
}

//Here's the app calculation
//The inputs are just the names provided - their order in the curly brackets is unimportant!
//By convention the input values are provided with the correct units within Main
function CalcIt({A,gamma,v,S,Vol,TK,tMax,lS2Max}) {
    const kT=TK*1.38e-23
    const B=16*Math.PI*Math.pow(gamma,3)*v*v/(3*Math.pow(kT,3))
    const nstar=2*B/Math.pow(Math.log(S),2)
    let J=A*Math.exp(-B/Math.log(S*S))
    let PCurve=[],JCurve=[]
    for (let t=0;t<=tMax;t++){
        P=1-Math.exp(-J*Vol*t)
        PCurve.push({x:t,y:P})
    }
    for (let lS2=1;lS2<=lS2Max;lS2++){
        J=A*Math.exp(-B*lS2)
        JCurve.push({x:lS2,y:Math.log(J)})
    }


    //Now set up all the graphing data detail by detail.
    const prmap = {
        plotData: [PCurve], //An array of 1 or more datasets
        lineLabels: ["Probability"], //An array of labels for each dataset
        hideLegend: true,
       xLabel: "t&s", //Label for the x axis, with an & to separate the units
        yLabel: "P& ", //Label for the y axis, with an & to separate the units
        y2Label: null, //Label for the y2 axis, null if not needed
        yAxisL1R2: [], //Array to say which axis each dataset goes on. Blank=Left=1
        logX: false, //Is the x-axis in log form?
        xTicks: undefined, //We can define a tick function if we're being fancy
        logY: false, //Is the y-axis in log form?
        yTicks: undefined, //We can define a tick function if we're being fancy
        legendPosition: 'top', //Where we want the legend - top, bottom, left, right
        xMinMax: [,], //Set min and max, e.g. [-10,100], leave one or both blank for auto
        yMinMax: [,], //Set min and max, e.g. [-10,100], leave one or both blank for auto
        y2MinMax: [,], //Set min and max, e.g. [-10,100], leave one or both blank for auto
        xSigFigs: 'F0', //These are the sig figs for the Tooltip readout. A wide choice!
        ySigFigs: 'F3', //F for Fixed, P for Precision, E for exponential
    };
    const prmap1 = {
        plotData: [JCurve], //An array of 1 or more datasets
        lineLabels: ["J"], //An array of labels for each dataset
        hideLegend: true,
        xLabel: "1/ln²S& ", //Label for the x axis, with an & to separate the units
        yLabel: "lnJ& ", //Label for the y axis, with an & to separate the units
        y2Label: null, //Label for the y2 axis, null if not needed
        yAxisL1R2: [], //Array to say which axis each dataset goes on. Blank=Left=1
        logX: false, //Is the x-axis in log form?
        xTicks: undefined, //We can define a tick function if we're being fancy
        logY: false, //Is the y-axis in log form?
        yTicks: undefined, //We can define a tick function if we're being fancy
        legendPosition: 'top', //Where we want the legend - top, bottom, left, right
        xMinMax: [,], //Set min and max, e.g. [-10,100], leave one or both blank for auto
        yMinMax: [,], //Set min and max, e.g. [-10,100], leave one or both blank for auto
        y2MinMax: [,], //Set min and max, e.g. [-10,100], leave one or both blank for auto
        xSigFigs: 'F1', //These are the sig figs for the Tooltip readout. A wide choice!
        ySigFigs: 'F2', //F for Fixed, P for Precision, E for exponential
    };

   //Now we return everything - text boxes, plot and the name of the canvas, which is 'canvas' for a single plot

    return {
        B:B.toPrecision(3),
        nstar:nstar.toPrecision(3),
        plots: [prmap,prmap1],
        canvas: ['canvas','canvas1'],
    };

}
                        

Crystallization is frustrating to study because it can take anywhere from 1 to 10,000s for a first crystal to appear; it's a matter of probability. We have had CNT, Crystal Nucleation Theory, for decades that relates the probabilities to a few key parameters. For generations, researchers have measured the probability of a crystal appearing after a given time then calculated those parameters. Here we input the parameters and see what the experimental results would look like. This "reverse" app technique is not as strange as it sounds. You still get to see the theory working - it's just that you get an idealised version of what the experiments would have found. For those who have never heard of Ideal Solubility and how it impacts the possible saturated solubility (and its temperature behaviour), the Ideal Solubility app provides an instant guide.

The key experimental data are gathered using modern high-throughput crystallization techniques. Tubes of volume V are heated to a temperature that dissolves all the solute then cooled to a temperature where the solution is satured by a value S, typically 1.1 to 1.2, i.e 10 or 20% saturated. An optical device detects the first crystal, the time is recorded, the tube re-heated, re-cooled, re-measured till enough data points are gathered for any give S. The whole thing is then repeated at a few different values of S. For each S we end up with a curve of the probability P that a crystal has appeared at time t. In the app you see the curve for your chosen S and, because the timescales for these curves vary dramatically you need to set t_Max which is the longest time you are prepared to wait to determine the P curve.

Note that this cycling technique can be used to measure the temperature dependent saturated solubility by noting the temperature when crystals appear and disappear (rather than waiting for appearance at a T that corresponds to an S). Any single value is likely to be misleading, but with enough cycles, a good estimate can be obtained.

The P curve is analysed via an equation that involves a steady state nucleation rate, J, the time t and the tube volume V via:

`P=1-exp(-JVt)`

So now we need to analyse J. CNT tells us that it depends on S (of course) and on a kinetic parameter A and a thermodynamic parameter B. In the app we just set A as a number because it's rather hard to relate it to experimental values (though higher viscosity tends to lead to lower A). The key equation is:

`J=Aexp(-B/(ln^2S))`

Now you see why we need curves at different values of S. By plotting lnJ versus 1/ln²S we can find A and B.

Finally, B can be related to the molecular volume v of the solute, the temperature T, the surface energy of the nucleus, γ and the Boltzmann constant k via:

`B=(16πγ^3v^2)/(3(kT)^3)`

Finally, n*, the number of molecules in a critical nucleus (one big enough to grow spontaneously into a large crystal), can be calculated via:

`n"*"=(2B)/(ln^2S)`

Has it been worth it?

After all those hours waiting for crystals to appear then a bit of data fitting, you know A, B, γ and n*. What do you do with those values? The answer seems to be "not much other than writing papers that speculate on the meaning of those values for a specific solute in a few solvents". A lot of people (including myself) find many of the claimed γ values to be meaningless, and who knows if the n* values are correct, especially when they often turn out to be close to 1.

That sounds a bit pessimistic, but it's largely true. There are countless papers on correlations of those parameters with various properties of the solvents, but there is not a single reliable rule that anyone can apply to a new solute in its own solvents.

A core reason for the failure is that CNT itself is not correct. More recently the "two-step" formalism has gained in popularity because it explains a number of puzzles that CNT cannot explain. Standard CNT simply has proto-crystals identical to the real ones, fluctuating in size till the thermodynamics say that they are large enough to grow spontaneously. This is simply not true. Molecules cluster more or less in a different arrangements that may include specific interactions (e.g. hydrogen bonding donor/acceptor) with the solvent. It turns out to be perfectly possible for these clusters to grow in size till there is a liquid-liquid phase separation, what we all know as "oiling out". Oiling out isn't some annoying mistake on the part of the solute - it's a respectable thermodynamic outcome when molecular clusters fluctuate naturally.

I also have a problem with the fact that although J is based on the time of first appearance and linked to the theory of nucleus size, there seems general agreement that nucleation takes place on bits of junk (or maybe stray seeds from past attempts) in the solution. A counter to this is that by the time was see the "first" crystals, we're well past the primal seeding event and indeed, we're seeing larger crystals growing on secondary nuclei from the breakup of primary ones. It seems to me to be a handwaving mess - or maybe I just don't get it.

The future

The idea that CNT is flawed and that measuring lots more P curves is not a great idea has been stated by experts in the field in authoritative reviews so is not just my personal view. With some smart use of Raman, FTIR, neutron & x-ray scattering and with some ingenious experiments where the re-heating is done to just above full dissolution, it's becoming possible to get a better idea of the pre-clustering effects (or their absence in those systems where crystallization seems to be easy). The specific effects of the solvents on those pre-clusters are also starting to be identified.

My personal view is that we need to switch over to a Kirkwood-Buff style of thinking about solubility and crystallization. After all, KB theory is the "Fluctuation Theory of Solubility" and fluctuations are core to crystal formation and/or oiling out. Watch this space for an update on whether my view is justified or not.