7 use hecmw,
only : kint, kreal
17 subroutine stf_c1( etype,nn,ecoord,area,gausses,stiff, u ,temperature )
19 integer(kind=kint),
intent(in) :: etype
20 integer(kind=kint),
intent(in) :: nn
21 real(kind=kreal),
intent(in) :: ecoord(3,nn)
22 real(kind=kreal),
intent(in) :: area
24 real(kind=kreal),
intent(out) :: stiff(:,:)
25 real(kind=kreal),
intent(in),
optional :: u(:,:)
26 real(kind=kreal),
intent(in),
optional :: temperature(nn)
28 real(kind=kreal) llen, llen0, elem(3,nn)
30 real(kind=kreal) ina(1), outa(1), direc(3), direc0(3), coeff, strain
31 integer(kind=kint) :: i,j
35 elem(:,:) = ecoord(:,:) + u(:,:)
37 elem(:,:) = ecoord(:,:)
40 direc = elem(:,2)-elem(:,1)
41 llen = dsqrt( dot_product(direc, direc) )
43 direc0 = ecoord(:,2)-ecoord(:,1)
44 llen0 = dsqrt( dot_product(direc0, direc0) )
46 if(
present(temperature) )
then
47 ina(1) = 0.5d0*(temperature(1)+temperature(2))
48 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa, ierr, ina )
50 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa, ierr )
52 if( ierr ) outa(1) = gausses(1)%pMaterial%variables(m_youngs)
53 coeff = outa(1)*area*llen0/(llen*llen)
54 strain = gausses(1)%strain(1)
58 stiff(i,i) = coeff*strain
60 stiff(i,j) = stiff(i,j) + coeff*(1.d0-2.d0*strain)*direc(i)*direc(j)
64 stiff(4:6,1:3) = -stiff(1:3,1:3)
65 stiff(1:3,4:6) = transpose(stiff(4:6,1:3))
66 stiff(4:6,4:6) = stiff(1:3,1:3)
73 subroutine update_c1( etype, nn, ecoord, area, u, du, qf ,gausses, TT, T0 )
78 integer(kind=kint),
intent(in) :: etype
79 integer(kind=kint),
intent(in) :: nn
80 real(kind=kreal),
intent(in) :: ecoord(3,nn)
81 real(kind=kreal),
intent(in) :: area
82 real(kind=kreal),
intent(in) :: u(3,nn)
83 real(kind=kreal),
intent(in) :: du(3,nn)
84 real(kind=kreal),
intent(out) :: qf(nn*3)
86 real(kind=kreal),
intent(in),
optional :: tt(nn)
87 real(kind=kreal),
intent(in),
optional :: t0(nn)
90 real(kind=kreal) :: direc(3), direc0(3)
91 real(kind=kreal) :: llen, llen0, ina(1), outa(1)
92 real(kind=kreal) :: elem(3,nn)
93 real(kind=kreal) :: young
94 real(kind=kreal) :: ttc, tt0, alp, alp0, epsth
99 elem(:,:) = ecoord(:,:) + u(:,:) + du(:,:)
101 direc = elem(:,2)-elem(:,1)
102 llen = dsqrt( dot_product(direc, direc) )
104 direc0 = ecoord(:,2)-ecoord(:,1)
105 llen0 = dsqrt( dot_product(direc0, direc0) )
108 if(
present(tt) .and.
present(t0) )
then
109 ttc = 0.5d0*(tt(1)+tt(2))
110 tt0 = 0.5d0*(t0(1)+t0(2))
113 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa, ierr, ina )
114 if( ierr ) outa(1) = gausses(1)%pMaterial%variables(m_youngs)
117 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa(:), ierr, ina )
118 if( ierr ) outa(1) = gausses(1)%pMaterial%variables(m_exapnsion)
122 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa(:), ierr, ina )
123 if( ierr ) outa(1) = gausses(1)%pMaterial%variables(m_exapnsion)
128 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa, ierr )
129 if( ierr ) outa(1) = gausses(1)%pMaterial%variables(m_youngs)
133 gausses(1)%strain(1) = dlog(llen/llen0)
134 gausses(1)%stress(1) = young*(gausses(1)%strain(1)-epsth)
137 gausses(1)%strain_out(1) = gausses(1)%strain(1)
138 gausses(1)%stress_out(1) = gausses(1)%stress(1)
140 qf(1) = gausses(1)%stress(1)*area*llen0/llen
141 qf(1:3) = -qf(1)*direc
153 integer(kind=kint),
intent(in) :: etype,nn
155 real(kind=kreal),
intent(out) :: ndstrain(nn,6)
156 real(kind=kreal),
intent(out) :: ndstress(nn,6)
158 ndstrain(1,1:6) = gausses(1)%strain(1:6)
159 ndstress(1,1:6) = gausses(1)%stress(1:6)
160 ndstrain(2,1:6) = gausses(1)%strain(1:6)
161 ndstress(2,1:6) = gausses(1)%stress(1:6)
174 integer(kind=kint),
intent(in) :: ETYPE
176 real(kind=kreal),
intent(out) :: strain(6)
177 real(kind=kreal),
intent(out) :: stress(6)
180 strain(:) = gausses(1)%strain(1:6)
181 stress(:) = gausses(1)%stress(1:6)
187 subroutine dl_c1(etype, nn, xx, yy, zz, rho, thick, ltype, params, &
193 integer(kind = kint),
intent(in) :: etype, nn
194 real(kind = kreal),
intent(in) :: xx(:), yy(:), zz(:)
195 real(kind = kreal),
intent(in) :: params(0:6)
196 real(kind = kreal),
intent(inout) :: vect(:)
197 real(kind = kreal) :: rho, thick
198 integer(kind = kint) :: ltype, nsize
200 integer(kind = kint) :: ndof = 3
201 integer(kind = kint) :: ivol, isuf, i
202 real(kind = kreal) :: vx, vy, vz, val, a, aa
209 if( ltype .LT. 10 )
then
211 else if( ltype .GE. 10 )
then
218 vect(1:nsize) = 0.0d0
221 if( ivol .EQ. 1 )
then
222 if( ltype .EQ. 4 )
then
223 aa = dsqrt( ( xx(2)-xx(1) )*( xx(2)-xx(1) ) &
224 +( yy(2)-yy(1) )*( yy(2)-yy(1) ) &
225 +( zz(2)-zz(1) )*( zz(2)-zz(1) ) )
231 vx = vx/dsqrt( params(1)**2+params(2)**2+params(3)**2 )
232 vy = vy/dsqrt( params(1)**2+params(2)**2+params(3)**2 )
233 vz = vz/dsqrt( params(1)**2+params(2)**2+params(3)**2 )
236 vect(3*i-2) = val*rho*a*0.5d0*aa*vx
237 vect(3*i-1) = val*rho*a*0.5d0*aa*vy
238 vect(3*i ) = val*rho*a*0.5d0*aa*vz
255 type (hecmwST_matrix) :: hecMAT
256 type (hecmwST_local_mesh) :: hecMESH
257 integer(kind=kint) :: itype, is, iE, ic_type, icel, jS, j, n
259 do itype = 1, hecmesh%n_elem_type
260 ic_type = hecmesh%elem_type_item(itype)
261 if(ic_type == 301)
then
262 is = hecmesh%elem_type_index(itype-1) + 1
263 ie = hecmesh%elem_type_index(itype )
265 js = hecmesh%elem_node_index(icel-1)
267 n = hecmesh%elem_node_item(js+j)
268 if( hecmat%D(9*n-8) == 0.0d0)
then
269 hecmat%D(9*n-8) = 1.0d0
272 if( hecmat%D(9*n-4) == 0.0d0)
then
273 hecmat%D(9*n-4) = 1.0d0
276 if( hecmat%D(9*n ) == 0.0d0)
then
277 hecmat%D(9*n ) = 1.0d0
294 type (hecmwST_matrix) :: hecMAT
295 type (hecmwST_local_mesh) :: hecMESH
296 integer :: n, nn, is, iE, i, a
300 is = hecmat%IndexL(n-1)+1
301 ie = hecmat%IndexL(n )
303 if(hecmat%AL(9*i-8) /= 0.0d0) a = 1
304 if(hecmat%AL(9*i-7) /= 0.0d0) a = 1
305 if(hecmat%AL(9*i-6) /= 0.0d0) a = 1
307 is = hecmat%IndexU(n-1)+1
308 ie = hecmat%IndexU(n )
310 if(hecmat%AU(9*i-8) /= 0.0d0) a = 1
311 if(hecmat%AU(9*i-7) /= 0.0d0) a = 1
312 if(hecmat%AU(9*i-6) /= 0.0d0) a = 1
314 if(hecmat%D(9*n-7) /= 0.0d0) a = 1
315 if(hecmat%D(9*n-6) /= 0.0d0) a = 1
317 hecmat%D(9*n-8) = 1.0d0
323 is = hecmat%IndexL(n-1)+1
324 ie = hecmat%IndexL(n )
326 if(hecmat%AL(9*i-5) /= 0.0d0) a = 1
327 if(hecmat%AL(9*i-4) /= 0.0d0) a = 1
328 if(hecmat%AL(9*i-3) /= 0.0d0) a = 1
330 is = hecmat%IndexU(n-1)+1
331 ie = hecmat%IndexU(n )
333 if(hecmat%AU(9*i-5) /= 0.0d0) a = 1
334 if(hecmat%AU(9*i-4) /= 0.0d0) a = 1
335 if(hecmat%AU(9*i-3) /= 0.0d0) a = 1
337 if(hecmat%D(9*n-5) /= 0.0d0) a = 1
338 if(hecmat%D(9*n-3) /= 0.0d0) a = 1
340 hecmat%D(9*n-4) = 1.0d0
346 is = hecmat%IndexL(n-1)+1
347 ie = hecmat%IndexL(n )
349 if(hecmat%AL(9*i-2) /= 0.0d0) a = 1
350 if(hecmat%AL(9*i-1) /= 0.0d0) a = 1
351 if(hecmat%AL(9*i ) /= 0.0d0) a = 1
353 is = hecmat%IndexU(n-1)+1
354 ie = hecmat%IndexU(n )
356 if(hecmat%AU(9*i-2) /= 0.0d0) a = 1
357 if(hecmat%AU(9*i-1) /= 0.0d0) a = 1
358 if(hecmat%AU(9*i ) /= 0.0d0) a = 1
360 if(hecmat%D(9*n-2) /= 0.0d0) a = 1
361 if(hecmat%D(9*n-1) /= 0.0d0) a = 1
363 hecmat%D(9*n ) = 1.0d0
This module encapsulate the basic functions of all elements provide by this software.
This module defined coomon data and basic structures for analysis.
real(kind=kreal), pointer ref_temp
REFTEMP.
This module provide common functions of 3D truss elements.
subroutine elementstress_c1(ETYPE, gausses, strain, stress)
subroutine nodalstress_c1(ETYPE, NN, gausses, ndstrain, ndstress)
subroutine search_diag_modify(n, nn, hecMAT, hecMESH)
subroutine update_c1(etype, nn, ecoord, area, u, du, qf, gausses, TT, T0)
Update strain and stress inside element.
subroutine stf_c1(etype, nn, ecoord, area, gausses, stiff, u, temperature)
This subroutine calculate stiff matrix of 2-nodes truss element.
subroutine truss_diag_modify(hecMAT, hecMESH)
subroutine dl_c1(etype, nn, xx, yy, zz, rho, thick, ltype, params, vect, nsize)
This modules defines a structure to record history dependent parameter in static analysis.
All data should be recorded in every quadrature points.