-
Notifications
You must be signed in to change notification settings - Fork 0
/
run_kclee_li.m
77 lines (65 loc) · 1.79 KB
/
run_kclee_li.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
clear all; %close all;
load('exp_data.mat')
t_0 = 1e-12;
chi_0 = 0.25;
% elastic properties
nu = 0.381;
G_298 = 2.83e3;
%% STEADY STATE
% calculation at two strain rates, T = const
T = 298;
q_1 = t_0*4.e-5;
q_2 = t_0*3e-4;
sig_1 = 0.48;
sig_2 = 0.62;
rho_ss = exp(-1/chi_0);
A = log(log(sqrt(rho_ss)/q_2)) - log(log(sqrt(rho_ss)/q_1));
sig_T = (sig_1 - sig_2)/A;
% T_p = T*exp(sig_1/sig_T)*log(sqrt(rho_ss)/q_1);
% calculation at multiple T
q = t_0*3e-5;
T_p=1.05*398*(-log(q)-1/(2*chi_0));% if 1.05-->1, 398K-->0
T_p=18500;
T_all = [198 248 273 298 348 398];
sig_ss = [1.85 0.99 0.71 0.6 0.49 0.37];
sig_T = zeros(size(sig_ss));
mu_T = zeros(size(sig_ss));
for ii = 1:length(T_all)
sig_T(ii) = sig_ss(ii)/(log(T_p/T_all(ii))-log(log(sqrt(rho_ss)/q)));
mu_T(ii) = sig_T(ii)/sqrt(rho_ss);
end
% ratio between Taylor stress and shear modulus
r_param = mu_T(4)/G_298;
clear A ii q q_1 q_2 sig_1 sig_2 sig_ss T
%% calculation for T=298 and edot=3e-5
eps_dot = 3.e-5;
T = 298;
color=['b' 'r' 'm' 'g' 'c' 'k'];
hold on
espan = [0 0.01];
e = linspace(0,0.01,1000);
c0 =500;
c1 = 1;
Kx = c0*exp(T_all/c1);
K_pf = 5;
% K_p = 1e5;
rho_ini = 7e-3;
chi_ini = 0.15;
y0 = [0.0 rho_ini chi_ini];
j=1;
legend_list={'198K','248K','273K','298K','348K','398K'};
for ii=1:6
e_exp = Li(~isnan(Li(:,2*ii-1)),2*ii-1);
s_exp = Li(~isnan(Li(:,2*ii)),2*ii);
% mu_T(3) corresponds to T=298 K
param = struct('mu_T',mu_T(ii),'Kx',Kx(ii),'K_pf',K_pf,'t_0',t_0,'T_p',T_p,'chi_0',chi_0,'r',r_param,'nu',nu);
sol = ode15s(@(t,y)kclee(t,y,T_all(ii),eps_dot,param),espan,y0);
s = deval(sol,e,j);
% error(4)=mean((deval(sol,e4,1)-sig4).^2);
scatter(e_exp,s_exp,color(ii));
plot(e*100,s,color(ii));
% xlabel('e(%)')
% ylabel('Stress, MPa')
end
% legend
title('Different Temperature')