%%
%% ===================================
%%  OxygenExperiments.m
%%  by Bart Trzynadlowski, 2012-2013
%% ===================================
%%
%% Using OModel(), replicates all experimental results discussed in my doctoral 
%% dissertation ("Reduced Moment-Based Models for Oxygen Precipitates and 
%% Dislocation Loops in Silicon," University of Washington, 2013).
%%
%% Version History
%% ---------------
%%
%% 2013.03.20:
%%  - Initial release included in doctoral dissertation.
%%

%
% OxygenExperiments()
%
% Runs the RKPM oxygen model with the best-fit parameters for a series of
% oxygen precipitation experiments and plots the results of each. Experiments 
% are run one by one. The full set takes a long time to complete.
%
function OxygenExperiments
    global NetCI1200 IOCToASTM ppmaToConc;
    
    % CI*-CV* at T=1200 C is used to express the initial NetCI in a more
    % convenient form
    NetCI1200 = ((2980/51.4)*5e22*exp((-4.95+1.77)/(8.62e-5*(1200+273.15)))-(86/3.07)*5e22*exp((-4.56+2.12)/(8.62e-5*(1200+273.15))));
    
    % Conversion factor for IOC-88 -> new ASTM (1983) calibration standard
    IOCToASTM = 2.45 / 3.14;
    
    % Conversion factor for ppma -> cm^-3
    ppmaToConc  = 5e22 / 1e6;
    
    % Run all the experiments
    tic;
    ExperimentStewart();
    ExperimentChiou();
    ExperimentSwaroop();
    ExperimentKennel();
    time = toc;
    
    fprintf('All simulations finished. Run time: %1.1f minutes.\n', time/60);
end

%
% ExperimentStewart()
%
% SANS study of precipitation during a single-step 750 C anneal.
%
% Reference:
%
%   R. J. Stewart, S. Messoloras and S. Rycroft, "A small angle neutron 
%   scattering study of oxygen precipitation in silicon," in Proceedings of the 
%   NATO Advanced Research Workshop, Exeter, UK, 1996.
%
function ExperimentStewart
    global NetCI1200 IOCToASTM;
    
    % Best-fit initial conditions
    InitialM0       = (1e8)^0;
    InitialSize     = 72 + (1e4)^0;
    InitialNetCI    = NetCI1200 * (1e-3)^0.54;
    
    % 
    % Data from Stewart et al.
    %
    % The initial oxygen concentration was measured using FTIR with the IOC-88
    % calibration standard. Precipitated oxygen was measured using SANS, which
    % should not need to be re-calibrated. Coincidentally, however, I discovered
    % that scaling the SANS measurements does preduce a much better fit, but 
    % that is not done here.
    %
    
    COInit      = IOCToASTM * 9.6e17;               % FTIR measurement
    Times       = [ 24 72 216 ];                    % [hr]
    COData      = COInit - 1e17 * [ 1.5 1.3 6.2 ];  % SANS measurements
    
    % Simulate
    [n, t, temp, CO, m0, m1] = OModel('rkpm', 1.1, COInit, InitialM0, InitialSize, InitialNetCI, 400, 'steps', [ 750, 1000 ]);
    
    % Plot    
    figure;
    semilogx(t/3600, CO, '-', Times, COData, 's');
    axis([10 1000 0 1e18]);
    legend('Simulation', 'Data', 'Location', 'NorthEast');        
    title('SANS Experiment by Stewart et al. (1996)');
    xlabel('Time (hr)');
    ylabel('Interstitial Oxygen (cm^{-3})');
end

%
% ExperimentSwaroop()
%
% Two-step precipitation experiment: 750 C / 4 hr + 1050 C / 16 hr, 10 C/min 
% ramp rate. The data here is a statistical reduction of a larger set of 
% measurements that appeared in the original paper by Swaroop et al.
%
% Reference:
%
%   R. Swaroop, N. Kim, W. Lin, M. Bullis, L. Shive, A. Rice, E. Castel and M. 
%   Christ, "Testing for Oxygen Precipitation in Silicon Wafers," Solid State 
%   Technology, pp. 85-89, March 1987.
%
function ExperimentSwaroop
    global NetCI1200 ppmaToConc;
    
    % Best-fit initial conditions
    InitialM0       = (1e8)^1;
    InitialSize     = 72 + (1e4)^0.9939;
    InitialNetCI    = NetCI1200 * (1e-3)^1.4;

    % Statistically reduced data from Swaroop et al. Initial vs. precipitated
    % oxygen. FTIR data, new ASTM (1983) calibration factor.
    COInit = ppmaToConc * [ 10.2261 10.47 11.3236 11.8461 12.3515 12.8209 13.3938 13.7954 14.2295 14.8199 15.219 15.8443 16.2794 16.7318 17.3759 17.7934 18.2113 18.6116 ];
    COLoss = ppmaToConc * [ 0.392663 0.464778 0.503442 0.719044 0.471481 1.6131 4.03753 2.93464 4.61048 6.9281 8.81756 10.9572 11.5287 12.2072 12.7438 13.4221 13.7799 14.1375 ];

    % Run each point
    for i = 1:length(COInit)
        [n, t, temp, CO] = OModel('rkpm', 1.1, COInit(i), InitialM0, InitialSize, InitialNetCI, 400, 'steps', [ 750, 4, 10, 1050, 16 ]);
        SimLoss(i) = COInit(i)-CO(end);
    end
    
    % Plot
    figure;
    axis([5e17 10e17 0 10e17]);
    plot(COInit, SimLoss, '-', COInit, COLoss, 's');
    legend('Simulation (750 ^oC / 4 hr + 1050 ^oC /16 hr)', 'Data', 'Location', 'NorthWest');
    title('Two-Step Precipitation Experiment by Swaroop et al. (1987)');
    xlabel('Initial Oxygen (cm^{-3})');
    ylabel('Precipitated Oxygen (cm^{-3})');
end

%
% ExperimentChiou()
%
% Two-step precipitation experiments by Chiou & Shive: 
%
%   1. Set A: 800 C / 2hr + 1050 C / 16 hr, 40 C/min ramp rate
%   2. Set B: 800 C / 1hr + 1050 C / 8 hr, 40 C/min ramp rate
%
% Reference:
%
%   H.-D. Chiou and L. W. Shive, "Test Methods for Oxygen Precipitation in
%   Silicon," VLSI Science and Technology/1985, pp. 429-435, 1985.
%
function ExperimentChiou
    global NetCI1200 ppmaToConc;
    
    % Best-fit initial conditions
    InitialM0       = (1e8)^1.1234;
    InitialSize     = 72 + (1e4)^0.972;
    InitialNetCI    = NetCI1200 * (1e-3)^1.036;
    
    %
    % Statistically reduced FTIR data from Chiou & Shive (Fig. 5), new ASTM
    %
    
    % Set A: 800 C / 2 hr + 1050 C / 16 hr
    COInitA = ppmaToConc * [ 13.42363 14.10402 14.49531 15.0197 15.5011 16.0044 16.5077 17.0547 17.558 18.0613 18.5646 19.0678 ];
    COLossA = ppmaToConc * [ 0.586008 1.296908 1.9181 3.862636572 5.15996605 6.973240236 9.101962048 11.08297337 11.880519 12.30062409 13.03411303 13.33950589 ];
    
    % Set B: 800 C / 1 hr + 1050 C / 8 hr
    COInitB = ppmaToConc * [ 12.7002 12.8753 13.6411 13.9256 14.4945 14.9759 15.9825 16.5295 17.0328 17.4486 18.5646 19.046 19.5711 20.0525 20.4683 ];
    COLossB = ppmaToConc * [ 0.331459384 0.422737441 0.463782308 0.68762132 0.887194537 1.459914686 2.267171456 3.480235626 4.287054927 5.929382049 9.98575922 10.95526322 12.73558444 12.88391971 13.33950589 ];

    % Run set A    
    for i = 1:length(COInitA)
        [n, t, temp, CO] = OModel('rkpm', 1.1, COInitA(i), InitialM0, InitialSize, InitialNetCI, 400, 'steps', [ 800, 2, 40, 1050, 16 ]);
        SimLossA(i) = COInitA(i)-CO(end);
    end
    
    % Run set B
    for i = 1:length(COInitB)
        [n, t, temp, CO] = OModel('rkpm', 1.1, COInitB(i), InitialM0, InitialSize, InitialNetCI, 400, 'steps', [ 800, 1, 40, 1050, 8 ]);
        SimLossB(i) = COInitB(i)-CO(end);
    end
 
    % Plot
    figure;
    axis([5e17 10e17 0 10e17]);
    plot(COInitA, SimLossA, '-', COInitA, COLossA, 's', COInitB, SimLossB, '-', COInitB, COLossB, 's');
    legend('Simulation (800 ^oC / 2 hr + 1050 ^oC /16 hr)', 'Data', 'Simulation (800 ^oC / 1 hr + 1050 ^oC / 8 hr)', 'Data', 'Location', 'NorthWest');
    title('Two-Step Precipitation Experiments by Chiou & Shive (1985)');
    xlabel('Initial Oxygen (cm^{-3})');
    ylabel('Precipitated Oxygen (cm^{-3})');
end

%
% ExperimentKennel()
%
% Precipitation experiments with variable-length nucleation (750 C) and growth
% (1100 C) anneals.
%
% Reference:
%
%   H. W. Kennel, "Physical Modeling of Oxygen Precipitation, Defect Formation, 
%   and Diffusion in Silicon," doctoral dissertation, Stanford University, 1991.
%
function ExperimentKennel
    global NetCI1200 IOCToASTM;
    
    % Best-fit initial conditions
    InitialM0       = (1e8)^1;
    InitialSize     = 72 + (1e4)^0.9;
    InitialNetCI    = NetCI1200 * (1e-3)^1.1;
    
    %
    % Data from Kennel Fig. 3.5
    %
    
    COInit =        IOCToASTM * 9.8e17; % initial oxygen
    AnnealTime =    [ 0 4 8 20 ];       % growth step anneal times [hr]
    
    % CO remaining after anneal for 0 hr and 2 hr nucleation times
    COData0 =   IOCToASTM * [ 9.5522E+17 9.1429E+17 8.8331E+17 7.87055E+17  ];
    COData2 =   IOCToASTM * [ 9.6019E+17 8.74485E+17 7.3074E+17 4.93526E+17 ];
    
    % Run 0 hr nucleation data
    for i = 1:length(AnnealTime)        
        if AnnealTime(i) ~= 0
            [n, t, temp, CO] = OModel('rkpm', 1.1, COInit, InitialM0, InitialSize, InitialNetCI, 400, 'steps', [ 750, 1/3600, 20, 800, 10/3600, 7, 1100, AnnealTime(i) ]);
        else
            [n, t, temp, CO] = OModel('rkpm', 1.1, COInit, InitialM0, InitialSize, InitialNetCI, 400, 'steps', [ 750, 1/3600, 20, 800, 10/3600 ]);
        end
        COSim0(i) = CO(end);
    end
    
    % Run 2 hr nucleation data
    for i = 1:length(AnnealTime)        
        if AnnealTime(i) ~= 0
            [n, t, temp, CO] = OModel('rkpm', 1.1, COInit, InitialM0, InitialSize, InitialNetCI, 400, 'steps', [ 750, 2, 20, 800, 10/3600, 7, 1100, AnnealTime(i) ]);
        else
            [n, t, temp, CO] = OModel('rkpm', 1.1, COInit, InitialM0, InitialSize, InitialNetCI, 400, 'steps', [ 750, 2, 20, 800, 10/3600 ]);
        end
        COSim2(i) = CO(end);
    end
    
    % Plot
    figure;
    axis([0 10 2.5e17 1e18]);
    plot(AnnealTime, COSim0, '-', AnnealTime, COData0, 's', AnnealTime, COSim2, '-', AnnealTime, COData2, 's');
    legend('Simulation (No nucleation step)', 'Data', 'Simulation (2 hr 750 ^oC nucleation)', 'Data', 'Location', 'SouthWest');
    title('1100 ^oC Annealing Experiments by Kennel (1991)');
    xlabel('Anneal Time (hr)');
    ylabel('Interstitial Oxygen (cm^{-3})');
end