10 private :: update_abort
29 type (hecmwST_matrix) :: hecMAT
30 type (hecmwST_local_mesh) :: hecMESH
31 type (fstr_solid) :: fstrSOLID
32 real(kind=kreal),
intent(in) :: time
33 real(kind=kreal),
intent(in) :: tincr
34 integer,
intent(in) :: iter
36 integer(kind=kint) :: nodLOCAL(fstrSOLID%max_ncon)
37 real(kind=kreal) :: ecoord(3, fstrsolid%max_ncon)
38 real(kind=kreal) :: thick
39 integer(kind=kint) :: ndof, itype, is, iE, ic_type, nn, icel, iiS, i, j
41 real(kind=kreal) :: total_disp(6, fstrsolid%max_ncon), du(6, fstrsolid%max_ncon), ddu(6, fstrsolid%max_ncon)
42 real(kind=kreal) :: tt(fstrsolid%max_ncon), tt0(fstrsolid%max_ncon), ttn(fstrsolid%max_ncon)
43 real(kind=kreal) :: qf(fstrsolid%max_ncon*6), coords(3, 3)
44 integer :: isect, ihead, cdsys_ID
45 integer :: ndim, initt
47 real(kind=kreal),
optional :: strainenergy
48 real(kind=kreal) :: tmp
49 real(kind=kreal) :: ddaux(3,3)
52 fstrsolid%QFORCE=0.0d0
79 do itype = 1, hecmesh%n_elem_type
80 is = hecmesh%elem_type_index(itype-1)+1
81 ie = hecmesh%elem_type_index(itype )
82 ic_type= hecmesh%elem_type_item(itype)
83 if (hecmw_is_etype_link(ic_type)) cycle
84 if (hecmw_is_etype_patch(ic_type)) cycle
97 iis = hecmesh%elem_node_index(icel-1)
98 nn = hecmesh%elem_node_index(icel)-iis
102 nodlocal(j) = hecmesh%elem_node_item (iis+j)
104 ecoord(i,j) = hecmesh%node(3*nodlocal(j)+i-3)
107 ddu(i,j) = hecmat%X(ndof*nodlocal(j)+i-ndof)
108 du(i,j) = fstrsolid%dunode(ndof*nodlocal(j)+i-ndof)
109 total_disp(i,j) = fstrsolid%unode(ndof*nodlocal(j)+i-ndof)
112 if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres > 0 )
then
113 if( iselastoplastic(fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype) .or. &
114 fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype == norton )
then
115 tt0(j)=fstrsolid%last_temp( nodlocal(j) )
118 if( hecmesh%hecmw_flag_initcon == 1 ) tt0(j) = hecmesh%node_init_val_item(nodlocal(j))
119 if( initt>0 ) tt0(j) =
g_initialcnd(initt)%realval(nodlocal(j))
121 ttn(j) = fstrsolid%last_temp( nodlocal(j) )
122 tt(j) = fstrsolid%temperature( nodlocal(j) )
126 isect = hecmesh%section_ID(icel)
127 ihead = hecmesh%section%sect_R_index(isect-1)
128 thick = hecmesh%section%sect_R_item(ihead+1)
129 cdsys_id = hecmesh%section%sect_orien_ID(isect)
130 if( cdsys_id > 0 )
call get_coordsys(cdsys_id, hecmesh, fstrsolid, coords)
133 if( getspacedimension( ic_type ) == 2 ) thick = 1.0d0
135 if( ic_type == 241 .or. ic_type == 242 .or. ic_type == 231 .or. ic_type == 232 .or. ic_type == 2322 )
then
136 call update_c2( ic_type,nn,ecoord(1:3,1:nn),fstrsolid%elements(icel)%gausses(:), &
137 thick,fstrsolid%elements(icel)%iset, &
138 total_disp(1:2,1:nn), ddu(1:2,1:nn), qf(1:nn*ndof), &
139 tt(1:nn), tt0(1:nn), ttn(1:nn) )
141 else if( ic_type == 301 )
then
142 call update_c1( ic_type,nn,ecoord(:,1:nn), thick, total_disp(1:3,1:nn), du(1:3,1:nn), &
143 qf(1:nn*ndof),fstrsolid%elements(icel)%gausses(:) )
145 else if( ic_type == 361 )
then
146 if( fstrsolid%sections(isect)%elemopt361 ==
kel361fi )
then
147 call update_c3( ic_type, nn, ecoord(:,1:nn), total_disp(1:3,1:nn), du(1:3,1:nn), cdsys_id, coords, &
148 qf(1:nn*ndof), fstrsolid%elements(icel)%gausses(:), iter, time, tincr, tt(1:nn), tt0(1:nn), ttn(1:nn) )
149 else if( fstrsolid%sections(isect)%elemopt361 ==
kel361bbar )
then
150 call update_c3d8bbar( ic_type, nn, ecoord(:,1:nn), total_disp(1:3,1:nn), du(1:3,1:nn), cdsys_id, coords, &
151 qf(1:nn*ndof), fstrsolid%elements(icel)%gausses(:), iter, time, tincr, tt(1:nn), tt0(1:nn), ttn(1:nn) )
152 else if( fstrsolid%sections(isect)%elemopt361 ==
kel361ic )
then
153 call update_c3d8ic( ic_type,nn,ecoord(:,1:nn), total_disp(1:3,1:nn), du(1:3,1:nn), ddu(1:3,1:nn), cdsys_id, coords,&
154 qf(1:nn*ndof), fstrsolid%elements(icel)%gausses(:), iter, time, tincr, &
155 fstrsolid%elements(icel)%aux, ddaux(1:3,1:3), tt(1:nn), tt0(1:nn), ttn(1:nn) )
156 fstrsolid%elements(icel)%aux(1:3,1:3) = fstrsolid%elements(icel)%aux(1:3,1:3) + ddaux(1:3,1:3)
157 else if( fstrsolid%sections(isect)%elemopt361 ==
kel361fbar )
then
158 call update_c3d8fbar( ic_type, nn, ecoord(:,1:nn), total_disp(1:3,1:nn), du(1:3,1:nn), cdsys_id, coords, &
159 qf(1:nn*ndof), fstrsolid%elements(icel)%gausses(:), iter, time, tincr, tt(1:nn), tt0(1:nn), ttn(1:nn) )
162 else if (ic_type == 341 .or. ic_type == 351 .or. ic_type == 342 .or. ic_type == 352 .or. ic_type == 362 )
then
163 if( ic_type==341 .and. fstrsolid%sections(isect)%elemopt341 ==
kel341sesns ) cycle
164 call update_c3( ic_type, nn, ecoord(:,1:nn), total_disp(1:3,1:nn), du(1:3,1:nn), cdsys_id, coords, &
165 qf(1:nn*ndof), fstrsolid%elements(icel)%gausses(:), iter, time, tincr, tt(1:nn), tt0(1:nn), ttn(1:nn) )
167 else if( ic_type == 511)
then
168 call update_connector( ic_type,nn,ecoord(:,1:nn), total_disp(1:3,1:nn), du(1:3,1:nn), &
169 qf(1:nn*ndof),fstrsolid%elements(icel)%gausses(:), tincr )
171 else if( ic_type == 611)
then
172 if(
fstrpr%nlgeom )
call update_abort( ic_type, 2 )
173 CALL updatest_beam(ic_type, nn, ecoord, total_disp(1:6,1:nn), du(1:6,1:nn), &
174 & hecmesh%section%sect_R_item(ihead+1:), fstrsolid%elements(icel)%gausses(:), qf(1:nn*ndof))
176 else if( ic_type == 641 )
then
177 if(
fstrpr%nlgeom )
call update_abort( ic_type, 2 )
178 call updatest_beam_641(ic_type, nn, ecoord, total_disp(1:ndof,1:nn), du(1:ndof,1:nn), &
179 & fstrsolid%elements(icel)%gausses(:), hecmesh%section%sect_R_item(ihead+1:), qf(1:nn*ndof))
181 else if( ( ic_type == 741 ) .or. ( ic_type == 743 ) .or. ( ic_type == 731 ) )
then
182 if(
fstrpr%nlgeom )
call update_abort( ic_type, 2 )
183 call updatest_shell_mitc(ic_type, nn, ndof, ecoord(1:3, 1:nn), total_disp(1:ndof,1:nn), du(1:ndof,1:nn), &
184 & fstrsolid%elements(icel)%gausses(:), qf(1:nn*ndof), thick, 0)
186 else if( ic_type == 761 )
then
187 if(
fstrpr%nlgeom )
call update_abort( ic_type, 2 )
188 call updatest_shell_mitc33(731, 3, 6, ecoord(1:3, 1:3), total_disp(1:ndof,1:nn), du(1:ndof,1:nn), &
189 & fstrsolid%elements(icel)%gausses(:), qf(1:nn*ndof), thick, 2)
191 else if( ic_type == 781 )
then
192 if(
fstrpr%nlgeom )
call update_abort( ic_type, 2 )
193 call updatest_shell_mitc33(741, 4, 6, ecoord(1:3, 1:4), total_disp(1:ndof,1:nn), du(1:ndof,1:nn), &
194 & fstrsolid%elements(icel)%gausses(:), qf(1:nn*ndof), thick, 1)
196 else if ( ic_type == 3414 )
then
197 if(fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype /= incomp_newtonian) &
198 &
call update_abort( ic_type, 3, fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype )
200 ( ic_type, nn, ecoord(:,1:nn), total_disp(1:4,1:nn), du(1:4,1:nn), &
201 fstrsolid%elements(icel)%gausses(:) )
204 else if ( ic_type == 881 .or. ic_type == 891 )
then
205 call update_c3_sesns( ic_type, nn, nodlocal, ecoord(:,1:nn), total_disp(1:3,1:nn), du(1:3,1:nn), cdsys_id, coords, &
206 qf(1:nn*ndof), fstrsolid%elements(icel)%gausses(:), time, tincr, tt(1:nn), tt0(1:nn), ttn(1:nn) )
209 write(*, *)
'###ERROR### : Element type not supported for nonlinear static analysis'
210 write(*, *)
' ic_type = ', ic_type
211 call hecmw_abort(hecmw_comm_get_comm())
219 fstrsolid%QFORCE(ndof*(nodlocal(j)-1)+i) = fstrsolid%QFORCE(ndof*(nodlocal(j)-1)+i)+qf(ndof*(j-1)+i)
224 if(
present(strainenergy))
then
225 ndim = getspacedimension( fstrsolid%elements(icel)%etype )
228 tmp = 0.5d0*( fstrsolid%elements(icel)%equiForces(ndim*(j-1)+i)+qf(ndim*(j-1)+i) )*ddu(i,j)
230 strainenergy = strainenergy+tmp
231 fstrsolid%elements(icel)%equiForces(ndim*(j-1)+i) = qf(ndim*(j-1)+i)
244 call hecmw_update_r(hecmesh,fstrsolid%QFORCE,hecmesh%n_node, ndof)
255 type(hecmwst_local_mesh) :: hecmesh
257 real(kind=kreal) :: tincr
258 integer(kind=kint) :: itype, is, iE, ic_type, icel, ngauss, i
260 if(
associated( fstrsolid%temperature ) )
then
261 do i = 1, hecmesh%n_node
262 fstrsolid%last_temp(i) = fstrsolid%temperature(i)
266 do itype = 1, hecmesh%n_elem_type
267 is = hecmesh%elem_type_index(itype-1) + 1
268 ie = hecmesh%elem_type_index(itype )
269 ic_type= hecmesh%elem_type_item(itype)
270 if( ic_type == 301 ) ic_type = 111
271 if( hecmw_is_etype_link(ic_type) ) cycle
272 if( hecmw_is_etype_patch(ic_type) ) cycle
274 ngauss = numofquadpoints( ic_type )
276 if( iselastoplastic( fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype ) )
then
280 elseif( fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype == norton )
then
281 if( tincr>0.d0 )
then
286 elseif( isviscoelastic( fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype ) )
then
287 if( tincr > 0.d0 )
then
295 fstrsolid%elements(icel)%gausses(i)%strain_bak = fstrsolid%elements(icel)%gausses(i)%strain
296 fstrsolid%elements(icel)%gausses(i)%stress_bak = fstrsolid%elements(icel)%gausses(i)%stress
297 fstrsolid%elements(icel)%gausses(i)%strain_energy_bak = fstrsolid%elements(icel)%gausses(i)%strain_energy
302 do i = 1, hecmesh%n_node
303 fstrsolid%QFORCE_bak(i) = fstrsolid%QFORCE(i)
308 subroutine update_abort( ic_type, flag, mtype )
309 integer(kind=kint),
intent(in) :: ic_type
310 integer(kind=kint),
intent(in) :: flag
311 integer(kind=kint),
intent(in),
optional :: mtype
314 write(*,*)
'###ERROR### : Element type not supported for static analysis'
315 else if( flag == 2 )
then
316 write(*,*)
'###ERROR### : Element type not supported for nonlinear static analysis'
317 else if( flag == 3 )
then
318 write(*,*)
'###ERROR### : This element is not supported for this material'
320 write(*,*)
' ic_type = ', ic_type
321 if(
present(mtype) )
write(*,*)
' mtype = ', mtype
322 call hecmw_abort(hecmw_comm_get_comm())
This module provide functions for elastoplastic calculation.
subroutine, public updateepstate(gauss)
Clear elatoplastic state.
This module provides function to calculate to do updates.
subroutine fstr_updatestate(hecMESH, fstrSOLID, tincr)
Update elastiplastic status.
subroutine fstr_updatenewton(hecMESH, hecMAT, fstrSOLID, time, tincr, iter, strainEnergy)
Update displacement, stress, strain and internal forces.
This module defines common data and basic structures for analysis.
integer(kind=kint), parameter kel361bbar
subroutine get_coordsys(cdsys_ID, hecMESH, fstrSOLID, coords)
This subroutine fetch coords defined by local coordinate system.
integer(kind=kint), parameter kel341sesns
integer(kind=kint), parameter kel361fi
integer(kind=kint), parameter kel361ic
type(fstr_param), target fstrpr
GLOBAL VARIABLE INITIALIZED IN FSTR_SETUP.
integer(kind=kint), parameter kel361fbar
type(tinitialcondition), dimension(:), pointer, save g_initialcnd
This modules just summarizes all modules used in static analysis.
This module provides functions for creep calculation.
subroutine updateviscostate(gauss)
Update viscoplastic state.
This module provides functions for viscoelastic calculation.
subroutine updateviscoelasticstate(gauss)
Update viscoplastic state.