29 subroutine matlmatrix( gauss, sectType, matrix, time, dtime, cdsys, temperature, isEp )
31 integer,
intent(in) :: sectType
32 real(kind=kreal),
intent(out) :: matrix(:,:)
33 real(kind=kreal),
intent(in) :: time
34 real(kind=kreal),
intent(in) :: dtime
35 real(kind=kreal),
intent(in) :: cdsys(3,3)
36 real(kind=kreal),
intent(in),
optional :: temperature
37 integer(kind=kint),
intent(in),
optional :: isEp
41 real(kind=kreal) :: cijkl(3,3,3,3)
46 if(
present(isep) )
then
47 if( isep == 1 )flag = 1
53 if(
present(temperature) )
then
58 elseif(
iselastic(matl%mtype) .or. flag==1 )
then
66 if(
present(temperature) )
then
72 if(
present(temperature) )
then
78 print *,
"Elasticity type", matl%mtype,
"not supported"
84 call mat_c2d( cijkl, matrix, secttype )
87 call mat_c2d( cijkl, matrix, secttype )
90 call mat_c2d( cijkl, matrix, secttype )
94 if(
present( temperature ) )
then
96 gauss%istatus(1), gauss%fstatus, gauss%plstrain, matrix, temperature )
99 gauss%istatus(1), gauss%fstatus, gauss%plstrain, matrix )
102 call umatlmatrix( matl%name, matl%variables(101:), gauss%strain, &
103 gauss%stress, gauss%fstatus, matrix, dtime, time )
104 elseif( matl%mtype==
norton )
then
105 if(
present( temperature ) )
then
106 call iso_creep( matl, secttype, gauss%stress, gauss%strain, gauss%fstatus, &
107 gauss%plstrain, dtime, time, matrix, temperature )
109 call iso_creep( matl, secttype, gauss%stress, gauss%strain, gauss%fstatus, &
110 gauss%plstrain, dtime, time, matrix )
113 stop
"Material type not supported!"
120 subroutine stressupdate( gauss, sectType, strain, stress, cdsys, time, dtime, temp, tempn )
122 integer,
intent(in) :: sectType
123 real(kind=kreal),
intent(in) :: strain(6)
124 real(kind=kreal),
intent(out) :: stress(6)
125 real(kind=kreal),
intent(in) :: cdsys(3,3)
126 real(kind=kreal),
intent(in),
optional :: time
127 real(kind=kreal),
intent(in),
optional :: dtime
128 real(kind=kreal),
optional :: temp
129 real(kind=kreal),
optional :: tempn
138 call uelasticupdate( gauss%pMaterial%variables(101:), strain, stress )
140 if( .not.
present(dtime) ) stop
"error in viscoelastic update!"
141 if(
present(temp) .and.
present(tempn) )
then
142 call updateviscoelastic( gauss%pMaterial, secttype, strain, stress, gauss%fstatus, dtime, temp, tempn )
144 call updateviscoelastic( gauss%pMaterial, secttype, strain, stress, gauss%fstatus, dtime )
146 elseif ( gauss%pMaterial%mtype==
norton )
then
147 if( .not.
present(dtime) ) stop
"error in viscoelastic update!"
148 if(
present(temp) )
then
149 call update_iso_creep( gauss%pMaterial, secttype, strain, stress, gauss%fstatus, gauss%plstrain, dtime, time, temp )
151 call update_iso_creep( gauss%pMaterial, secttype, strain, stress, gauss%fstatus, gauss%plstrain, dtime, time )
154 call uupdate( gauss%pMaterial%name, gauss%pMaterial%variables(101:), &
155 strain, stress, gauss%fstatus, dtime, time )
162 real(kind=kreal),
intent(in) :: cijkl(3,3,3,3)
163 real(kind=kreal),
intent(out) :: dij(6,6)
164 integer,
intent(in) :: itype
169 dij(1,1) = cijkl(1,1,1,1)
170 dij(1,2) = cijkl(1,1,2,2)
171 dij(1,3) = cijkl(1,1,3,3)
172 dij(1,4) = cijkl(1,1,1,2)
173 dij(1,5) = cijkl(1,1,2,3)
174 dij(1,6) = cijkl(1,1,3,1)
175 dij(2,1) = cijkl(2,2,1,1)
176 dij(2,2) = cijkl(2,2,2,2)
177 dij(2,3) = cijkl(2,2,3,3)
178 dij(2,4) = cijkl(2,2,1,2)
179 dij(2,5) = cijkl(2,2,2,3)
180 dij(2,6) = cijkl(2,2,3,1)
181 dij(3,1) = cijkl(3,3,1,1)
182 dij(3,2) = cijkl(3,3,2,2)
183 dij(3,3) = cijkl(3,3,3,3)
184 dij(3,4) = cijkl(3,3,1,2)
185 dij(3,5) = cijkl(3,3,2,3)
186 dij(3,6) = cijkl(3,3,3,1)
187 dij(4,1) = cijkl(1,2,1,1)
188 dij(4,2) = cijkl(1,2,2,2)
189 dij(4,3) = cijkl(1,2,3,3)
190 dij(4,4) = cijkl(1,2,1,2)
191 dij(4,5) = cijkl(1,2,2,3)
192 dij(4,6) = cijkl(1,2,3,1)
193 dij(5,1) = cijkl(2,3,1,1)
194 dij(5,2) = cijkl(2,3,2,2)
195 dij(5,3) = cijkl(2,3,3,3)
196 dij(5,4) = cijkl(2,3,1,2)
197 dij(5,5) = cijkl(2,3,2,3)
198 dij(5,6) = cijkl(2,3,3,1)
199 dij(6,1) = cijkl(3,1,1,1)
200 dij(6,2) = cijkl(3,1,2,2)
201 dij(6,3) = cijkl(3,1,3,3)
202 dij(6,4) = cijkl(3,1,1,2)
203 dij(6,5) = cijkl(3,1,2,3)
204 dij(6,6) = cijkl(3,1,3,1)
207 dij(1,1) = cijkl(1,1,1,1)
208 dij(1,2) = cijkl(1,1,2,2)
209 dij(1,3) = cijkl(1,1,1,2)
210 dij(2,1) = cijkl(2,2,1,1)
211 dij(2,2) = cijkl(2,2,2,2)
212 dij(2,3) = cijkl(2,2,1,2)
213 dij(3,1) = cijkl(1,2,1,1)
214 dij(3,2) = cijkl(1,2,2,2)
215 dij(3,3) = cijkl(1,2,1,2)
217 dij(1,1) = cijkl(1,1,1,1)
218 dij(1,2) = cijkl(1,1,2,2)
219 dij(1,3) = cijkl(1,1,1,2)
220 dij(1,4) = cijkl(1,1,3,3)
221 dij(2,1) = cijkl(2,2,1,1)
222 dij(2,2) = cijkl(2,2,2,2)
223 dij(2,3) = cijkl(2,2,1,2)
224 dij(2,4) = cijkl(2,2,3,3)
225 dij(3,1) = cijkl(1,2,1,1)
226 dij(3,2) = cijkl(1,2,2,2)
227 dij(3,3) = cijkl(1,2,1,2)
228 dij(3,4) = cijkl(1,2,3,3)
229 dij(4,1) = cijkl(3,3,1,1)
230 dij(4,2) = cijkl(3,3,2,2)
231 dij(4,3) = cijkl(3,3,1,2)
232 dij(4,4) = cijkl(3,3,3,3)
241 (gauss, secttype, d, &
242 e1_hat, e2_hat, e3_hat, cg1, cg2, cg3, &
247 integer,
intent(in) :: sectType, n_layer
248 real(kind = kreal),
intent(out) :: d(:, :)
249 real(kind = kreal),
intent(in) :: e1_hat(3), e2_hat(3), e3_hat(3)
250 real(kind = kreal),
intent(in) :: cg1(3), cg2(3), cg3(3)
251 real(kind = kreal),
intent(out) :: alpha
255 real(kind = kreal) :: c(3, 3, 3, 3)
260 matl => gauss%pMaterial
266 (matl, secttype, c, &
267 e1_hat, e2_hat, e3_hat, cg1, cg2, cg3, &
272 stop
"Material type not supported!"
290 real(kind = kreal),
intent(in) :: c(:, :, :, :)
291 real(kind = kreal),
intent(out) :: d(:, :)
292 integer,
intent(in) :: itype
296 integer :: index_i(5), index_j(5), &
297 index_k(5), index_l(5)
298 integer :: i, j, k, l
345 d(is, js) = c(i, j, k, l)
This module provides functions for elastic material.
subroutine calelasticmatrix_ortho(matl, secttype, bij, dmat, temp)
Calculate orthotropic elastic matrix.
subroutine linearelastic_shell(matl, secttype, c, e1_hat, e2_hat, e3_hat, cg1, cg2, cg3, alpha, n_layer)
subroutine calelasticmatrix(matl, secttype, d, temp)
Calculate isotropic elastic matrix.
This module provide functions for elastoplastic calculation.
subroutine calelastoplasticmatrix(matl, secttype, stress, istat, extval, plstrain, d, temperature)
This subroutine calculates elastoplastic constitutive relation.
This module manages calculation relates with materials.
subroutine stressupdate(gauss, secttype, strain, stress, cdsys, time, dtime, temp, tempn)
Update strain and stress for elastic and hyperelastic materials.
subroutine matlmatrix_shell(gauss, secttype, d, e1_hat, e2_hat, e3_hat, cg1, cg2, cg3, alpha, n_layer)
subroutine matlmatrix(gauss, secttype, matrix, time, dtime, cdsys, temperature, isep)
Calculate constituive matrix.
integer function getnlgeomflag(gauss)
Fetch the nlgeom flag of the material.
subroutine mat_c2d_shell(c, d, itype)
subroutine mat_c2d(cijkl, dij, itype)
Transfer rank 4 constituive matrix to rank 2 form.
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.
This module provides functions for hyperelastic calculation.
subroutine calupdateelasticarrudaboyce(matl, secttype, dstrain, dstress)
This subroutine provides to update stress and strain for Arrude-Royce material.
subroutine calelasticmooneyrivlin(matl, secttype, cijkl, strain)
This subroutine provides elastic tangent coefficient for Mooney-Rivlin hyperelastic material.
subroutine calelasticarrudaboyce(matl, secttype, cijkl, strain)
This subroutine provides elastic tangent coefficient for Arruda-Boyce hyperelastic material.
subroutine calelasticmooneyrivlinaniso(matl, secttype, cijkl, strain, cdsys)
This subroutine provides elastic tangent coefficient for anisotropic Mooney-Rivlin hyperelastic mater...
subroutine calupdateelasticmooneyrivlinaniso(matl, secttype, strain, stress, cdsys)
This subroutine provides to update stress and strain for anisotropic Mooney-Rivlin material.
subroutine calupdateelasticmooneyrivlin(matl, secttype, strain, stress)
This subroutine provides to update stress and strain for Mooney-Rivlin material.
This module summarizes all infomation of material properties.
integer(kind=kint), parameter mooneyrivlin
integer(kind=kint), parameter mooneyrivlin_aniso
integer(kind=kint), parameter planestress
integer function getelastictype(mtype)
Get elastic type.
integer(kind=kint), parameter d3
integer(kind=kint), parameter planestrain
integer(kind=kint), parameter elastic
integer(kind=kint), parameter arrudaboyce
integer(kind=kint), parameter shell
integer(kind=kint), parameter norton
integer(kind=kint), parameter axissymetric
integer(kind=kint), parameter userelastic
integer(kind=kint), parameter neohooke
integer(kind=kint), parameter usermaterial
logical function iselastic(mtype)
If it is an elastic material?
integer(kind=kint), parameter userhyperelastic
logical function isviscoelastic(mtype)
If it is an viscoelastic material?
logical function iselastoplastic(mtype)
If it is an elastoplastic material?
This modules defines a structure to record history dependent parameter in static analysis.
This module provides interface for elastic or hyperelastic calculation.
subroutine uelasticmatrix(matl, strain, d)
This subroutine calculates constitutive relation.
subroutine uelasticupdate(matl, strain, stress)
This subroutine calculate updated strain and stress.
This subroutine read in used-defined material properties tangent.
subroutine uupdate(mname, matl, strain, stress, fstat, dtime, ttime, temperature)
This subroutine calculate strain and stress increment.
subroutine umatlmatrix(mname, matl, strain, stress, fstat, d, dtime, ttime, temperature)
This subroutine calculates constitutive matrix.
This module provides functions for viscoelastic calculation.
subroutine calviscoelasticmatrix(matl, secttype, dt, d, temp)
This subroutine calculates tangent moduli for isotropic viscoelastic material.
subroutine updateviscoelastic(matl, secttype, eps, sig, vsig, dt, temp, tempn)
This subroutine provides to update stress for viscoelastic material.
Stucture to management all material relates data.
All data should be recorded in every quadrature points.