17 subroutine iso_creep(matl, sectType, stress, strain, extval,plstrain, &
18 dtime,ttime,stiffness, temp)
20 integer,
intent(in) :: sectType
21 real(kind=
kreal),
intent(in) :: stress(6)
22 real(kind=
kreal),
intent(in) :: strain(6)
23 real(kind=
kreal),
intent(in) :: extval(:)
24 real(kind=
kreal),
intent(in) :: plstrain
25 real(kind=
kreal),
intent(in) :: ttime
26 real(kind=
kreal),
intent(in) :: dtime
27 real(kind=
kreal),
intent(out) :: stiffness(6,6)
28 real(kind=
kreal),
optional :: temp
32 real(kind=
kreal) :: ina(1), outa(3)
33 real(kind=
kreal) :: xxn, aa
35 real(kind=
kreal) :: c3,e,un,g,stri(6),p,dstri,c4,c5,f,df, eqvs
43 if( dtime==0.d0 .or. all(stress==0.d0) )
return
47 if(
present(temp) )
then
49 call fetch_tabledata(
mc_isoelastic, matl%dict, outa(1:2), ierr, ina )
51 call fetch_tabledata(
mc_isoelastic, matl%dict, outa(1:2), ierr )
54 stop
"error in isotropic elasticity definition"
61 if( matl%mtype==
norton )
then
62 if(
present( temp ) )
then
64 call fetch_tabledata(
mc_norton, matl%dict, outa, ierr, ina )
66 call fetch_tabledata(
mc_norton, matl%dict, outa, ierr )
69 aa=outa(1)*((ttime+dtime)**(outa(3)+1.d0)-ttime**(outa(3)+1.d0))/(outa(3)+1.d0)
77 p=-(stri(1)+stri(2)+stri(3))/3.d0
82 dstri=dsqrt(1.5d0*(stri(1)*stri(1)+stri(2)*stri(2)+stri(3)*stri(3)+ &
83 2.d0*(stri(4)*stri(4)+stri(5)*stri(5)+stri(6)*stri(6))) )
90 if( eqvs<1.d-10 ) eqvs=1.d-10
97 c4=c3*plstrain/(dstri+3.d0*g*plstrain)
98 c3=c4-c3*df/(3.d0*g*df+1.d0)
103 stiffness(i,j) = stiffness(i,j) +c3*stri(i)*stri(j)
107 stiffness(i,i) = stiffness(i,i) - c4
109 stiffness(i,j) = stiffness(i,j) +c5
113 stiffness(i,i) = stiffness(i,i) - c4/2.d0
123 integer,
intent(in) :: secttype
124 real(kind=
kreal),
intent(in) :: strain(6)
125 real(kind=
kreal),
intent(inout) :: stress(6)
126 real(kind=
kreal),
intent(inout) :: extval(:)
127 real(kind=
kreal),
intent(out) :: plstrain
128 real(kind=
kreal),
intent(in) :: ttime
129 real(kind=
kreal),
intent(in) :: dtime
130 real(kind=
kreal),
optional :: temp
134 real(kind=
kreal) :: ina(1), outa(3)
135 real(kind=
kreal) :: xxn, aa
137 real(kind=
kreal) :: e,un,g,dg,ddg,stri(6),p,dstri,f,df, eqvs
142 if( dtime==0.d0 )
return
146 if(
present(temp) )
then
148 call fetch_tabledata(
mc_isoelastic, matl%dict, outa(1:2), ierr, ina )
150 call fetch_tabledata(
mc_isoelastic, matl%dict, outa(1:2), ierr )
153 stop
"error in isotropic elasticity definition"
160 if( matl%mtype==
norton )
then
161 if(
present( temp ) )
then
163 call fetch_tabledata(
mc_norton, matl%dict, outa, ierr, ina )
165 call fetch_tabledata(
mc_norton, matl%dict, outa, ierr )
168 stop
"error in isotropic elasticity definition"
171 aa=outa(1)*((ttime+dtime)**(outa(3)+1.d0)-ttime**(outa(3)+1.d0))/(outa(3)+1.d0)
181 p=-(stri(1)+stri(2)+stri(3))/3.d0
186 dstri=dsqrt(1.5d0*(stri(1)*stri(1)+stri(2)*stri(2)+stri(3)*stri(3)+ &
187 2.d0*(stri(4)*stri(4)+stri(5)*stri(5)+stri(6)*stri(6))) )
193 if( matl%mtype==
norton )
then
194 eqvs = dstri-3.d0*g*dg
197 ddg = (f-dg)/(3.d0*g*df+1.d0)
199 if((ddg<dg*1.d-6).or.(ddg<1.d-12))
exit
203 stri(:) = stri(:)-3.d0*g*dg*stri(:)/dstri
204 stress(1:3) = stri(1:3)-p
205 stress(4:6) = stri(4:6)
220 gauss%fstatus(2) = gauss%fstatus(2)+gauss%plstrain
integer(kind=4), parameter kreal
This module provides functions for elastic material.
subroutine calelasticmatrix(matl, secttype, d, temp)
Calculate isotropic elastic matrix.
This module provides functions for creep calculation.
subroutine update_iso_creep(matl, secttype, strain, stress, extval, plstrain, dtime, ttime, temp)
This subroutine calculates stresses and creep status for an elastically isotropic material with isotr...
subroutine iso_creep(matl, secttype, stress, strain, extval, plstrain, dtime, ttime, stiffness, temp)
This subroutine calculates stiffness for elastically isotropic materials with isotropic creep.
subroutine updateviscostate(gauss)
Update viscoplastic state.
This module summarizes all infomation of material properties.
character(len=dict_key_length) mc_norton
integer(kind=kint), parameter norton
character(len=dict_key_length) mc_isoelastic
This modules defines a structure to record history dependent parameter in static analysis.
Stucture to management all material relates data.
All data should be recorded in every quadrature points.