23 type (hecmwst_local_mesh) :: hecmesh
24 type (hecmwst_matrix) :: hecmat
26 real(kind=kreal),
intent(in) :: time
27 real(kind=kreal),
intent(in) :: tincr
29 type( tmaterial ),
pointer :: material
31 real(kind=kreal) :: stiffness(20*6, 20*6)
32 integer(kind=kint) :: nodlocal(20)
33 real(kind=kreal) :: tt(20), ecoord(3,20)
34 real(kind=kreal) :: thick
35 integer(kind=kint) :: ndof, itype, is, ie, ic_type, nn, icel, iis, i, j
36 real(kind=kreal) :: u(6,20), du(6,20), coords(3,3), u_prev(6,20)
37 integer :: isect, ihead, cdsys_id
40 call hecmw_mat_clear( hecmat )
43 do itype= 1, hecmesh%n_elem_type
44 is= hecmesh%elem_type_index(itype-1) + 1
45 ie= hecmesh%elem_type_index(itype )
46 ic_type= hecmesh%elem_type_item(itype)
48 if (hecmw_is_etype_link(ic_type)) cycle
49 if (hecmw_is_etype_patch(ic_type)) cycle
51 nn = hecmw_get_max_node(ic_type)
62 iis= hecmesh%elem_node_index(icel-1)
64 nodlocal(j)= hecmesh%elem_node_item (iis+j)
66 ecoord(i,j) = hecmesh%node(3*nodlocal(j)+i-3)
69 du(i,j) = fstrsolid%dunode(ndof*nodlocal(j)+i-ndof)
70 u(i,j) = fstrsolid%unode(ndof*nodlocal(j)+i-ndof) + du(i,j)
71 u_prev(i,j) = fstrsolid%unode(ndof*nodlocal(j)+i-ndof)
73 if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres >0 ) &
74 tt(j)=fstrsolid%temperature( nodlocal(j) )
77 isect = hecmesh%section_ID(icel)
78 cdsys_id = hecmesh%section%sect_orien_ID(isect)
79 if( cdsys_id > 0 )
call get_coordsys(cdsys_id, hecmesh, fstrsolid, coords)
81 material => fstrsolid%elements(icel)%gausses(1)%pMaterial
82 thick = material%variables(m_thick)
83 if( getspacedimension( ic_type )==2 ) thick =1.d0
84 if( ic_type==241 .or. ic_type==242 .or. ic_type==231 .or. ic_type==232 .or. ic_type==2322)
then
85 if( material%nlgeom_flag /= infinitesimal )
call stiffmat_abort( ic_type, 2 )
86 call stf_c2( ic_type,nn,ecoord(1:2,1:nn),fstrsolid%elements(icel)%gausses(:),thick, &
87 stiffness(1:nn*ndof,1:nn*ndof), fstrsolid%elements(icel)%iset, &
90 elseif ( ic_type==301 )
then
91 isect= hecmesh%section_ID(icel)
92 ihead = hecmesh%section%sect_R_index(isect-1)
93 thick = hecmesh%section%sect_R_item(ihead+1)
94 call stf_c1( ic_type,nn,ecoord(:,1:nn),thick,fstrsolid%elements(icel)%gausses(:), &
95 stiffness(1:nn*ndof,1:nn*ndof), u(1:3,1:nn) )
97 elseif ( ic_type==361 )
then
99 if( fstrsolid%sections(isect)%elemopt361 ==
kel361fi )
then
100 if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres >0 )
then
102 ( ic_type, nn, ecoord(:, 1:nn), fstrsolid%elements(icel)%gausses(:), &
103 stiffness(1:nn*ndof, 1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3,1:nn), tt(1:nn) )
106 ( ic_type,nn,ecoord(:, 1:nn),fstrsolid%elements(icel)%gausses(:), &
107 stiffness(1:nn*ndof, 1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3, 1:nn) )
109 else if( fstrsolid%sections(isect)%elemopt361 ==
kel361bbar )
then
110 if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres >0 )
then
112 ( ic_type, nn, ecoord(:, 1:nn), fstrsolid%elements(icel)%gausses(:), &
113 stiffness(1:nn*ndof,1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3, 1:nn), tt(1:nn) )
116 ( ic_type, nn, ecoord(:, 1:nn),fstrsolid%elements(icel)%gausses(:), &
117 stiffness(1:nn*ndof, 1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3, 1:nn) )
119 else if( fstrsolid%sections(isect)%elemopt361 ==
kel361ic )
then
120 if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres >0 )
then
122 ( ic_type, nn, ecoord(:,1:nn), fstrsolid%elements(icel)%gausses(:), &
123 stiffness(1:nn*ndof, 1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3,1:nn), &
124 fstrsolid%elements(icel)%aux, tt(1:nn) )
127 ( ic_type, nn, ecoord(:,1:nn), fstrsolid%elements(icel)%gausses(:), &
128 stiffness(1:nn*ndof, 1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3,1:nn), &
129 fstrsolid%elements(icel)%aux )
131 else if( fstrsolid%sections(isect)%elemopt361 ==
kel361fbar )
then
132 if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres >0 )
then
134 ( ic_type, nn, ecoord(:, 1:nn), fstrsolid%elements(icel)%gausses(:), &
135 stiffness(1:nn*ndof,1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3, 1:nn), tt(1:nn) )
138 ( ic_type, nn, ecoord(:, 1:nn),fstrsolid%elements(icel)%gausses(:), &
139 stiffness(1:nn*ndof, 1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3, 1:nn) )
143 elseif (ic_type==341 .or. ic_type==351 .or. &
144 ic_type==342 .or. ic_type==352 .or. ic_type==362 )
then
145 if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres >0 )
then
147 ( ic_type, nn, ecoord(:, 1:nn), fstrsolid%elements(icel)%gausses(:), &
148 stiffness(1:nn*ndof, 1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3,1:nn), tt(1:nn) )
151 ( ic_type,nn,ecoord(:, 1:nn),fstrsolid%elements(icel)%gausses(:), &
152 stiffness(1:nn*ndof, 1:nn*ndof), cdsys_id, coords, time, tincr, u(1:3, 1:nn) )
155 else if( ic_type == 611)
then
156 if( material%nlgeom_flag /= infinitesimal )
call stiffmat_abort( ic_type, 2 )
157 isect = hecmesh%section_ID(icel)
158 ihead = hecmesh%section%sect_R_index(isect-1)
159 call stf_beam(ic_type, nn, ecoord, hecmesh%section%sect_R_item(ihead+1:), &
160 & material%variables(m_youngs), material%variables(m_poisson), stiffness(1:nn*ndof,1:nn*ndof))
162 else if( ic_type == 641 )
then
163 if( material%nlgeom_flag /= infinitesimal )
call stiffmat_abort( ic_type, 2 )
164 isect = hecmesh%section_ID(icel)
165 ihead = hecmesh%section%sect_R_index(isect-1)
166 call stf_beam_641(ic_type, nn, ecoord, fstrsolid%elements(icel)%gausses(:), &
167 & hecmesh%section%sect_R_item(ihead+1:), stiffness(1:nn*ndof,1:nn*ndof))
169 else if( ( ic_type == 741 ) .or. ( ic_type == 743 ) .or. ( ic_type == 731 ) )
then
170 if( material%nlgeom_flag /= infinitesimal )
call stiffmat_abort( ic_type, 2 )
171 isect = hecmesh%section_ID(icel)
172 ihead = hecmesh%section%sect_R_index(isect-1)
173 thick = hecmesh%section%sect_R_item(ihead+1)
174 call stf_shell_mitc(ic_type, nn, ndof, ecoord(1:3, 1:nn), fstrsolid%elements(icel)%gausses(:), &
175 & stiffness(1:nn*ndof, 1:nn*ndof), thick, 0)
177 else if( ic_type == 761 )
then
178 if( material%nlgeom_flag /= infinitesimal )
call stiffmat_abort( ic_type, 2 )
179 isect = hecmesh%section_ID(icel)
180 ihead = hecmesh%section%sect_R_index(isect-1)
181 thick = hecmesh%section%sect_R_item(ihead+1)
182 call stf_shell_mitc(731, 3, 6, ecoord(1:3, 1:3), fstrsolid%elements(icel)%gausses(:), &
183 & stiffness(1:nn*ndof, 1:nn*ndof), thick, 2)
185 else if( ic_type == 781 )
then
186 if( material%nlgeom_flag /= infinitesimal )
call stiffmat_abort( ic_type, 2 )
187 isect = hecmesh%section_ID(icel)
188 ihead = hecmesh%section%sect_R_index(isect-1)
189 thick = hecmesh%section%sect_R_item(ihead+1)
190 call stf_shell_mitc(741, 4, 6, ecoord(1:3, 1:4), fstrsolid%elements(icel)%gausses(:), &
191 & stiffness(1:nn*ndof, 1:nn*ndof), thick, 1)
193 elseif ( ic_type==3414 )
then
194 if(fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype /= incomp_newtonian)
then
195 write(*, *)
'###ERROR### : This element is not supported for this material'
196 write(*, *)
'ic_type = ', ic_type,
', mtype = ', fstrsolid%elements(icel)%gausses(1)%pMaterial%mtype
197 call hecmw_abort(hecmw_comm_get_comm())
200 ( ic_type, nn, ecoord(:, 1:nn),fstrsolid%elements(icel)%gausses(:), &
201 stiffness(1:nn*ndof, 1:nn*ndof), tincr, u_prev(1:4, 1:nn) )
214 call hecmw_mat_ass_elem(hecmat, nn, nodlocal, stiffness)
224 integer(kind=kint),
intent(in) :: ic_type
225 integer(kind=kint),
intent(in) :: flag
228 write(*,*)
'###ERROR### : Element type not supported for static analysis'
229 else if( flag == 2 )
then
230 write(*,*)
'###ERROR### : Element type not supported for nonlinear static analysis'
232 write(*,*)
' ic_type = ', ic_type
233 call hecmw_abort(hecmw_comm_get_comm())
This module provides function to calcualte tangent stiffness matrix.
subroutine, public fstr_stiffmatrix(hecMESH, hecMAT, fstrSOLID, time, tincr)
接線剛性マトリックスを作成するサブルーチン
subroutine stiffmat_abort(ic_type, flag)
This module defined coomon 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 kel361fi
section control
integer(kind=kint), parameter kel361ic
integer(kind=kint), parameter kel361fbar
This modules just summarizes all modules used in static analysis.
This modules defines a structure to record history dependent parameter in static analysis.