
% ------------------------------------------------------------------------
% Function: CRSS( x_1,x_2 )
% ------------------------------------------------------------------------
% Aim: 
% To solve the zero-temperature energy barrier and corresponding CRSS
% ------------------------------------------------------------------------
% Ref:
% 1/ Leyson, G. P. M., L. G. Hector, and W. A. Curtin. Acta Materialia 60.13 (2012): 5197-5203.
% 2/ Leyson, G. P. M., L. G. Hector, and W. A. Curtin. Acta Materialia 60.9 (2012): 3873-3884.
% ------------------------------------------------------------------------

function [ crss0,Eb0,wc ] = CRSS( x_1,x_2 )

    % define global variables

        global eb es;
        global solute_concentration LT;
        global shear_modulus poisson_ratio;
        global mistype;
        global BurVect;
        global wc_range;

    % define the volume misfit parameter (eb) and slip misfit parameter (es)

        eb=x_1;
        es=x_2;

    % dislocation line tension (unit: 1E-10*J/m )
    % Ref: 
    % 1/ U.F. Kocks, A.S. Argon, M.F. Ashby, Thermodynamic and Kinetics of 
    %    Slip, vol. 19, Pergamon Press, Oxford, 1975.

        if mistype=='e'
            LT=(1-2*poisson_ratio)/(1-poisson_ratio)*shear_modulus*BurVect^2/2/10;   % for edge dislocation
        end

        if mistype=='s'
            LT=(1+poisson_ratio)/(1-poisson_ratio)*shear_modulus*BurVect^2/2/10;   % for screw dislocation
        end

    % solve the characteristic roughtnening amplitude wc by minimizing Etot 

        options=optimset('MaxIter',500,'TolX',1E-8);
        [wc,~]=fminbnd(@Etotal,wc_range(1),wc_range(2),options);
        
        % [wc,~]=fminbnd(@Etotal,0,20,options);
        % for i=10:10:30
        %     [wc_tmp,~]=fminbnd(@Etotal,i,i+10,options);
        %     if Etotal(wc_tmp) < Etotal(wc)
        %         wc=wc_tmp;
        %     end
        % end

    % zero-temperature energy barrier( ev ) and corresponding CRSS ( MPa )

        Eb0=0.0625*((4*sqrt(2)-1)/3)*(3^(5/6)/2^(5/3))*(solute_concentration*wc^2*LT*Epoten(wc)^2/BurVect)^(1/3);

        crss0=10^4*pi*(4*sqrt(2)-1)/(8*6^(1/3))*((solute_concentration^2*Epoten(wc)^4)/(LT*BurVect^5*wc^5))^(1/3);

end

