FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_dynamic_nlexplicit.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!-------------------------------------------------------------------------------
6
8 use m_fstr
21
22contains
23
24 !C================================================================C
25 !C-- subroutine fstr_solve_LINEAR_DYNAMIC
26 !C================================================================C
27 subroutine fstr_solve_dynamic_nlexplicit(hecMESH,hecMAT,fstrSOLID,fstrEIG &
28 ,fstrDYN,fstrRESULT,fstrPARAM,infoCTChange &
29 ,fstrCPL, restrt_step_num )
30 implicit none
31 type(hecmwst_local_mesh) :: hecMESH
32 type(hecmwst_matrix) :: hecMAT
33 type(fstr_eigen) :: fstrEIG
34 type(fstr_solid) :: fstrSOLID
35 type(hecmwst_result_data) :: fstrRESULT
36 type(fstr_param) :: fstrPARAM
37 type(fstr_dynamic) :: fstrDYN
38 type(fstrst_matrix_contact_lagrange) :: fstrMAT
39 type(fstr_info_contactchange) :: infoCTChange
40 type(fstr_couple) :: fstrCPL !for COUPLE
41 type(hecmwst_matrix), pointer :: hecMATmpc
42 integer(kind=kint), allocatable :: mark(:)
43 integer(kind=kint) :: nnod, ndof, nn, numnp
44 integer(kind=kint) :: i, j, ids, ide, kk
45 integer(kind=kint) :: kkk0, kkk1
46 integer(kind=kint) :: ierror
47 integer(kind=kint) :: iiii5, iexit
48 integer(kind=kint) :: revocap_flag
49 real(kind=kreal), allocatable :: prevb(:)
50 real(kind=kreal) :: a1, a2, a3, b1, b2, b3, c1, c2
51 real(kind=kreal) :: bsize, res
52 real(kind=kreal) :: time_1, time_2
53 integer(kind=kint) :: restrt_step_num
54 real(kind=kreal), parameter :: pi = 3.14159265358979323846d0
55
56 a1 = 0.0d0; a2 = 0.0d0; a3 = 0.0d0; b1 = 0.0d0; b2 = 0.0d0; b3 = 0.0d0
57 c1 = 0.0d0; c2 = 0.0d0
58
59 call hecmw_mpc_mat_init_explicit(hecmesh, hecmat, hecmatmpc)
60
61 hecmat%NDOF=hecmesh%n_dof
62 nnod=hecmesh%n_node
63 ndof=hecmat%NDOF
64 nn=ndof*ndof
65
66 if( fstrparam%fg_couple == 1) then
67 if( fstrparam%fg_couple_type==5 .or. &
68 fstrparam%fg_couple_type==6 ) then
69 allocate( prevb(hecmat%NP*ndof) ,stat=ierror )
70 prevb = 0.0d0
71 if( ierror /= 0 ) then
72 write(idbg,*) 'stop due to allocation error <fstr_solve_NONLINEAR_DYNAMIC, prevB>'
73 write(idbg,*) ' rank = ', hecmesh%my_rank,' ierror = ',ierror
74 call flush(idbg)
75 call hecmw_abort( hecmw_comm_get_comm())
76 endif
77 endif
78 endif
79
80 fstrsolid%dunode(:) =0.d0
81
82 a1 = 1.d0/fstrdyn%t_delta**2
83 a2 = 1.d0/(2.d0*fstrdyn%t_delta)
84
85 call setmass(fstrsolid,hecmesh,hecmat,fstreig)
86 call hecmw_mpc_trans_mass(hecmesh, hecmat, fstreig%mass)
87
88 allocate(mark(hecmat%NP * hecmat%NDOF))
89 call hecmw_mpc_mark_slave(hecmesh, hecmat, mark)
90
91 do j = 1 ,ndof*nnod
92 fstrdyn%VEC1(j) = (a1 + a2 *fstrdyn%ray_m) * fstreig%mass(j)
93 if(mark(j) == 1) fstrdyn%VEC1(j) = 1.d0
94 if(dabs(fstrdyn%VEC1(j)) < 1.0e-20) then
95 if( hecmesh%my_rank == 0 ) then
96 write(*,*) 'stop due to fstrDYN%VEC(j) = 0 , j = ', j
97 write(imsg,*) 'stop due to fstrDYN%VEC(j) = 0 , j = ', j
98 end if
99 call hecmw_abort( hecmw_comm_get_comm())
100 endif
101 end do
102
103 deallocate(mark)
104
105 !C-- output of initial state
106 if( restrt_step_num == 1 ) then
107 do j = 1 ,ndof*nnod
108 fstrdyn%DISP(j,3) = fstrdyn%DISP(j,1) - fstrdyn%VEL (j,1)/(2.d0*a2) + fstrdyn%ACC (j,1)/ (2.d0*a1)
109 fstrdyn%DISP(j,2) = fstrdyn%DISP(j,1) - fstrdyn%VEL (j,1)/ a2 + fstrdyn%ACC (j,1)/ (2.d0*a1) * 4.d0
110 end do
111
112 call fstr_dynamic_output(hecmesh, fstrsolid, fstrdyn, fstrparam)
113 call dynamic_output_monit(hecmesh, fstrparam, fstrdyn, fstreig, fstrsolid)
114 end if
115
116 if( associated( fstrsolid%contacts ) ) then
117 call initialize_contact_output_vectors(fstrsolid,hecmat)
118 call forward_increment_lagrange(1,ndof,fstrdyn%VEC1,hecmesh,fstrsolid,infoctchange,&
119 & fstrdyn%DISP(:,2),fstrsolid%ddunode)
120 endif
121
122 do i= restrt_step_num, fstrdyn%n_step
123
124 fstrdyn%i_step = i
125 fstrdyn%t_curr = fstrdyn%t_delta * i
126
127 !C-- mechanical boundary condition
128 call dynamic_mat_ass_load (hecmesh, hecmat, fstrsolid, fstrdyn, fstrparam)
129 do j=1, hecmesh%n_node* hecmesh%n_dof
130 hecmat%B(j)=hecmat%B(j)-fstrsolid%QFORCE(j)
131 end do
132
133 !C ********************************************************************************
134 !C for couple analysis
135 if( fstrparam%fg_couple == 1 ) then
136 if( fstrparam%fg_couple_type==5 .or. &
137 fstrparam%fg_couple_type==6 ) then
138 do j = 1, hecmat%NP * ndof
139 prevb(j) = hecmat%B(j)
140 enddo
141 endif
142 endif
143 do
144 if( fstrparam%fg_couple == 1 ) then
145 if( fstrparam%fg_couple_type==1 .or. &
146 fstrparam%fg_couple_type==3 .or. &
147 fstrparam%fg_couple_type==5 ) call fstr_rcap_get( fstrcpl )
148 if( fstrparam%fg_couple_first /= 0 ) then
149 bsize = dfloat( i ) / dfloat( fstrparam%fg_couple_first )
150 if( bsize > 1.0 ) bsize = 1.0
151 do kkk0 = 1, fstrcpl%coupled_node_n
152 kkk1 = 3 * kkk0
153 fstrcpl%trac(kkk1-2) = bsize * fstrcpl%trac(kkk1-2)
154 fstrcpl%trac(kkk1-1) = bsize * fstrcpl%trac(kkk1-1)
155 fstrcpl%trac(kkk1 ) = bsize * fstrcpl%trac(kkk1 )
156 enddo
157 endif
158 if( fstrparam%fg_couple_window > 0 ) then
159 j = i - restrt_step_num + 1
160 kk = fstrdyn%n_step - restrt_step_num + 1
161 bsize = 0.5*(1.0-cos(2.0*pi*dfloat(j)/dfloat(kk)))
162 do kkk0 = 1, fstrcpl%coupled_node_n
163 kkk1 = 3 * kkk0
164 fstrcpl%trac(kkk1-2) = bsize * fstrcpl%trac(kkk1-2)
165 fstrcpl%trac(kkk1-1) = bsize * fstrcpl%trac(kkk1-1)
166 fstrcpl%trac(kkk1 ) = bsize * fstrcpl%trac(kkk1 )
167 enddo
168 endif
169 call dynamic_mat_ass_couple( hecmesh, hecmat, fstrsolid, fstrcpl )
170 endif
171 !C ********************************************************************************
172
173 call hecmw_mpc_trans_rhs(hecmesh, hecmat, hecmatmpc)
174
175 do j = 1 ,ndof*nnod
176 hecmatmpc%B(j) = hecmatmpc%B(j) + 2.d0*a1* fstreig%mass(j) * fstrdyn%DISP(j,1) &
177 + (- a1 + a2 * fstrdyn%ray_m) * fstreig%mass(j) * fstrdyn%DISP(j,3)
178 end do
179
180 !C
181 !C-- geometrical boundary condition
182
183 call dynamic_explicit_ass_bc(hecmesh, hecmatmpc, fstrsolid, fstrdyn)
184 call dynamic_explicit_ass_vl(hecmesh, hecmatmpc, fstrsolid, fstrdyn)
185 call dynamic_explicit_ass_ac(hecmesh, hecmatmpc, fstrsolid, fstrdyn)
186 !call dynamic_mat_ass_bc (hecMESH, hecMATmpc, fstrSOLID, fstrDYN, fstrPARAM, fstrMAT)
187 !call dynamic_mat_ass_bc_vl(hecMESH, hecMATmpc, fstrSOLID, fstrDYN, fstrPARAM, fstrMAT)
188 !call dynamic_mat_ass_bc_ac(hecMESH, hecMATmpc, fstrSOLID, fstrDYN, fstrPARAM, fstrMAT)
189
190 ! Finish the calculation
191 do j = 1 ,ndof*nnod
192 hecmatmpc%X(j) = hecmatmpc%B(j) / fstrdyn%VEC1(j)
193 if(dabs(hecmatmpc%X(j)) > 1.0d+5) then
194 if( hecmesh%my_rank == 0 ) then
195 print *, 'Displacement increment too large, please adjust your step size!',i,hecmatmpc%X(j)
196 write(imsg,*) 'Displacement increment too large, please adjust your step size!',i,hecmatmpc%B(j),fstrdyn%VEC1(j)
197 end if
198 call hecmw_abort( hecmw_comm_get_comm())
199 end if
200 end do
201 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
202
203 !C *****************************************************
204 !C for couple analysis
205 if( fstrparam%fg_couple == 1 ) then
206 if( fstrparam%fg_couple_type>1 ) then
207 do j=1, fstrcpl%coupled_node_n
208 if( fstrcpl%dof == 3 ) then
209 kkk0 = j*3
210 kkk1 = fstrcpl%coupled_node(j)*3
211
212 fstrcpl%disp (kkk0-2) = hecmat%X(kkk1-2)
213 fstrcpl%disp (kkk0-1) = hecmat%X(kkk1-1)
214 fstrcpl%disp (kkk0 ) = hecmat%X(kkk1 )
215
216 fstrcpl%velo (kkk0-2) = -b1*fstrdyn%ACC(kkk1-2,1) - b2*fstrdyn%VEL(kkk1-2,1) + &
217 b3*( hecmat%X(kkk1-2) - fstrdyn%DISP(kkk1-2,1) )
218 fstrcpl%velo (kkk0-1) = -b1*fstrdyn%ACC(kkk1-1,1) - b2*fstrdyn%VEL(kkk1-1,1) + &
219 b3*( hecmat%X(kkk1-1) - fstrdyn%DISP(kkk1-1,1) )
220 fstrcpl%velo (kkk0 ) = -b1*fstrdyn%ACC(kkk1,1) - b2*fstrdyn%VEL(kkk1,1) + &
221 b3*( hecmat%X(kkk1) - fstrdyn%DISP(kkk1,1) )
222 fstrcpl%accel(kkk0-2) = -a1*fstrdyn%ACC(kkk1-2,1) - a2*fstrdyn%VEL(kkk1-2,1) + &
223 a3*( hecmat%X(kkk1-2) - fstrdyn%DISP(kkk1-2,1) )
224 fstrcpl%accel(kkk0-1) = -a1*fstrdyn%ACC(kkk1-1,1) - a2*fstrdyn%VEL(kkk1-1,1) + &
225 a3*( hecmat%X(kkk1-1) - fstrdyn%DISP(kkk1-1,1) )
226 fstrcpl%accel(kkk0 ) = -a1*fstrdyn%ACC(kkk1,1) - a2*fstrdyn%VEL(kkk1,1) + &
227 a3*( hecmat%X(kkk1) - fstrdyn%DISP(kkk1,1) )
228 else
229 kkk0 = j*2
230 kkk1 = fstrcpl%coupled_node(j)*2
231
232 fstrcpl%disp (kkk0-1) = hecmat%X(kkk1-1)
233 fstrcpl%disp (kkk0 ) = hecmat%X(kkk1 )
234
235 fstrcpl%velo (kkk0-1) = -b1*fstrdyn%ACC(kkk1-1,1) - b2*fstrdyn%VEL(kkk1-1,1) + &
236 b3*( hecmat%X(kkk1-1) - fstrdyn%DISP(kkk1-1,1) )
237 fstrcpl%velo (kkk0 ) = -b1*fstrdyn%ACC(kkk1,1) - b2*fstrdyn%VEL(kkk1,1) + &
238 b3*( hecmat%X(kkk1) - fstrdyn%DISP(kkk1,1) )
239 fstrcpl%accel(kkk0-1) = -a1*fstrdyn%ACC(kkk1-1,1) - a2*fstrdyn%VEL(kkk1-1,1) + &
240 a3*( hecmat%X(kkk1-1) - fstrdyn%DISP(kkk1-1,1) )
241 fstrcpl%accel(kkk0 ) = -a1*fstrdyn%ACC(kkk1,1) - a2*fstrdyn%VEL(kkk1,1) + &
242 a3*( hecmat%X(kkk1) - fstrdyn%DISP(kkk1,1) )
243 endif
244 end do
245 call fstr_rcap_send( fstrcpl )
246 endif
247
248 select case ( fstrparam%fg_couple_type )
249 case (4)
250 call fstr_rcap_get( fstrcpl )
251 case (5)
252 call fstr_get_convergence( revocap_flag )
253 if( revocap_flag==0 ) then
254 do j = 1, hecmat%NP * ndof
255 hecmat%B(j) = prevb(j)
256 enddo
257 cycle
258 endif
259 case (6)
260 call fstr_get_convergence( revocap_flag )
261 if( revocap_flag==0 ) then
262 do j = 1, hecmat%NP * ndof
263 hecmat%B(j) = prevb(j)
264 enddo
265 call fstr_rcap_get( fstrcpl )
266 cycle
267 else
268 if( i /= fstrdyn%n_step ) call fstr_rcap_get( fstrcpl )
269 endif
270 end select
271 endif
272 exit
273 enddo
274
275 !C *****************************************************
276 !C-- contact corrector
277 !C
278 do j = 1 ,ndof*nnod
279 fstrsolid%unode(j) = fstrdyn%DISP(j,1)
280 fstrsolid%dunode(j) = hecmat%X(j)-fstrdyn%DISP(j,1)
281 enddo
282 if( associated( fstrsolid%contacts ) ) then
283 !call fstr_scan_contact_state( 1, fstrDYN%t_delta, kcaSLAGRANGE, hecMESH, fstrSOLID, infoCTChange )
284 call forward_increment_lagrange(1,ndof,fstrdyn%VEC1,hecmesh,fstrsolid,infoctchange,&
285 & fstrdyn%DISP(:,2),fstrsolid%ddunode)
286 do j = 1 ,ndof*nnod
287 hecmat%X(j) = hecmat%X(j) + fstrsolid%ddunode(j)
288 enddo
289 endif
290
291 !C-- new displacement, velocity and accelaration
292 do j = 1 ,ndof*nnod
293 fstrdyn%ACC (j,1) = a1*(hecmat%X(j) - 2.d0*fstrdyn%DISP(j,1) + fstrdyn%DISP(j,3))
294 fstrdyn%VEL (j,1) = a2*(hecmat%X(j) - fstrdyn%DISP(j,3))
295 fstrsolid%unode(j) = fstrdyn%DISP(j,1)
296 fstrsolid%dunode(j) = hecmat%X(j)-fstrdyn%DISP(j,1)
297 fstrdyn%DISP(j,3) = fstrdyn%DISP(j,1)
298 fstrdyn%DISP(j,1) = hecmat%X(j)
299 hecmat%X(j) = fstrsolid%dunode(j)
300 end do
301
302 ! ----- update strain, stress, and internal force
303 call fstr_updatenewton( hecmesh, hecmat, fstrsolid, fstrdyn%t_curr, fstrdyn%t_delta, 1 )
304
305 do j = 1 ,ndof*nnod
306 fstrsolid%unode(j) = fstrsolid%unode(j) + fstrsolid%dunode(j)
307 end do
308 call fstr_updatestate( hecmesh, fstrsolid, fstrdyn%t_delta )
309
310 if( fstrdyn%restart_nout > 0 .and. &
311 (mod(i,fstrdyn%restart_nout).eq.0 .or. i.eq.fstrdyn%n_step) ) then
312 call fstr_write_restart_dyna_nl(i,hecmesh,fstrsolid,fstrdyn,fstrparam)
313 end if
314 !
315 !C-- output new displacement, velocity and accelaration
316 call fstr_dynamic_output(hecmesh, fstrsolid, fstrdyn, fstrparam)
317 call dynamic_output_monit(hecmesh, fstrparam, fstrdyn, fstreig, fstrsolid)
318
319 enddo
320
321 if( fstrparam%fg_couple == 1) then
322 if( fstrparam%fg_couple_type==5 .or. &
323 fstrparam%fg_couple_type==6 ) then
324 deallocate( prevb ,stat=ierror )
325 if( ierror /= 0 ) then
326 write(idbg,*) 'stop due to deallocation error <fstr_solve_NONLINEAR_DYNAMIC, prevB>'
327 write(idbg,*) ' rank = ', hecmesh%my_rank,' ierror = ',ierror
328 call flush(idbg)
329 call hecmw_abort( hecmw_comm_get_comm())
330 endif
331 endif
332 endif
333
334 call hecmw_mpc_mat_finalize_explicit(hecmesh, hecmat, hecmatmpc)
335
336 end subroutine fstr_solve_dynamic_nlexplicit
337
338
339 subroutine forward_increment_lagrange(cstep,ndof,mmat,hecMESH,fstrSOLID,infoCTChange,wkarray,uc)
340 integer, intent(in) :: cstep
341 integer, intent(in) :: ndof
342 real(kind=kreal), intent(in) :: mmat(:)
343 type( hecmwst_local_mesh ), intent(in) :: hecmesh
344 type(fstr_solid), intent(inout) :: fstrSOLID
345 type(fstr_info_contactchange) :: infoCTChange
346 real(kind=kreal), intent(out) :: wkarray(:)
347 real(kind=kreal), intent(out) :: uc(:)
348 integer :: i, j, k, m, grpid, slave, nn, iSS, sid, etype, iter
349 real(kind=kreal) :: fdum, conv, dlambda, shapefunc(l_max_surface_node), lambda(3)
350
351 call fstr_scan_contact_state_exp( cstep, hecmesh, fstrsolid, infoctchange )
352 if( .not. infoctchange%active ) return
353
354 uc = 0.0d0
355
356 iter = 0
357 do
358 wkarray = 0.0d0
359 do i=1,size(fstrsolid%contacts)
360 do j= 1, size(fstrsolid%contacts(i)%slave)
361 if( fstrsolid%contacts(i)%states(j)%state == contactfree ) cycle
362 if( fstrsolid%contacts(i)%states(j)%distance>epsilon(1.d0) ) then
363 fstrsolid%contacts(i)%states(j)%state = contactfree
364 cycle
365 endif
366 if( iter==0 ) then
367 fstrsolid%contacts(i)%states(j)%multiplier(:) =0.d0
368 fstrsolid%contacts(i)%states(j)%wkdist =0.d0
369 cycle
370 endif
371 slave = fstrsolid%contacts(i)%slave(j)
372
373 sid = fstrsolid%contacts(i)%states(j)%surface
374 nn = size( fstrsolid%contacts(i)%master(sid)%nodes )
375 etype = fstrsolid%contacts(i)%master(sid)%etype
376 call getshapefunc( etype, fstrsolid%contacts(i)%states(j)%lpos(:), shapefunc )
377 wkarray( slave ) = -fstrsolid%contacts(i)%states(j)%multiplier(1)
378 do k=1,nn
379 iss = fstrsolid%contacts(i)%master(sid)%nodes(k)
380 wkarray( iss ) = wkarray( iss ) + shapefunc(k) * fstrsolid%contacts(i)%states(j)%multiplier(1)
381 enddo
382 enddo
383 enddo
384
385 if(iter > 0)then
386 do i=1,size(fstrsolid%contacts)
387 do j= 1, size(fstrsolid%contacts(i)%slave)
388 if( fstrsolid%contacts(i)%states(j)%state == contactfree ) cycle
389 slave = fstrsolid%contacts(i)%slave(j)
390 sid = fstrsolid%contacts(i)%states(j)%surface
391 nn = size( fstrsolid%contacts(i)%master(sid)%nodes )
392 etype = fstrsolid%contacts(i)%master(sid)%etype
393 call getshapefunc( etype, fstrsolid%contacts(i)%states(j)%lpos(:), shapefunc )
394 fstrsolid%contacts(i)%states(j)%wkdist = -wkarray( slave )/mmat( (slave-1)*ndof+1 )
395 do k=1,nn
396 iss = fstrsolid%contacts(i)%master(sid)%nodes(k)
397 fstrsolid%contacts(i)%states(j)%wkdist = fstrsolid%contacts(i)%states(j)%wkdist &
398 + shapefunc(k) * wkarray(iss) / mmat( (iss-1)*ndof+1 )
399 enddo
400 enddo
401 enddo
402 endif
403
404 conv = 0.d0
405 wkarray = 0.d0
406 do i=1,size(fstrsolid%contacts)
407 do j= 1, size(fstrsolid%contacts(i)%slave)
408 if( fstrsolid%contacts(i)%states(j)%state == contactfree ) cycle
409 slave = fstrsolid%contacts(i)%slave(j)
410 sid = fstrsolid%contacts(i)%states(j)%surface
411 nn = size( fstrsolid%contacts(i)%master(sid)%nodes )
412 etype = fstrsolid%contacts(i)%master(sid)%etype
413 call getshapefunc( etype, fstrsolid%contacts(i)%states(j)%lpos(:), shapefunc )
414 fdum = 1.d0/mmat( (slave-1)*ndof+1 )
415 do k=1,nn
416 iss = fstrsolid%contacts(i)%master(sid)%nodes(k)
417 fdum = fdum + shapefunc(k)*shapefunc(k)/mmat( (iss-1)*ndof+1 )
418 enddo
419 dlambda= (fstrsolid%contacts(i)%states(j)%distance-fstrsolid%contacts(i)%states(j)%wkdist) /fdum
420 conv = conv + dlambda*dlambda;
421 fstrsolid%contacts(i)%states(j)%multiplier(1) = fstrsolid%contacts(i)%states(j)%multiplier(1) + dlambda
422 if( fstrsolid%contacts(i)%fcoeff>0.d0 ) then
423 if( fstrsolid%contacts(i)%states(j)%state == contactslip ) then
424 fstrsolid%contacts(i)%states(j)%multiplier(2) = &
425 fstrsolid%contacts(i)%fcoeff * fstrsolid%contacts(i)%states(j)%multiplier(1)
426 else ! stick
427 ! fstrSOLID%contacts(i)%states(j)%multiplier(2) =
428 endif
429 endif
430 lambda = fstrsolid%contacts(i)%states(j)%multiplier(1)* fstrsolid%contacts(i)%states(j)%direction
431 wkarray((slave-1)*ndof+1:(slave-1)*ndof+3) = lambda(:)
432 do k=1,nn
433 iss = fstrsolid%contacts(i)%master(sid)%nodes(k)
434 wkarray((iss-1)*ndof+1:(iss-1)*ndof+3) = wkarray((iss-1)*ndof+1:(iss-1)*ndof+3) -lambda(:)*shapefunc(k)
435 enddo
436 enddo
437 enddo
438 if( dsqrt(conv)<1.d-8 ) exit
439 iter = iter+1
440 enddo
441
442 do i=1,hecmesh%n_node*ndof
443 uc(i) = wkarray(i)/mmat(i)
444 enddo
445 end subroutine forward_increment_lagrange
446
This module contains subroutines for nonlinear explicit dynamic analysis.
subroutine forward_increment_lagrange(cstep, ndof, mmat, hecmesh, fstrsolid, infoctchange, wkarray, uc)
subroutine fstr_solve_dynamic_nlexplicit(hecmesh, hecmat, fstrsolid, fstreig, fstrdyn, fstrresult, fstrparam, infoctchange, fstrcpl, restrt_step_num)
This module provides functions of reconstructing.
This module contains functions to set acceleration boundary condition in dynamic analysis.
subroutine dynamic_explicit_ass_ac(hecmesh, hecmat, fstrsolid, fstrdynamic, iter)
This module contains functions to set velocity boundary condition in dynamic analysis.
subroutine dynamic_explicit_ass_vl(hecmesh, hecmat, fstrsolid, fstrdynamic, iter)
This module contains functions to set displacement boundary condition in dynamic analysis.
subroutine dynamic_explicit_ass_bc(hecmesh, hecmat, fstrsolid, fstrdynamic, iter)
This subroutine setup disp boundary condition.
This module contains functions relates to coupling analysis.
subroutine dynamic_mat_ass_couple(hecmesh, hecmat, fstrsolid, fstrcpl)
This module contains function to set boundary condition of external load in dynamic analysis.
subroutine dynamic_mat_ass_load(hecmesh, hecmat, fstrsolid, fstrdynamic, fstrparam, iter)
This function sets boundary condition of external load.
This module provides functions to output result.
subroutine fstr_dynamic_output(hecmesh, fstrsolid, fstrdynamic, fstrparam)
Output result.
subroutine dynamic_output_monit(hecmesh, fstrparam, fstrdynamic, fstreig, fstrsolid)
Set up lumped mass matrix.
subroutine setmass(fstrsolid, hecmesh, hecmat, fstreig)
subroutine, public fstr_rcap_send(fstrcpl)
subroutine, public fstr_rcap_get(fstrcpl)
subroutine fstr_get_convergence(revocap_flag)
This module provides functions to read in and write out restart fiels.
Definition: fstr_Restart.f90:8
subroutine fstr_write_restart_dyna_nl(cstep, hecmesh, fstrsolid, fstrdynamic, fstrparam, contactnode)
write out restart file for nonlinear dynamic analysis
This module provides function to calcualte to do updates.
Definition: fstr_Update.f90:6
subroutine fstr_updatenewton(hecmesh, hecmat, fstrsolid, time, tincr, iter, strainenergy)
変位/応力・ひずみ/内力のアップデート
Definition: fstr_Update.f90:26
subroutine fstr_updatestate(hecmesh, fstrsolid, tincr)
Update elastiplastic status.
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
integer(kind=kint), parameter idbg
Definition: m_fstr.f90:95
This modules just summarizes all modules used in static analysis.
Definition: static_LIB.f90:6
Structure for Lagrange multiplier-related part of stiffness matrix (Lagrange multiplier-related matri...
Data for coupling analysis.
Definition: m_fstr.f90:580
Data for DYNAMIC ANSLYSIS (fstrDYNAMIC)
Definition: m_fstr.f90:473
Package of data used by Lanczos eigenvalue solver.
Definition: m_fstr.f90:562
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138