% Try to reproduce figure 20 in Lacmann_1999
% try to reproduce Zacher_1996

% clearly given are:
% sigma = 0.05
% T = 25°C
% G_s = 6.7 * 10^9 N/m^2
% nu_p = 0.33
% l_m = 0.2 mu-m = 0.2 * 10^-6 m
%
% the length of the Burger vector is not 100% clear, as it is written b in
% the equations and text, but B in the figure. B_0 only exists as
% nucleation rate in the list of symbols. Best to assume now is:
% b = 1.22 * 10^-9 m (MISPRINT IN LACMANN)
%
%
% Also unclear is whether the angle is to be set in the equation as degree
% or as radian. So either:
% eta = 1
% or
% eta = 1*pi/180

% Equation 11 states:
%  R*T*sigma = 2*V_m * gamma_c / l_m
% If we assume that the interface tension should equal the grain boundary
% tension from equation 9, we have:
%  sigma = 2*V_m/l_m * G_s * b * eta *sqrt(l_m) / 4 / pi / (1-nu_p) /
%  sqrt(L) * log(0.43 / eta * sqrt(L / l_m)
% Let's try this one:

sigma = 0.05;              % - (supersaturation)
temperature = 25 + 273.15; % K

gasConstant = 8314*10^-3 % J / mol /K
%boltzmannConstant = 1.380648*10^-23; % J/K
shearModule = 6.7*10^9;         % N/m^2
poissonRatio = 0.33;            % -
mosaicLength = 0.2*10^-6;       % m
burgerLength = 1.22*10^-9;      % m
mosaicAngle = 0.5*pi/180;       % degree
% not provided:
%molarVolume = 474.39/1725*10^-3         % m^3/mol
%1/molarVolume
molarConstant = 3.7*10^3; % mol/m^3
% moleculeVolume = molarVolume/6.022141/10^23 % m^3/molecule

L = (1:100) * 10^-6;
% gamma_gb = shearModule  *burgerLength * mosaicAngle * sqrt(mosaicLength) / ...
%     (4*pi*(1-poissonRatio)) ./ sqrt(L) .* log(0.43/mosaicAngle *sqrt(L/mosaicLength));
% sigma = 2*molarVolume / mosaicLength / boltzmannConstant / temperature .* gamma_gb;

c1 = 0.24 * shearModule * burgerLength / (1-poissonRatio) / molarConstant / sqrt(mosaicLength)
c2 = 0.43 / sqrt(mosaicLength)
sigma = sigma - c1 ./ gasConstant ./ temperature .* mosaicAngle ./ ...
    sqrt(L) .* log(c2 .* sqrt(L) ./ mosaicAngle);

plot(L, sigma)
xlim([min(L) max(L)])
ylim([-0.05 0.05])
set(gca, 'xscale', 'log')