-
Notifications
You must be signed in to change notification settings - Fork 5
/
kHardening.f
143 lines (102 loc) · 4.88 KB
/
kHardening.f
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
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
********************************************
** KHARDENING updates the state variables **
** that determine mechanical hardening **
********************************************
SUBROUTINE kHardening(pdot,p,plasStrainrate,dtime,slip,gammaDot,
+ nSys,irradiate,tauSolute,gammast,rhossd,iphase,rhofor,
+ Temperature,rhosub,twinon,twinvolfrac,nTwin,
+ nTwinStart,nTwinEnd,gammaTwinDot)
INCLUDE 'ABA_PARAM.INC'
! number of slip systems
INTEGER,intent(in) :: nSys
! plastic strain rate
REAL*8,intent(in) :: plasStrainRate(3,3)
! time increment
REAL*8,intent(in) :: dtime
! plastic shear rate on slip systems
! and absolute value
REAL*8,intent(in) :: gammadot(nSys)
! activate irradiation effect
INTEGER,intent(in) :: irradiate
! prefactor for SSD evolution equation
REAL*8,intent(in) :: gammast
! crystal type
INTEGER,intent(in) :: iphase
! Current temperature
REAL*8,intent(in) :: Temperature
! twin systems activation flag
INTEGER,intent(in) :: twinon
! total number of twin systems
INTEGER,intent(in) :: nTwin
! the active twins are the ones in the
! interval [nTwinStart,nTwinEnd] in the
! twin system file
INTEGER,intent(in) :: nTwinStart,nTwinEnd
! twin rate
REAL*8,intent(in) :: gammaTwinDot(nTwin)
! Von Mises invariant plastic strain rate
REAL*8,intent(inout) :: p
! crystallographic slip
! needed for model with irradiation
REAL*8,intent(inout) :: slip
! increase in tauc due to solute force
REAL*8,intent(inout) :: tauSolute
! total SSD density
REAL*8,intent(inout) :: rhossd
! forest dislocation density
REAL*8,intent(inout) :: rhofor(nSys)
! substructure dislocation density
REAL*8,intent(inout) :: rhosub
! twin volume fraction
REAL*8,intent(inout) :: twinvolfrac(nTwin)
! Von Mises invariant plastic strain rate
REAL*8,intent(out) :: pdot
! absolute value of the plastic shear rate on slip systems
REAL*8,dimension(nSys) :: tgammadot
INTEGER :: I
! calculate Von Mises invariant plastic strain rate
pdot=sqrt(2./3.*sum(plasStrainrate*plasStrainrate))
p = p + pdot*dtime
! update crystallographic slip
slip = slip + sum(abs(gammaDot))*dtime
! update solute force
if (irradiate == 1) then
tauSolute = 750.0*exp(-slip/0.025)
end if
!=========================================================================
! SSD Evolution
rhossd = rhossd + (gammast*pdot*dtime)
!=========================================================================
!=========================================================================
! Dislocation density evolution (explicit Euler time integration)
! alpha-Uranium model
if (iphase == 5) then
! slip independent of twin fraction
DO I=1,nSys
tgammaDot(I) = abs(gammaDot(I))
END DO
! forest dislocations evolution
! using constants from calibration of tensile bar 3
! using twin-slip interaction model
rhofor(1) = rhofor(1) + 43.2*max(sqrt(rhofor(1))-(0.17100+2.6093e-03*Temperature)*rhofor(1),0.0)*tgammaDot(1)*dtime
rhofor(2) = rhofor(2) + 6320.0*max(sqrt(rhofor(2))-(0.25650+5.8708e-04*Temperature)*rhofor(2),0.0)*tgammaDot(2)*dtime
rhofor(3) = rhofor(3) + 0.24*max(sqrt(rhofor(3))-(0.11718+1.9289e-04*Temperature)*rhofor(3),0.0)*tgammaDot(3)*dtime
rhofor(4) = rhofor(4) + 0.24*max(sqrt(rhofor(4))-(0.11718+1.9289e-04*Temperature)*rhofor(4),0.0)*tgammaDot(4)*dtime
rhofor(5) = rhofor(5) + 800.0*max(sqrt(rhofor(5))-(0.12+1.5e-05*Temperature)*rhofor(5),0.0)*tgammaDot(5)*dtime
rhofor(6) = rhofor(6) + 800.0*max(sqrt(rhofor(6))-(0.12+1.5e-05*Temperature)*rhofor(6),0.0)*tgammaDot(6)*dtime
rhofor(7) = rhofor(7) + 800.0*max(sqrt(rhofor(7))-(0.12+1.5e-05*Temperature)*rhofor(7),0.0)*tgammaDot(7)*dtime
rhofor(8) = rhofor(8) + 800.0*max(sqrt(rhofor(8))-(0.12+1.5e-05*Temperature)*rhofor(8),0.0)*tgammaDot(8)*dtime
! substructure dislocations evolution
rhosub = rhosub + 0.216*(17.545+0.26771*Temperature)*rhofor(1)*sqrt(rhosub)*tgammaDot(1)*dtime
! twin volume fraction evolution
! dot(f)^beta = dot(gamma)^beta / gamma^twin
! see Kalidindi JMPS 1998 (equations 32 and 33a)
! if twinvolfractot = 1.0 => all gammaTwinDot are zero => no evolution
! McCabe 2010: (130) twin induces a total strain 0.299
if (twinon == 1) then ! twin active
twinvolfrac(nTwinStart) = twinvolfrac(nTwinStart) + gammaTwinDot(nTwinStart)*dtime/0.299
twinvolfrac(nTwinEnd) = twinvolfrac(nTwinEnd) + gammaTwinDot(nTwinEnd)*dtime/0.299
end if
end if
RETURN
END