FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
static_output.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
7 implicit none
8contains
9
11 !----------------------------------------------------------------------*
12 subroutine fstr_static_output( cstep, istep, time, hecMESH, fstrSOLID, fstrPARAM, flag, outflag )
13 !----------------------------------------------------------------------*
14 use m_fstr
18 use mcontact
19 integer, intent(in) :: cstep
20 integer, intent(in) :: istep
21 real(kind=kreal), intent(in) :: time
22 type (hecmwST_local_mesh), intent(in) :: hecMESH
23 type (fstr_solid), intent(inout) :: fstrSOLID
24 type (fstr_param ) :: fstrPARAM
25 integer, intent(in) :: flag
26 logical, intent(in) :: outflag
27
28 type ( hecmwST_result_data ) :: fstrRESULT
29 integer(kind=kint) :: i, j, ndof, fnum, is, iE, gid
30
31 ndof = hecmesh%n_dof
32
33 if( fstrsolid%TEMP_ngrp_tot>0 .or. fstrsolid%TEMP_irres>0 ) then
34 if( ndof==3 ) then
35 if( .not. associated(fstrsolid%tnstrain) ) allocate( fstrsolid%tnstrain(hecmesh%n_node*6) )
36 if( .not. associated(fstrsolid%testrain) ) allocate( fstrsolid%testrain(hecmesh%n_elem*6) )
37 else if( ndof==2 ) then
38 if( .not. associated(fstrsolid%tnstrain) ) allocate( fstrsolid%tnstrain(hecmesh%n_node*3) )
39 if( .not. associated(fstrsolid%testrain) ) allocate( fstrsolid%testrain(hecmesh%n_elem*3) )
40 endif
41 endif
42
43 if( ndof==3 ) then
44 call fstr_nodalstress3d( hecmesh, fstrsolid )
45 else if( ndof==2 ) then
46 call fstr_nodalstress2d( hecmesh, fstrsolid )
47 else if( ndof==6 ) then
48 call fstr_nodalstress6d( hecmesh, fstrsolid )
49 endif
50
51 if( associated( fstrsolid%contacts ) ) then
52 call setup_contact_output_variables( hecmesh, fstrsolid, -1 )
53 endif
54
55 if( flag==kststaticeigen ) then
56 if( iresult==1 .and. &
57 (mod(istep,fstrsolid%output_ctrl(3)%freqency)==0 .or. outflag) ) then
58 call fstr_write_result( hecmesh, fstrsolid, fstrparam, istep, time, 1 )
59 endif
60 return
61 endif
62
63 if( iresult==1 .and. &
64 (mod(istep,fstrsolid%output_ctrl(3)%freqency)==0 .or. outflag) ) then
65 call setup_contact_output_variables( hecmesh, fstrsolid, 3 )
66 call fstr_write_result( hecmesh, fstrsolid, fstrparam, istep, time, 0 )
67 endif
68
69 if( ivisual==1 .and. &
70 (mod(istep,fstrsolid%output_ctrl(4)%freqency)==0 .or. outflag) ) then
71
72 call setup_contact_output_variables( hecmesh, fstrsolid, 4 )
73 call fstr_make_result( hecmesh, fstrsolid, fstrresult, istep, time )
74 call fstr2hecmw_mesh_conv( hecmesh )
75 call hecmw_visualize_init
76 call hecmw_visualize( hecmesh, fstrresult, istep )
77 call hecmw_visualize_finalize
78 call hecmw2fstr_mesh_conv( hecmesh )
79 call hecmw_result_free( fstrresult )
80 endif
81
82 if( (mod(istep,fstrsolid%output_ctrl(1)%freqency)==0 .or. outflag) ) then
83 fnum = fstrsolid%output_ctrl(1)%filenum
84 call fstr_static_post( fnum, hecmesh, fstrsolid, istep )
85 endif
86
87 if( fstrsolid%output_ctrl(2)%outinfo%grp_id>0 .and. &
88 (mod(istep,fstrsolid%output_ctrl(2)%freqency)==0 .or. outflag) ) then
89 is = fstrsolid%output_ctrl(2)%outinfo%grp_id
90 fnum = fstrsolid%output_ctrl(2)%filenum
91 do i = hecmesh%node_group%grp_index(is-1)+1, hecmesh%node_group%grp_index(is)
92 ie = hecmesh%node_group%grp_item(i)
93 gid = hecmesh%global_node_ID(ie)
94 write(fnum,'(2i10,1p6e15.7)') istep,gid,(fstrsolid%unode(ndof*(ie-1)+j),j=1,ndof)
95 enddo
96 endif
97
98 end subroutine fstr_static_output
99
101 !----------------------------------------------------------------------*
102 subroutine fstr_static_post( fnum, hecMESH, fstrSOLID, istep )
103 !----------------------------------------------------------------------*
104 use m_fstr
105 integer, intent(in) :: fnum, istep
106 type (hecmwST_local_mesh), intent(in) :: hecMESH
107 type (fstr_solid), intent(in) :: fstrSOLID
108
109 real(kind=kreal) :: umax(6), umin(6), emax(14), emin(14), smax(14), smin(14)
110 real(kind=kreal) :: mmax(1), mmin(1), emmax(1), emmin(1)
111 real(kind=kreal) :: eemax(14), eemin(14), esmax(14), esmin(14)
112
113 real(kind=kreal) :: gumax(6), gumin(6), gemax(14), gemin(14), gsmax(14), gsmin(14)
114 real(kind=kreal) :: gmmax(1), gmmin(1), gemmax(1), gemmin(1)
115 real(kind=kreal) :: geemax(14), geemin(14), gesmax(14), gesmin(14)
116
117 integer(kind=kint) :: IUmax(6), IUmin(6), IEmax(14), IEmin(14), ISmax(14), ISmin(14)
118 integer(kind=kint) :: IMmax(1), IMmin(1), IEMmax(1), IEMmin(1)
119 integer(kind=kint) :: IEEmax(14), IEEmin(14), IESmax(14), IESmin(14)
120 integer(kind=kint) :: i, j, k, ndof, mdof, ID_area
121 integer(kind=kint) :: label(6)
122
123 ndof = hecmesh%n_dof
124
125 write( fnum, '(''#### Result step='',I6)') istep
126 select case (ndof)
127 case (2)
128 mdof = 3
129 label(1)=11;label(2)=22
130 label(3)=12
131 case (3,6)
132 mdof = 6
133 label(1)=11;label(2)=22;label(3)=33;
134 label(4)=12;label(5)=23;label(6)=31;
135 end select
136
137 j = hecmesh%global_node_ID(1)
138 do k = 1, ndof
139 umax(k) = fstrsolid%unode(k); umin(k) = fstrsolid%unode(k)
140 iumax(k)= j; iumin(k)= j
141 enddo
142 do k = 1, mdof
143 emax(k) = fstrsolid%STRAIN(k); emin(k) = fstrsolid%STRAIN(k)
144 smax(k) = fstrsolid%STRESS(k); smin(k) = fstrsolid%STRESS(k)
145 eemax(k) = fstrsolid%ESTRAIN(k); eemin(k) = fstrsolid%ESTRAIN(k)
146 esmax(k) = fstrsolid%ESTRESS(k); esmin(k) = fstrsolid%ESTRESS(k)
147 iemax(k)= j; iemin(k)= j; ismax(k)= j; ismin(k)= j
148 ieemax(k)= j; ieemin(k)= j; iesmax(k)= j; iesmin(k)= j
149 enddo
150 mmax(1) = fstrsolid%MISES(1); mmin(1) = fstrsolid%MISES(1)
151 emmax(1) = fstrsolid%EMISES(1); emmin(1) = fstrsolid%EMISES(1)
152 immax(1)= j; immin(1)= j; iemmax(1)= j; iemmin(1)= j
153
154 !C*** Show Displacement
155 do i = 1, hecmesh%nn_internal
156 if(fstrsolid%is_rot(i)==1)cycle
157 j = hecmesh%global_node_ID(i)
158 do k = 1, ndof
159 if ( fstrsolid%unode(ndof*(i-1)+k) > umax(k) ) then
160 umax(k) = fstrsolid%unode(ndof*(i-1)+k)
161 iumax(k)= j
162 else if( fstrsolid%unode(ndof*(i-1)+k) < umin(k) ) then
163 umin(k) = fstrsolid%unode(ndof*(i-1)+k)
164 iumin(k)= j
165 endif
166 enddo
167 enddo
168 !C*** Nodal Strain / Stress / MISES
169 !C @node
170 do i = 1, hecmesh%nn_internal
171 if(fstrsolid%is_rot(i)==1)cycle
172 j = hecmesh%global_node_ID(i)
173 do k = 1, mdof
174 if ( fstrsolid%STRAIN(mdof*(i-1)+k) > emax(k) ) then
175 emax(k) = fstrsolid%STRAIN(mdof*(i-1)+k)
176 iemax(k)= j
177 else if( fstrsolid%STRAIN(mdof*(i-1)+k) < emin(k) ) then
178 emin(k) = fstrsolid%STRAIN(mdof*(i-1)+k)
179 iemin(k)= j
180 endif
181 if ( fstrsolid%STRESS(mdof*(i-1)+k) > smax(k) ) then
182 smax(k) = fstrsolid%STRESS(mdof*(i-1)+k)
183 ismax(k)= j
184 else if( fstrsolid%STRESS(mdof*(i-1)+k) < smin(k) ) then
185 smin(k) = fstrsolid%STRESS(mdof*(i-1)+k)
186 ismin(k)= j
187 endif
188 enddo
189 if ( fstrsolid%MISES(i) > mmax(1) ) then
190 mmax(1) = fstrsolid%MISES(i)
191 immax(1)= j
192 else if( fstrsolid%MISES(i) < mmin(1) ) then
193 mmin(1) = fstrsolid%MISES(i)
194 immin(1)= j
195 endif
196 enddo
197 !C*** Elemental Strain / STRESS
198 do i = 1, hecmesh%n_elem
199 id_area = hecmesh%elem_ID(i*2)
200 if( id_area==hecmesh%my_rank ) then
201 j = hecmesh%global_elem_ID(i)
202 do k = 1, mdof
203 if( fstrsolid%ESTRAIN(mdof*(i-1)+k) > eemax(k) ) then
204 eemax(k) = fstrsolid%ESTRAIN(mdof*(i-1)+k)
205 ieemax(k)= j
206 else if( fstrsolid%ESTRAIN(mdof*(i-1)+k) < eemin(k) ) then
207 eemin(k) = fstrsolid%ESTRAIN(mdof*(i-1)+k)
208 ieemin(k)= j
209 endif
210 if( fstrsolid%ESTRESS(mdof*(i-1)+k) > esmax(k) ) then
211 esmax(k) = fstrsolid%ESTRESS(mdof*(i-1)+k)
212 iesmax(k)= j
213 else if( fstrsolid%ESTRESS(mdof*(i-1)+k) < esmin(k) ) then
214 esmin(k) = fstrsolid%ESTRESS(mdof*(i-1)+k)
215 iesmin(k)= j
216 endif
217 enddo
218 if( fstrsolid%EMISES(i) > emmax(1) ) then
219 emmax(1) = fstrsolid%EMISES(i)
220 iemmax(1)= j
221 else if( fstrsolid%EMISES(i) < emmin(1) ) then
222 emmin(1) = fstrsolid%EMISES(i)
223 iemmin(1)= j
224 endif
225 endif
226 enddo
227
228
229 write(ilog,*) '##### Local Summary @Node :Max/IdMax/Min/IdMin####'
230 do i = 1, ndof; write(ilog,1029) ' //U',i, ' ',umax(i),iumax(i),umin(i),iumin(i); end do
231 do i = 1, mdof; write(ilog,1029) ' //E',label(i),' ',emax(i),iemax(i),emin(i),iemin(i); end do
232 do i = 1, mdof; write(ilog,1029) ' //S',label(i),' ',smax(i),ismax(i),smin(i),ismin(i); end do
233 write(ilog,1009) '//SMS ' ,mmax(1),immax(1),mmin(1),immin(1)
234 write(ilog,*) '##### Local Summary @Element :Max/IdMax/Min/IdMin####'
235 do i = 1, mdof; write(ilog,1029) ' //E',label(i),' ',eemax(i),ieemax(i),eemin(i),ieemin(i); end do
236 do i = 1, mdof; write(ilog,1029) ' //S',label(i),' ',esmax(i),iesmax(i),esmin(i),iesmin(i); end do
237 write(ilog,1009) '//SMS ' ,emmax(1),iemmax(1),emmin(1),iemmin(1)
238
239 !C*** Show Summary
240 gumax = umax; gumin = umin;
241 gemax = emax; gemin = emin; geemax = eemax; geemin = eemin;
242 gsmax = smax; gsmin = smin; gesmax = esmax; gesmin = esmin;
243 gmmax = mmax; gmmin = mmin; gemmax = emmax; gemmin = emmin;
244
245 call hecmw_allreduce_r(hecmesh,gumax,ndof,hecmw_max)
246 call hecmw_allreduce_r(hecmesh,gumin,ndof,hecmw_min)
247 call hecmw_allreduce_r(hecmesh,gemax,mdof,hecmw_max)
248 call hecmw_allreduce_r(hecmesh,gemin,mdof,hecmw_min)
249 call hecmw_allreduce_r(hecmesh,gsmax,mdof,hecmw_max)
250 call hecmw_allreduce_r(hecmesh,gsmin,mdof,hecmw_min)
251 call hecmw_allreduce_r(hecmesh,gmmax,1,hecmw_max)
252 call hecmw_allreduce_r(hecmesh,gmmin,1,hecmw_min)
253 call hecmw_allreduce_r(hecmesh,geemax,mdof,hecmw_max)
254 call hecmw_allreduce_r(hecmesh,geemin,mdof,hecmw_min)
255 call hecmw_allreduce_r(hecmesh,gesmax,mdof,hecmw_max)
256 call hecmw_allreduce_r(hecmesh,gesmin,mdof,hecmw_min)
257 call hecmw_allreduce_r(hecmesh,gemmax,1,hecmw_max)
258 call hecmw_allreduce_r(hecmesh,gemmin,1,hecmw_min)
259
260 do i=1,ndof
261 if(gumax(i) > umax(i)) iumax(i) = 0
262 if(gumin(i) < umin(i)) iumin(i) = 0
263 enddo
264 do i=1,mdof
265 if(gemax(i) > emax(i)) iemax(i) = 0
266 if(gsmax(i) > smax(i)) ismax(i) = 0
267 if(geemax(i) > eemax(i)) ieemax(i) = 0
268 if(gesmax(i) > esmax(i)) iesmax(i) = 0
269 if(gemin(i) < emin(i)) iemin(i) = 0
270 if(gsmin(i) < smin(i)) ismin(i) = 0
271 if(geemin(i) < eemin(i)) ieemin(i) = 0
272 if(gesmin(i) < esmin(i)) iesmin(i) = 0
273 enddo
274 do i=1,1
275 if(gmmax(i) > mmax(i)) immax(i) = 0
276 if(gemmax(i) > emmax(i)) iemmax(i) = 0
277 if(gmmin(i) < mmin(i)) immin(i) = 0
278 if(gemmin(i) < emmin(i)) iemmin(i) = 0
279 enddo
280
281 call hecmw_allreduce_i(hecmesh,iumax,ndof,hecmw_max)
282 call hecmw_allreduce_i(hecmesh,iumin,ndof,hecmw_max)
283 call hecmw_allreduce_i(hecmesh,iemax,mdof,hecmw_max)
284 call hecmw_allreduce_i(hecmesh,iemin,mdof,hecmw_max)
285 call hecmw_allreduce_i(hecmesh,ismax,mdof,hecmw_max)
286 call hecmw_allreduce_i(hecmesh,ismin,mdof,hecmw_max)
287 call hecmw_allreduce_i(hecmesh,immax,1,hecmw_max)
288 call hecmw_allreduce_i(hecmesh,immin,1,hecmw_max)
289 call hecmw_allreduce_i(hecmesh,ieemax,mdof,hecmw_max)
290 call hecmw_allreduce_i(hecmesh,ieemin,mdof,hecmw_max)
291 call hecmw_allreduce_i(hecmesh,iesmax,mdof,hecmw_max)
292 call hecmw_allreduce_i(hecmesh,iesmin,mdof,hecmw_max)
293 call hecmw_allreduce_i(hecmesh,iemmax,1,hecmw_max)
294 call hecmw_allreduce_i(hecmesh,iemmin,1,hecmw_max)
295
296 if( hecmesh%my_rank==0 ) then
297 write(ilog,*) '##### Global Summary @Node :Max/IdMax/Min/IdMin####'
298 do i = 1, ndof; write(ilog,1029) ' //U',i, ' ',gumax(i),iumax(i),gumin(i),iumin(i); end do
299 do i = 1, mdof; write(ilog,1029) ' //E',label(i),' ',gemax(i),iemax(i),gemin(i),iemin(i); end do
300 do i = 1, mdof; write(ilog,1029) ' //S',label(i),' ',gsmax(i),ismax(i),gsmin(i),ismin(i); end do
301 write(ilog,1009) '//SMS ' ,gmmax(1),immax(1),gmmin(1),immin(1)
302 write(ilog,*) '##### Global Summary @Element :Max/IdMax/Min/IdMin####'
303 do i = 1, mdof; write(ilog,1029) ' //E',label(i),' ',geemax(i),ieemax(i),geemin(i),ieemin(i); end do
304 do i = 1, mdof; write(ilog,1029) ' //S',label(i),' ',gesmax(i),iesmax(i),gesmin(i),iesmin(i); end do
305 write(ilog,1009) '//SMS ' ,gemmax(1),iemmax(1),gemmin(1),iemmin(1)
306 endif
307
308 1009 format(a7,1pe12.4,i10,1pe12.4,i10)
309 1029 format(a,i0,a,1pe12.4,i10,1pe12.4,i10)
310 end subroutine fstr_static_post
311
312end module m_static_output
This module provides functions to caluclation nodal stress.
subroutine fstr_nodalstress2d(hecmesh, fstrsolid)
Calculate NODAL STRESS of plane elements.
subroutine fstr_nodalstress6d(hecmesh, fstrsolid)
Calculate NODAL STRESS of shell elements.
subroutine fstr_nodalstress3d(hecmesh, fstrsolid)
Calculate NODAL STRESS of solid elements.
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint), pointer iresult
Definition: m_fstr.f90:106
integer(kind=kint), parameter ilog
FILE HANDLER.
Definition: m_fstr.f90:91
integer(kind=kint), parameter kststaticeigen
Definition: m_fstr.f90:42
integer(kind=kint), pointer ivisual
Definition: m_fstr.f90:107
HECMW to FSTR Mesh Data Converter. Convering Conectivity of Element Type 232, 342 and 352.
subroutine hecmw2fstr_mesh_conv(hecmesh)
subroutine fstr2hecmw_mesh_conv(hecmesh)
This module provide a function to prepare output of static analysis.
Definition: make_result.f90:6
subroutine, public fstr_make_result(hecmesh, fstrsolid, fstrresult, istep, time, fstrdynamic)
MAKE RESULT for static and dynamic analysis (WITHOUT ELEMENTAL RESULTS) -----------------------------...
subroutine, public setup_contact_output_variables(hecmesh, fstrsolid, phase)
subroutine, public fstr_write_result(hecmesh, fstrsolid, fstrparam, istep, time, flag, fstrdynamic)
OUTPUT result file for static and dynamic analysis.
Definition: make_result.f90:23
This module provides functions to output result.
subroutine fstr_static_output(cstep, istep, time, hecmesh, fstrsolid, fstrparam, flag, outflag)
Output result.
subroutine fstr_static_post(fnum, hecmesh, fstrsolid, istep)
Summarizer of output data which prints out max and min output values.
This module provides functions to calcualte contact stiff matrix.
Definition: fstr_contact.f90:6