clc
clear
close all
%% MATERIAL PARAMETERS %%
A = 3*10^-10;

% ADJUST Q AND n FOR EACH LOADING CONDITION
n = 1.89;
Q = 123;
R = 8.314;

%% EXPERIMENTAL DATA %%
sigma_cr = 47*10^6;

temp = [523 523 523 523 548 548 548 548 558 558 558 558 573 573 573 573];
stress = [14.4 21.6 28.8 33.6 14.1 21.15 28.2 32.9 13.95 20.925 27.9 32.55 13.5 20.25 27 31.5]*10^6;

eps_ss = [4.1406 6.0938 9.01961 14.6299 6.0548 7.29631 16.8293 28.521 5.22214 6.57895 23.3066 59.6639 5.46939 27.0492 135.593 1102.22]/10^8;

temp1 = [523 523 523 523];
stress1 = [14.4 21.6 28.8 33.6]*10^6;

eps_ss1 = [4.1406 6.0938 9.01961 14.6299]/10^8;

temp2 = [548 548 548 548];
stress2 = [14.1 21.15 28.2 32.9]*10^6;

eps_ss2 = [6.0548 7.29631 16.8293 28.521]/10^8;

temp3 = [558 558 558 558];
stress3 = [13.95 20.925 27.9 32.55]*10^6;

eps_ss3 = [5.22214 6.57895 25.8963 72.2689]/10^8;

temp4 = [573 573 573 573];
stress4 = [13.5 20.25 27 31.5]*10^6;

eps_ss4 = [5.46939 27.0492 135.593 1102.22]/10^8;

Z1 = log10(eps_ss1./exp(-10^3./(R*temp1)));
Z2 = log10(eps_ss2./exp(-10^3./(R*temp2)));
Z3 = log10(eps_ss3./exp(-10^3./(R*temp3)));
Z4 = log10(eps_ss4./exp(-10^3./(R*temp4)));

log_stress1 = log10(stress1);
log_stress2 = log10(stress2);
log_stress3 = log10(stress3);
log_stress4 = log10(stress4);


i = 0;
for sigma = 10*10^6:10^5:35*10^6;
    i = i+1;
    
    Temp = 573;

if Temp>548 && Temp<=558
    if sigma<=0.45*sigma_cr
        Q = 123*10^3;
    end
    if sigma>0.45*sigma_cr && sigma<=0.6*sigma_cr
        a = 123-(sigma/sigma_cr-0.45)/0.15*3;
        Q = (123-(Temp-548)/10*(123-a))*10^3;
    end
    if sigma>0.6*sigma_cr
        a = 120-(sigma/sigma_cr-0.6)/0.1*4;
        Q = (123-(Temp-548)/10*(123-a))*10^3;
    end
end

if Temp>558
    if sigma<=0.3*sigma_cr
        Q = 123*10^3;
    end
    if sigma>0.3*sigma_cr && sigma<=0.45*sigma_cr
        a = 123-(sigma/sigma_cr-0.3)/0.15*4;
        Q = (123-(Temp-558)/15*(123-a))*10^3;
    end
    if sigma>0.45*sigma_cr && sigma<=0.6*sigma_cr
        a = 119-(sigma/sigma_cr-0.45)/0.15*5;
        b = 123-(sigma/sigma_cr-0.45)/0.15*3;
        Q = (b-(Temp-558)/15*(b-a))*10^3;
    end
    if sigma>0.6*sigma_cr
        a = 114-(sigma/sigma_cr-0.6)/0.1*8;
        b = 120-(sigma/sigma_cr-0.6)/0.1*4;
        Q = (b-(Temp-558)/15*(b-a))*10^3;
    end
end

eps_model1(i) = A.*sigma.^n.*exp(-Q./(R.*temp4(1)));
sigma1(i) = sigma;
% eps_model2 = A.*sigma.^n.*exp(-Q*10^3./(R.*temp2));
% eps_model3 = A.*sigma.^n.*exp(-Q*10^3./(R.*temp3));
% eps_model4 = A.*sigma.^n.*exp(-Q*10^3./(R.*temp4));
end

%% PLOT FIGURE %%
figure
plot(stress4/10^6,eps_ss4,'*r','LineWidth',2)
hold on
plot(sigma1/10^6,eps_model1,'LineWidth',2)
xlabel('stress (MPa)')
ylabel('$\dot{\epsilon}_{ss}$ (1/s)','Interpreter','latex','FontSize',12)
legend('Experiment','Model')
text(20,2.8*10^-5,'Tempeature = 573K')

