FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_dynamic_nlimplicit.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
20
21 !-------- for couple -------
24
25contains
26
28
29 subroutine fstr_solve_dynamic_nlimplicit(cstep, hecMESH,hecMAT,fstrSOLID,fstrEIG, &
30 fstrDYNAMIC,fstrRESULT,fstrPARAM,fstrCPL, restrt_step_num )
31 implicit none
32 !C-- global variable
33 integer, intent(in) :: cstep
34 type(hecmwst_local_mesh) :: hecMESH
35 type(hecmwst_matrix) :: hecMAT
36 type(fstr_eigen) :: fstrEIG
37 type(fstr_solid) :: fstrSOLID
38 type(hecmwst_result_data) :: fstrRESULT
39 type(fstr_param) :: fstrPARAM
40 type(fstr_dynamic) :: fstrDYNAMIC
41 type(fstrst_matrix_contact_lagrange) :: fstrMAT
42 type(fstr_couple) :: fstrCPL !for COUPLE
43
44 !C-- local variable
45 type(hecmwst_local_mesh), pointer :: hecMESHmpc
46 type(hecmwst_matrix), pointer :: hecMATmpc
47 type(hecmwst_matrix), pointer :: hecMAT0
48 integer(kind=kint) :: nnod, ndof, numnp, nn
49 integer(kind=kint) :: i, j, ids, ide, ims, ime, kk, idm, imm
50 integer(kind=kint) :: iter
51 integer(kind=kint) :: iiii5, iexit
52 integer(kind=kint) :: revocap_flag
53 integer(kind=kint) :: kkk0, kkk1
54 integer(kind=kint) :: restrt_step_num
55 integer(kind=kint) :: n_node_global
56 integer(kind=kint) :: ierr
57
58 real(kind=kreal) :: a1, a2, a3, b1, b2, b3, c1, c2
59 real(kind=kreal) :: bsize, res, resb
60 real(kind=kreal) :: time_1, time_2
61 real(kind=kreal), parameter :: pi = 3.14159265358979323846d0
62 real(kind=kreal), allocatable :: coord(:)
63
64 iexit = 0
65 resb = 0.0d0
66
67 call hecmw_mpc_mat_init(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
68 nullify(hecmat0)
69
70 ! sum of n_node among all subdomains (to be used to calc res)
71 n_node_global = hecmesh%nn_internal
72 call hecmw_allreduce_i1(hecmesh,n_node_global,hecmw_sum)
73
74 hecmat%NDOF=hecmesh%n_dof
75 nnod=hecmesh%n_node
76 ndof=hecmat%NDOF
77 nn =ndof*ndof
78
79 allocate(coord(hecmesh%n_node*ndof))
80
81 !!-- initial value
82 time_1 = hecmw_wtime()
83
84 !C-- check parameters
85 if(dabs(fstrdynamic%beta) < 1.0e-20) then
86 if( hecmesh%my_rank == 0 ) then
87 write(imsg,*) 'stop due to Newmark-beta = 0'
88 endif
89 call hecmw_abort( hecmw_comm_get_comm())
90 endif
91
92 !C-- matrix [M] lumped mass matrix
93 if(fstrdynamic%idx_mas == 1) then
94 call setmass(fstrsolid,hecmesh,hecmat,fstreig)
95
96 !C-- consistent mass matrix
97 else if(fstrdynamic%idx_mas == 2) then
98 if( hecmesh%my_rank .eq. 0 ) then
99 write(imsg,*) 'stop: consistent mass matrix is not yet available !'
100 endif
101 call hecmw_abort( hecmw_comm_get_comm())
102 endif
103
104 hecmat%Iarray(98) = 1 !Assmebly complete
105 hecmat%Iarray(97) = 1 !Need numerical factorization
106
107 !C-- time step loop
108 a1 = 0.5d0/fstrdynamic%beta - 1.0d0
109 a2 = 1.0d0/(fstrdynamic%beta*fstrdynamic%t_delta)
110 a3 = 1.0d0/(fstrdynamic%beta*fstrdynamic%t_delta*fstrdynamic%t_delta)
111 b1 = (0.5d0*fstrdynamic%ganma/fstrdynamic%beta - 1.0d0 )*fstrdynamic%t_delta
112 b2 = fstrdynamic%ganma/fstrdynamic%beta - 1.0d0
113 b3 = fstrdynamic%ganma/(fstrdynamic%beta*fstrdynamic%t_delta)
114 c1 = 1.0d0 + fstrdynamic%ray_k*b3
115 c2 = a3 + fstrdynamic%ray_m*b3
116
117 !C-- output of initial state
118 if( restrt_step_num == 1 ) then
119 call fstr_dynamic_output(hecmesh, fstrsolid, fstrdynamic, fstrparam)
120 call dynamic_output_monit(hecmesh, fstrparam, fstrdynamic, fstreig, fstrsolid)
121 endif
122
123 fstrdynamic%VEC3(:) =0.0d0
124 hecmat%X(:) =0.0d0
125
126 !! step = 1,2,....,fstrDYNAMIC%n_step
127 do i= restrt_step_num, fstrdynamic%n_step
128 if(ndof == 4 .and. hecmesh%my_rank==0) write(*,'(a,i5)')"iter: ",i
129
130 fstrdynamic%i_step = i
131 fstrdynamic%t_curr = fstrdynamic%t_delta * i
132
133 if(hecmesh%my_rank==0) then
134 !write(*,'('' time step='',i10,'' time='',1pe13.4e3)') i,fstrDYNAMIC%t_curr
135 write(ista,'('' time step='',i10,'' time='',1pe13.4e3)') i,fstrdynamic%t_curr
136 endif
137
138 do j = 1 ,ndof*nnod
139 fstrdynamic%VEC1(j) = a1*fstrdynamic%ACC(j,1) + a2*fstrdynamic%VEL(j,1)
140 fstrdynamic%VEC2(j) = b1*fstrdynamic%ACC(j,1) + b2*fstrdynamic%VEL(j,1)
141 enddo
142
143 !C ********************************************************************************
144 !C for couple analysis
145 do
146 fstrsolid%dunode(:) =0.d0
147 ! call fstr_UpdateEPState( hecMESH, fstrSOLID )
148
149 if( fstrparam%fg_couple == 1) then
150 if( fstrparam%fg_couple_type==1 .or. &
151 fstrparam%fg_couple_type==3 .or. &
152 fstrparam%fg_couple_type==5 ) call fstr_rcap_get( fstrcpl )
153 endif
154
155 do iter = 1, fstrsolid%step_ctrl(cstep)%max_iter
156 if (fstrparam%nlgeom) then
157 call fstr_stiffmatrix( hecmesh, hecmat, fstrsolid, fstrdynamic%t_curr, fstrdynamic%t_delta )
158 else
159 if (.not. associated(hecmat0)) then
160 call fstr_stiffmatrix( hecmesh, hecmat, fstrsolid, fstrdynamic%t_curr, fstrdynamic%t_delta )
161 allocate(hecmat0)
162 call hecmw_mat_init(hecmat0)
163 call hecmw_mat_copy_profile(hecmat, hecmat0)
164 call hecmw_mat_copy_val(hecmat, hecmat0)
165 else
166 call hecmw_mat_copy_val(hecmat0, hecmat)
167 endif
168 endif
169
170 if( fstrdynamic%ray_k/=0.d0 .or. fstrdynamic%ray_m/=0.d0 ) then
171 do j = 1 ,ndof*nnod
172 hecmat%X(j) = fstrdynamic%VEC2(j) - b3*fstrsolid%dunode(j)
173 enddo
174 endif
175 if( fstrdynamic%ray_k/=0.d0 ) then
176 if( hecmesh%n_dof == 3 ) then
177 call hecmw_matvec (hecmesh, hecmat, hecmat%X, fstrdynamic%VEC3)
178 else if( hecmesh%n_dof == 2 ) then
179 call hecmw_matvec (hecmesh, hecmat, hecmat%X, fstrdynamic%VEC3)
180 else if( hecmesh%n_dof == 6 ) then
181 call matvec(fstrdynamic%VEC3, hecmat%X, hecmat, ndof, hecmat%D, hecmat%AU, hecmat%AL)
182 endif
183 endif
184
185 !C-- mechanical boundary condition
186 call dynamic_mat_ass_load (hecmesh, hecmat, fstrsolid, fstrdynamic, fstrparam, iter)
187 do j=1, hecmesh%n_node* hecmesh%n_dof
188 hecmat%B(j) = hecmat%B(j)- fstrsolid%QFORCE(j) + fstreig%mass(j)*( fstrdynamic%VEC1(j)-a3*fstrsolid%dunode(j) &
189 + fstrdynamic%ray_m* hecmat%X(j) ) + fstrdynamic%ray_k*fstrdynamic%VEC3(j)
190 enddo
191
192 !C ********************************************************************************
193 !C for couple analysis
194 if( fstrparam%fg_couple == 1) then
195 if( fstrparam%fg_couple_first /= 0 ) then
196 bsize = dfloat( i ) / dfloat( fstrparam%fg_couple_first )
197 if( bsize > 1.0 ) bsize = 1.0
198 do kkk0 = 1, fstrcpl%coupled_node_n
199 kkk1 = 3 * kkk0
200 fstrcpl%trac(kkk1-2) = bsize * fstrcpl%trac(kkk1-2)
201 fstrcpl%trac(kkk1-1) = bsize * fstrcpl%trac(kkk1-1)
202 fstrcpl%trac(kkk1 ) = bsize * fstrcpl%trac(kkk1 )
203 enddo
204 endif
205 if( fstrparam%fg_couple_window > 0 ) then
206 j = i - restrt_step_num + 1
207 kk = fstrdynamic%n_step - restrt_step_num + 1
208 bsize = 0.5*(1.0-cos(2.0*pi*dfloat(j)/dfloat(kk)))
209 do kkk0 = 1, fstrcpl%coupled_node_n
210 kkk1 = 3 * kkk0
211 fstrcpl%trac(kkk1-2) = bsize * fstrcpl%trac(kkk1-2)
212 fstrcpl%trac(kkk1-1) = bsize * fstrcpl%trac(kkk1-1)
213 fstrcpl%trac(kkk1 ) = bsize * fstrcpl%trac(kkk1 )
214 enddo
215 endif
216 call dynamic_mat_ass_couple( hecmesh, hecmat, fstrsolid, fstrcpl )
217 endif
218
219 do j = 1 ,nn*hecmat%NP
220 hecmat%D(j) = c1* hecmat%D(j)
221 enddo
222 do j = 1 ,nn*hecmat%NPU
223 hecmat%AU(j) = c1* hecmat%AU(j)
224 enddo
225 do j = 1 ,nn*hecmat%NPL
226 hecmat%AL(j) = c1*hecmat%AL(j)
227 enddo
228 do j=1,nnod
229 do kk=1,ndof
230 idm = nn*(j-1)+1 + (ndof+1)*(kk-1)
231 imm = ndof*(j-1) + kk
232 hecmat%D(idm) = hecmat%D(idm) + c2*fstreig%mass(imm)
233 enddo
234 enddo
235
236 !C-- geometrical boundary condition
237 call hecmw_mpc_mat_ass(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
238 call hecmw_mpc_trans_rhs(hecmesh, hecmat, hecmatmpc)
239 call dynamic_mat_ass_bc (hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, iter)
240 call dynamic_mat_ass_bc_vl(hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, iter)
241 call dynamic_mat_ass_bc_ac(hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, iter)
242
243 !C-- RHS LOAD VECTOR CHECK
244 numnp=hecmatmpc%NP
245 call hecmw_innerproduct_r(hecmesh, ndof, hecmatmpc%B, hecmatmpc%B, bsize)
246
247 if(iter == 1)then
248 resb = bsize
249 endif
250
251 !res = dsqrt(bsize)/n_node_global
252 res = dsqrt(bsize/resb)
253 if( fstrparam%nlgeom .and. ndof /= 4 ) then
254 if(hecmesh%my_rank==0) write(*,'(a,i5,a,1pe12.4)')"iter: ",iter,", res: ",res
255 if(hecmesh%my_rank==0) write(ista,'(''iter='',I5,''- Residual'',E15.7)')iter,res
256 if( res<fstrsolid%step_ctrl(cstep)%converg ) exit
257 endif
258
259 !C-- linear solver [A]{X} = {B}
260 hecmatmpc%X = 0.0d0
261 if( iexit .ne. 1 ) then
262 if( fstrparam%nlgeom ) then
263 if( iter == 1 ) then
264 hecmatmpc%Iarray(97) = 2 !Force numerical factorization
265 else
266 hecmatmpc%Iarray(97) = 1 !Need numerical factorization
267 endif
268 call fstr_set_current_config_to_mesh(hecmeshmpc,fstrsolid,coord)
269 endif
270 call solve_lineq(hecmeshmpc,hecmatmpc)
271 if( fstrparam%nlgeom ) then
272 call fstr_recover_initial_config_to_mesh(hecmeshmpc,fstrsolid,coord)
273 endif
274 endif
275 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
276
277 do j=1,hecmesh%n_node*ndof
278 fstrsolid%dunode(j) = fstrsolid%dunode(j)+hecmat%X(j)
279 enddo
280 ! ----- update the strain, stress, and internal force
281 call fstr_updatenewton( hecmesh, hecmat, fstrsolid, fstrdynamic%t_curr, &
282 & fstrdynamic%t_delta, iter, fstrdynamic%strainEnergy )
283
284 if(.not. fstrparam%nlgeom) exit
285 if(ndof == 4) exit
286 enddo
287
288 ! ----- not convergence
289 if( iter>fstrsolid%step_ctrl(cstep)%max_iter ) then
290 if( hecmesh%my_rank==0) then
291 write(ilog,*) '### Fail to Converge : at step=', i
292 write(ista,*) '### Fail to Converge : at step=', i
293 write( *,*) ' ### Fail to Converge : at step=', i
294 endif
295 stop
296 endif
297
298 !C *****************************************************
299 !C for couple analysis
300 if( fstrparam%fg_couple == 1 ) then
301 if( fstrparam%fg_couple_type>1 ) then
302 do j=1, fstrcpl%coupled_node_n
303 if( fstrcpl%dof == 3 ) then
304 kkk0 = j*3
305 kkk1 = fstrcpl%coupled_node(j)*3
306
307 fstrcpl%disp (kkk0-2) = fstrsolid%unode(kkk1-2) + fstrsolid%dunode(kkk1-2)
308 fstrcpl%disp (kkk0-1) = fstrsolid%unode(kkk1-1) + fstrsolid%dunode(kkk1-1)
309 fstrcpl%disp (kkk0 ) = fstrsolid%unode(kkk1 ) + fstrsolid%dunode(kkk1 )
310
311 fstrcpl%velo (kkk0-2) = -b1*fstrdynamic%ACC(kkk1-2,1) - b2*fstrdynamic%VEL(kkk1-2,1) + &
312 b3*fstrsolid%dunode(kkk1-2)
313 fstrcpl%velo (kkk0-1) = -b1*fstrdynamic%ACC(kkk1-1,1) - b2*fstrdynamic%VEL(kkk1-1,1) + &
314 b3*fstrsolid%dunode(kkk1-1)
315 fstrcpl%velo (kkk0 ) = -b1*fstrdynamic%ACC(kkk1,1) - b2*fstrdynamic%VEL(kkk1,1) + &
316 b3*fstrsolid%dunode(kkk1)
317 fstrcpl%accel(kkk0-2) = -a1*fstrdynamic%ACC(kkk1-2,1) - a2*fstrdynamic%VEL(kkk1-2,1) + &
318 a3*fstrsolid%dunode(kkk1-2)
319 fstrcpl%accel(kkk0-1) = -a1*fstrdynamic%ACC(kkk1-1,1) - a2*fstrdynamic%VEL(kkk1-1,1) + &
320 a3*fstrsolid%dunode(kkk1-1)
321 fstrcpl%accel(kkk0 ) = -a1*fstrdynamic%ACC(kkk1,1) - a2*fstrdynamic%VEL(kkk1,1) + &
322 a3*fstrsolid%dunode(kkk1)
323 else
324 kkk0 = j*2
325 kkk1 = fstrcpl%coupled_node(j)*2
326
327 fstrcpl%disp (kkk0-1) = fstrsolid%unode(kkk1-1) + fstrsolid%dunode(kkk1-1)
328 fstrcpl%disp (kkk0 ) = fstrsolid%unode(kkk1 ) + fstrsolid%dunode(kkk1 )
329
330 fstrcpl%velo (kkk0-1) = -b1*fstrdynamic%ACC(kkk1-1,1) - b2*fstrdynamic%VEL(kkk1-1,1) + &
331 b3*fstrsolid%dunode(kkk1-1)
332 fstrcpl%velo (kkk0 ) = -b1*fstrdynamic%ACC(kkk1,1) - b2*fstrdynamic%VEL(kkk1,1) + &
333 b3*fstrsolid%dunode(kkk1)
334 fstrcpl%accel(kkk0-1) = -a1*fstrdynamic%ACC(kkk1-1,1) - a2*fstrdynamic%VEL(kkk1-1,1) + &
335 a3*fstrsolid%dunode(kkk1-1)
336 fstrcpl%accel(kkk0 ) = -a1*fstrdynamic%ACC(kkk1,1) - a2*fstrdynamic%VEL(kkk1,1) + &
337 a3*fstrsolid%dunode(kkk1)
338 endif
339 enddo
340 call fstr_rcap_send( fstrcpl )
341 endif
342
343 select case ( fstrparam%fg_couple_type )
344 case (4)
345 call fstr_rcap_get( fstrcpl )
346 case (5)
347 call fstr_get_convergence( revocap_flag )
348 if( revocap_flag==0 ) cycle
349 case (6)
350 call fstr_get_convergence( revocap_flag )
351 if( revocap_flag==0 ) then
352 call fstr_rcap_get( fstrcpl )
353 cycle
354 else
355 if( i /= fstrdynamic%n_step ) call fstr_rcap_get( fstrcpl )
356 endif
357 end select
358 endif
359 exit
360 enddo
361 !C *****************************************************
362
363 !C-- new displacement, velocity and accelaration
364 fstrdynamic%kineticEnergy = 0.0d0
365 do j = 1 ,ndof*nnod
366 fstrdynamic%ACC (j,2) = -a1*fstrdynamic%ACC(j,1) - a2*fstrdynamic%VEL(j,1) + &
367 a3*fstrsolid%dunode(j)
368 fstrdynamic%VEL (j,2) = -b1*fstrdynamic%ACC(j,1) - b2*fstrdynamic%VEL(j,1) + &
369 b3*fstrsolid%dunode(j)
370 fstrdynamic%ACC (j,1) = fstrdynamic%ACC (j,2)
371 fstrdynamic%VEL (j,1) = fstrdynamic%VEL (j,2)
372
373 fstrsolid%unode(j) = fstrsolid%unode(j)+fstrsolid%dunode(j)
374 fstrdynamic%DISP(j,2) = fstrsolid%unode(j)
375
376 fstrdynamic%kineticEnergy = fstrdynamic%kineticEnergy + &
377 0.5d0*fstreig%mass(j)*fstrdynamic%VEL(j,2)*fstrdynamic%VEL(j,2)
378 enddo
379
380 !--- Restart info
381 if( fstrdynamic%restart_nout > 0 .and. &
382 (mod(i,fstrdynamic%restart_nout).eq.0 .or. i.eq.fstrdynamic%n_step) ) then
383 call fstr_write_restart_dyna_nl(i,hecmesh,fstrsolid,fstrdynamic,fstrparam)
384 endif
385
386 !C-- output new displacement, velocity and accelaration
387 call fstr_dynamic_output(hecmesh, fstrsolid, fstrdynamic, fstrparam)
388
389 !C-- output result of monitoring node
390 call dynamic_output_monit(hecmesh, fstrparam, fstrdynamic, fstreig, fstrsolid)
391 call fstr_updatestate( hecmesh, fstrsolid, fstrdynamic%t_delta )
392
393 enddo
394 !C-- end of time step loop
395 time_2 = hecmw_wtime()
396
397 if( hecmesh%my_rank == 0 ) then
398 write(ista,'(a,f10.2,a)') ' solve (sec) :', time_2 - time_1, 's'
399 endif
400
401 deallocate(coord)
402 call hecmw_mpc_mat_finalize(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
403 if (associated(hecmat0)) then
404 call hecmw_mat_finalize(hecmat0)
405 deallocate(hecmat0)
406 endif
407 end subroutine fstr_solve_dynamic_nlimplicit
408
411 subroutine fstr_solve_dynamic_nlimplicit_contactslag(cstep, hecMESH,hecMAT,fstrSOLID,fstrEIG &
412 ,fstrDYNAMIC,fstrRESULT,fstrPARAM &
413 ,fstrCPL,fstrMAT,restrt_step_num,infoCTChange &
414 ,conMAT )
415
416 use mcontact
420
421 implicit none
422 !C
423 !C-- global variable
424 !C
425 integer, intent(in) :: cstep
426 type(hecmwst_local_mesh) :: hecMESH
427 type(hecmwst_matrix) :: hecMAT
428 type(fstr_eigen) :: fstrEIG
429 type(fstr_solid) :: fstrSOLID
430 type(hecmwst_result_data) :: fstrRESULT
431 type(fstr_param) :: fstrPARAM
432 type(fstr_dynamic) :: fstrDYNAMIC
433 type(fstr_couple) :: fstrCPL !for COUPLE
434 type(fstrst_matrix_contact_lagrange) :: fstrMAT
435 type(fstr_info_contactchange) :: infoCTChange
436 type(hecmwst_matrix), optional :: conMAT
437
438 !C
439 !C-- local variable
440 !C
441
442 type(hecmwst_local_mesh), pointer :: hecMESHmpc
443 type(hecmwst_matrix), pointer :: hecMATmpc
444 integer(kind=kint) :: nnod, ndof, numnp, nn
445 integer(kind=kint) :: i, j, ids, ide, ims, ime, kk, idm, imm
446 integer(kind=kint) :: iter
447
448
449 real(kind=kreal) :: a1, a2, a3, b1, b2, b3, c1, c2
450 real(kind=kreal) :: bsize, res, res1, rf
451 real(kind=kreal) :: res0, relres
452 real :: time_1, time_2
453
454 integer(kind=kint) :: restrt_step_num
455
456 integer(kind=kint) :: ctAlgo
457 integer(kind=kint) :: max_iter_contact, count_step
458 integer(kind=kint) :: stepcnt
459 real(kind=kreal) :: maxdlag
460
461 logical :: is_mat_symmetric
462 integer(kind=kint) :: n_node_global
463 integer(kind=kint) :: contact_changed_global
464
465
466 integer(kind=kint) :: nndof,npdof
467 real(kind=kreal),allocatable :: tmp_conb(:)
468 integer :: istat
469 real(kind=kreal), allocatable :: coord(:)
470
471 call hecmw_mpc_mat_init(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
472
473 ! sum of n_node among all subdomains (to be used to calc res)
474 n_node_global = hecmesh%nn_internal
475 call hecmw_allreduce_i1(hecmesh,n_node_global,hecmw_sum)
476
477 ctalgo = fstrparam%contact_algo
478
479 if( hecmat%Iarray(99)==4 .and. .not.fstr_is_matrixstruct_symmetric(fstrsolid,hecmesh) ) then
480 write(*,*) ' This type of direct solver is not yet available in such case ! '
481 write(*,*) ' Please use intel MKL direct solver !'
482 call hecmw_abort(hecmw_comm_get_comm())
483 endif
484
485 hecmat%NDOF=hecmesh%n_dof
486
487 nnod=hecmesh%n_node
488 ndof=hecmat%NDOF
489 nn=ndof*ndof
490
491 allocate(coord(hecmesh%n_node*ndof))
492 if( associated( fstrsolid%contacts ) ) call initialize_contact_output_vectors(fstrsolid,hecmat)
493
494 !!
495 !!-- initial value
496 !!
497 time_1 = hecmw_wtime()
498
499 !C
500 !C-- check parameters
501 !C
502 if(dabs(fstrdynamic%beta) < 1.0e-20) then
503 if( hecmesh%my_rank == 0 ) then
504 write(imsg,*) 'stop due to Newmark-beta = 0'
505 endif
506 call hecmw_abort( hecmw_comm_get_comm())
507 endif
508
509
510 !C-- matrix [M]
511 !C-- lumped mass matrix
512 if(fstrdynamic%idx_mas == 1) then
513
514 call setmass(fstrsolid,hecmesh,hecmat,fstreig)
515
516 !C-- consistent mass matrix
517 else if(fstrdynamic%idx_mas == 2) then
518 if( hecmesh%my_rank .eq. 0 ) then
519 write(imsg,*) 'stop: consistent mass matrix is not yet available !'
520 endif
521 call hecmw_abort( hecmw_comm_get_comm())
522 endif
523 !C--
524 hecmat%Iarray(98) = 1 !Assmebly complete
525 hecmat%Iarray(97) = 1 !Need numerical factorization
526 !C
527 !C-- initialize variables
528 !C
529 if( restrt_step_num == 1 .and. fstrdynamic%VarInitialize .and. fstrdynamic%ray_m /= 0.0d0 ) &
530 call dynamic_init_varibles( hecmesh, hecmat, fstrsolid, fstreig, fstrdynamic, fstrparam )
531 !C
532 !C
533 !C-- time step loop
534 !C
535 a1 = .5d0/fstrdynamic%beta - 1.d0
536 a2 = 1.d0/(fstrdynamic%beta*fstrdynamic%t_delta)
537 a3 = 1.d0/(fstrdynamic%beta*fstrdynamic%t_delta*fstrdynamic%t_delta)
538 b1 = ( .5d0*fstrdynamic%ganma/fstrdynamic%beta - 1.d0 )*fstrdynamic%t_delta
539 b2 = fstrdynamic%ganma/fstrdynamic%beta - 1.d0
540 b3 = fstrdynamic%ganma/(fstrdynamic%beta*fstrdynamic%t_delta)
541 c1 = 1.d0 + fstrdynamic%ray_k*b3
542 c2 = a3 + fstrdynamic%ray_m*b3
543
544
545 !C-- output of initial state
546 if( restrt_step_num == 1 ) then
547 call fstr_dynamic_output(hecmesh, fstrsolid, fstrdynamic, fstrparam)
548 call dynamic_output_monit(hecmesh, fstrparam, fstrdynamic, fstreig, fstrsolid)
549 endif
550
551 fstrdynamic%VEC3(:) =0.d0
552 hecmat%X(:) =0.d0
553
555 call fstr_scan_contact_state(cstep, restrt_step_num, 0, fstrdynamic%t_delta, ctalgo, hecmesh, fstrsolid, infoctchange, hecmat%B)
556
557 if(paracontactflag.and.present(conmat)) then
558 call hecmw_mat_copy_profile( hecmat, conmat )
559 endif
560
561 if ( fstr_is_contact_active() ) then
562 ! ---- For Parallel Contact with Multi-Partition Domains
563 if(paracontactflag.and.present(conmat)) then
564 call fstr_mat_con_contact( cstep, hecmat, fstrsolid, fstrmat, infoctchange, conmat)
565 else
566 call fstr_mat_con_contact( cstep, hecmat, fstrsolid, fstrmat, infoctchange)
567 endif
568 elseif( hecmat%Iarray(99)==4 ) then
569 write(*,*) ' This type of direct solver is not yet available in such case ! '
570 write(*,*) ' Please change solver type to intel MKL direct solver !'
571 call hecmw_abort(hecmw_comm_get_comm())
572 endif
573 is_mat_symmetric = fstr_is_matrixstruct_symmetric(fstrsolid,hecmesh)
574 call solve_lineq_contact_init(hecmesh,hecmat,fstrmat,is_mat_symmetric)
575
576 !!
577 !! step = 1,2,....,fstrDYNAMIC%n_step
578 !!
579 do i= restrt_step_num, fstrdynamic%n_step
580
581 fstrdynamic%i_step = i
582 fstrdynamic%t_curr = fstrdynamic%t_delta * i
583
584 if(hecmesh%my_rank==0) then
585 write(ista,'('' time step='',i10,'' time='',1pe13.4e3)') i,fstrdynamic%t_curr
586 write(*,'(A)')'-------------------------------------------------'
587 write(*,'('' time step='',i10,'' time='',1pe13.4e3)') i,fstrdynamic%t_curr
588 endif
589
590 fstrsolid%dunode(:) =0.d0
591 ! call fstr_UpdateEPState( hecMESH, fstrSOLID )
592
593 do j = 1 ,ndof*nnod
594 fstrdynamic%VEC1(j) = a1*fstrdynamic%ACC(j,1) + a2*fstrdynamic%VEL(j,1)
595 fstrdynamic%VEC2(j) = b1*fstrdynamic%ACC(j,1) + b2*fstrdynamic%VEL(j,1)
596 enddo
597
598 max_iter_contact = 6 !1
599 count_step = 0
600 stepcnt = 0
601 loopforcontactanalysis: do while( .true. )
602 count_step = count_step + 1
603
604 ! ----- Inner Iteration
605 res0 = 0.d0
606 res1 = 0.d0
607 relres = 1.d0
608
609 do iter = 1, fstrsolid%step_ctrl(cstep)%max_iter
610 stepcnt=stepcnt+1
611 call fstr_stiffmatrix( hecmesh, hecmat, fstrsolid, fstrdynamic%t_curr, fstrdynamic%t_delta )
612 if( fstrdynamic%ray_k/=0.d0 .or. fstrdynamic%ray_m/=0.d0 ) then
613 do j = 1 ,ndof*nnod
614 hecmat%X(j) = fstrdynamic%VEC2(j) - b3*fstrsolid%dunode(j)
615 enddo
616 endif
617 if( fstrdynamic%ray_k/=0.d0 ) then
618 if( hecmesh%n_dof == 3 ) then
619 call hecmw_matvec (hecmesh, hecmat, hecmat%X, fstrdynamic%VEC3)
620 else if( hecmesh%n_dof == 2 ) then
621 call hecmw_matvec (hecmesh, hecmat, hecmat%X, fstrdynamic%VEC3)
622 else if( hecmesh%n_dof == 6 ) then
623 call matvec(fstrdynamic%VEC3, hecmat%X, hecmat, ndof, hecmat%D, hecmat%AU, hecmat%AL)
624 endif
625 endif
626 !C
627 !C-- mechanical boundary condition
628 call dynamic_mat_ass_load (hecmesh, hecmat, fstrsolid, fstrdynamic, fstrparam)
629 do j=1, hecmesh%n_node* hecmesh%n_dof
630 hecmat%B(j)=hecmat%B(j)- fstrsolid%QFORCE(j) + fstreig%mass(j)*( fstrdynamic%VEC1(j)-a3*fstrsolid%dunode(j) &
631 + fstrdynamic%ray_m* hecmat%X(j) ) + fstrdynamic%ray_k*fstrdynamic%VEC3(j)
632 enddo
633 do j = 1 ,nn*hecmat%NP
634 hecmat%D(j) = c1* hecmat%D(j)
635 enddo
636 do j = 1 ,nn*hecmat%NPU
637 hecmat%AU(j) = c1* hecmat%AU(j)
638 enddo
639 do j = 1 ,nn*hecmat%NPL
640 hecmat%AL(j) = c1*hecmat%AL(j)
641 enddo
642 do j=1,nnod
643 do kk=1,ndof
644 idm = nn*(j-1)+1 + (ndof+1)*(kk-1)
645 imm = ndof*(j-1) + kk
646 hecmat%D(idm) = hecmat%D(idm) + c2*fstreig%mass(imm)
647 enddo
648 enddo
649
650
651 if(paracontactflag.and.present(conmat)) then
652 call hecmw_mat_clear( conmat )
653 call hecmw_mat_clear_b( conmat )
654 conmat%X = 0.0d0
655 endif
656 if( fstr_is_contact_active() ) then
657 ! ---- For Parallel Contact with Multi-Partition Domains
658 if(paracontactflag.and.present(conmat)) then
659 call fstr_update_ndforce_contact(cstep,hecmesh,hecmat,fstrmat,fstrsolid,conmat)
660 call fstr_addcontactstiffness(cstep,iter,conmat,fstrmat,fstrsolid)
661 else
662 call fstr_update_ndforce_contact(cstep,hecmesh,hecmat,fstrmat,fstrsolid)
663 call fstr_addcontactstiffness(cstep,iter,hecmat,fstrmat,fstrsolid)
664 endif
665 endif
666 !
667 !C ********************************************************************************
668 !C for couple analysis
669 if( fstrparam%fg_couple == 1) then
670 if( fstrdynamic%i_step > 1 .or. &
671 (fstrdynamic%i_step==1 .and. fstrparam%fg_couple_first==1 )) then
672 call fstr_rcap_get( fstrcpl )
673 call dynamic_mat_ass_couple( hecmesh, hecmat, fstrsolid, fstrcpl )
674 endif
675 endif
676 !C ********************************************************************************
677
678 !C-- geometrical boundary condition
679 call hecmw_mpc_mat_ass(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
680 call hecmw_mpc_trans_rhs(hecmesh, hecmat, hecmatmpc)
681 if(paracontactflag.and.present(conmat)) then
682 call dynamic_mat_ass_bc (hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, stepcnt, conmat=conmat)
683 call dynamic_mat_ass_bc_vl(hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, stepcnt, conmat=conmat)
684 call dynamic_mat_ass_bc_ac(hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, stepcnt, conmat=conmat)
685 else
686 call dynamic_mat_ass_bc (hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, stepcnt)
687 call dynamic_mat_ass_bc_vl(hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, stepcnt)
688 call dynamic_mat_ass_bc_ac(hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, stepcnt)
689 endif
690
691 ! ----- check convergence
692 ! ---- For Parallel Contact with Multi-Partition Domains
693 if(paracontactflag.and.present(conmat)) then
694 res = fstr_get_norm_para_contact(hecmatmpc,fstrmat,conmat,hecmesh)
695 else
696 call hecmw_innerproduct_r(hecmesh,ndof,hecmatmpc%B,hecmatmpc%B,res)
697 endif
698 res = sqrt(res)/n_node_global
699
700 if( iter==1 ) res0=res
701 if( res0==0.d0 ) then
702 res0 =1.d0
703 else
704 relres = dabs(res1-res)/res0
705 endif
706
707 ! ----- check convergence
708 if( .not.fstr_is_contact_active() ) then
709 maxdlag = 0.0d0
710 elseif( maxdlag == 0.0d0) then
711 maxdlag = 1.0d0
712 endif
713 call hecmw_allreduce_r1(hecmesh, maxdlag, hecmw_max)
714
715 if( hecmesh%my_rank==0 ) then
716 ! if( mod(i,max(int(fstrDYNAMIC%nout/10),1)) == 0 ) &
717 write(*,'(a,i3,a,2e15.7)') ' - Residual(',iter,') =',res,relres
718 write(*,'(a,1e15.7)') ' - MaxDLag =',maxdlag
719 write(ista,'(''iter='',I5,''res/res0='',2E15.7)')iter,res,relres
720 write(ista,'(a,1e15.7)') ' - MaxDLag =',maxdlag
721 endif
722
723 ! ----- check convergence
724 ! if( .not.fstr_is_contact_active() ) maxDLag= 0.0d0
725 ! call hecmw_allreduce_R1(hecMESH, maxDlag, HECMW_MAX)
726 ! if( res<fstrSOLID%step_ctrl(cstep)%converg .and. maxDLag<1.0d-5 .and. iter>1 ) exit
727 if( (res<fstrsolid%step_ctrl(cstep)%converg .or. &
728 relres<fstrsolid%step_ctrl(cstep)%converg) .and. maxdlag<1.0d-4 ) exit
729 res1 = res
730 rf=1.0d0
731 if( iter>1 .and. res>res1 )rf=0.5d0*rf
732 res1=res
733
734 ! ---- For Parallel Contact with Multi-Partition Domains
735 hecmatmpc%X = 0.0d0
736 call fstr_set_current_config_to_mesh(hecmeshmpc,fstrsolid,coord)
737 if(paracontactflag.and.present(conmat)) then
738 call solve_lineq_contact(hecmeshmpc,hecmatmpc,fstrmat,istat,1.0d0,conmat)
739 else
740 call solve_lineq_contact(hecmeshmpc,hecmatmpc,fstrmat,istat,rf)
741 endif
742 call fstr_recover_initial_config_to_mesh(hecmeshmpc,fstrsolid,coord)
743 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
744
745 ! ----- update external nodal displacement increments
746 call hecmw_update_3_r (hecmesh, hecmat%X, hecmat%NP)
747
748 ! ----- update the strain, stress, and internal force
749 do j=1,hecmesh%n_node*ndof
750 fstrsolid%dunode(j) = fstrsolid%dunode(j)+hecmat%X(j)
751 enddo
752 call fstr_updatenewton( hecmesh, hecmat, fstrsolid, fstrdynamic%t_curr, &
753 & fstrdynamic%t_delta,iter, fstrdynamic%strainEnergy )
754
755
756 ! ----- update the Lagrange multipliers
757 if( fstr_is_contact_active() ) then
758 maxdlag = 0.0d0
759 do j=1,fstrmat%num_lagrange
760 fstrmat%lagrange(j) = fstrmat%lagrange(j) + hecmat%X(hecmesh%n_node*ndof+j)
761 if(dabs(hecmat%X(hecmesh%n_node*ndof+j))>maxdlag) maxdlag=dabs(hecmat%X(hecmesh%n_node*ndof+j))
762 ! write(*,*)'Lagrange:', j,fstrMAT%lagrange(j),hecMAT%X(hecMESH%n_node*ndof+j)
763 enddo
764 endif
765 enddo
766
767 ! ----- not convergence
768 if( iter>fstrsolid%step_ctrl(cstep)%max_iter ) then
769 if( hecmesh%my_rank==0) then
770 write(ilog,*) '### Fail to Converge : at step=', i
771 write(ista,*) '### Fail to Converge : at step=', i
772 write( *,*) ' ### Fail to Converge : at step=', i
773 endif
774 stop
775 endif
776
777 call fstr_scan_contact_state(cstep, i, count_step, fstrdynamic%t_delta, ctalgo, hecmesh, fstrsolid, infoctchange, hecmat%B)
778
779 if( hecmat%Iarray(99)==4 .and. .not. fstr_is_contact_active() ) then
780 write(*,*) ' This type of direct solver is not yet available in such case ! '
781 write(*,*) ' Please use intel MKL direct solver !'
782 call hecmw_abort(hecmw_comm_get_comm())
783 endif
784
785 is_mat_symmetric = fstr_is_matrixstruct_symmetric(fstrsolid,hecmesh)
786 contact_changed_global=0
787 if( fstr_is_contact_conv(ctalgo,infoctchange,hecmesh) ) then
788 exit loopforcontactanalysis
789 elseif( fstr_is_matrixstructure_changed(infoctchange) ) then
790 ! ---- For Parallel Contact with Multi-Partition Domains
791 if(paracontactflag.and.present(conmat)) then
792 call fstr_mat_con_contact( cstep, hecmat, fstrsolid, fstrmat, infoctchange, conmat)
793 else
794 call fstr_mat_con_contact( cstep, hecmat, fstrsolid, fstrmat, infoctchange)
795 endif
796 contact_changed_global=1
797 endif
798 call hecmw_allreduce_i1(hecmesh,contact_changed_global,hecmw_max)
799 if (contact_changed_global > 0) then
800 call hecmw_mat_clear_b( hecmat )
801 if(paracontactflag.and.present(conmat)) call hecmw_mat_clear_b( conmat )
802 call solve_lineq_contact_init(hecmesh,hecmat,fstrmat,is_mat_symmetric)
803 endif
804
805 if( count_step > max_iter_contact ) exit loopforcontactanalysis
806
807
808 enddo loopforcontactanalysis
809 !
810 !C-- new displacement, velocity and accelaration
811 !C
812 fstrdynamic%kineticEnergy = 0.0d0
813 do j = 1 ,ndof*nnod
814 fstrdynamic%ACC (j,2) = -a1*fstrdynamic%ACC(j,1) - a2*fstrdynamic%VEL(j,1) + &
815 a3*fstrsolid%dunode(j)
816 fstrdynamic%VEL (j,2) = -b1*fstrdynamic%ACC(j,1) - b2*fstrdynamic%VEL(j,1) + &
817 b3*fstrsolid%dunode(j)
818 fstrdynamic%ACC (j,1) = fstrdynamic%ACC (j,2)
819 fstrdynamic%VEL (j,1) = fstrdynamic%VEL (j,2)
820
821 fstrsolid%unode(j) = fstrsolid%unode(j)+fstrsolid%dunode(j)
822 fstrdynamic%DISP(j,2) = fstrsolid%unode(j)
823
824 fstrdynamic%kineticEnergy = fstrdynamic%kineticEnergy + &
825 0.5d0*fstreig%mass(j)*fstrdynamic%VEL(j,2)*fstrdynamic%VEL(j,2)
826 enddo
827
828 !--- Restart info
829 if( fstrdynamic%restart_nout > 0 .and. &
830 (mod(i,fstrdynamic%restart_nout).eq.0 .or. i.eq.fstrdynamic%n_step) ) then
831 call fstr_write_restart_dyna_nl(i,hecmesh,fstrsolid,fstrdynamic,fstrparam,&
832 infoctchange%contactNode_current)
833 endif
834
835 !C-- output new displacement, velocity and accelaration
836 call fstr_dynamic_output(hecmesh, fstrsolid, fstrdynamic, fstrparam)
837
838 !C-- output result of monitoring node
839 call dynamic_output_monit(hecmesh, fstrparam, fstrdynamic, fstreig, fstrsolid)
840
841 call fstr_updatestate( hecmesh, fstrsolid, fstrdynamic%t_delta )
842
843 enddo
844 !C
845 !C-- end of time step loop
846
847 time_2 = hecmw_wtime()
848 if( hecmesh%my_rank == 0 ) then
849 write(ista,'(a,f10.2,a)') ' solve (sec) :', time_2 - time_1, 's'
850 endif
851
852 deallocate(coord)
853 call hecmw_mpc_mat_finalize(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
855
This module contains subroutines for nonlinear implicit dynamic analysis.
subroutine fstr_solve_dynamic_nlimplicit(cstep, hecmesh, hecmat, fstrsolid, fstreig, fstrdynamic, fstrresult, fstrparam, fstrcpl, restrt_step_num)
This subroutine provides function of nonlinear implicit dynamic analysis using the Newmark method.
subroutine fstr_solve_dynamic_nlimplicit_contactslag(cstep, hecmesh, hecmat, fstrsolid, fstreig, fstrdynamic, fstrresult, fstrparam, fstrcpl, fstrmat, restrt_step_num, infoctchange, conmat)
This subroutine provides function of nonlinear implicit dynamic analysis using the Newmark method....
This module provides functions of reconstructing.
subroutine fstr_save_originalmatrixstructure(hecmat)
This subroutine saves original matrix structure constructed originally by hecMW_matrix.
logical function fstr_is_matrixstruct_symmetric(fstrsolid, hecmesh)
this function judges whether sitiffness matrix is symmetric or not
subroutine fstr_mat_con_contact(cstep, hecmat, fstrsolid, fstrmat, infoctchange, conmat)
this subroutine reconstructs node-based (stiffness) matrix structure \corresponding to contact state
This module provides functions: 1) obtain contact stiffness matrix of each contact pair and assemble ...
subroutine, public fstr_update_ndforce_contact(cstep, hecmesh, hecmat, fstrmat, fstrsolid, conmat)
This subroutine obtains contact nodal force vector of each contact pair and assembles it into right-h...
subroutine, public fstr_addcontactstiffness(cstep, iter, hecmat, fstrmat, fstrsolid)
This subroutine obtains contact stiffness matrix of each contact pair and assembles it into global st...
This module provides functions to initialize variables when initial velocity or acceleration boundary...
subroutine dynamic_init_varibles(hecmesh, hecmat, fstrsolid, fstreig, fstrdynamic, fstrparam)
This module contains functions to set acceleration boundary condition in dynamic analysis.
subroutine dynamic_mat_ass_bc_ac(hecmesh, hecmat, fstrsolid, fstrdynamic, fstrparam, fstrmat, iter, conmat)
This subrouitne set acceleration boundary condition in dynamic analysis.
This module contains functions to set velocity boundary condition in dynamic analysis.
subroutine dynamic_mat_ass_bc_vl(hecmesh, hecmat, fstrsolid, fstrdynamic, fstrparam, fstrmat, iter, conmat)
This subrouitne set velocity boundary condition in dynamic analysis.
This module contains functions to set displacement boundary condition in dynamic analysis.
subroutine dynamic_mat_ass_bc(hecmesh, hecmat, fstrsolid, fstrdynamic, fstrparam, fstrmat, iter, conmat)
This subroutine setup disp bundary 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 matvec(y, x, hecmat, ndof, d, au, al)
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 function to calcualte residual of nodal force.
real(kind=kreal) function, public fstr_get_norm_para_contact(hecmat, fstrmat, conmat, hecmesh)
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 tangent stiffness matrix.
subroutine, public fstr_stiffmatrix(hecmesh, hecmat, fstrsolid, time, tincr)
接線剛性マトリックスを作成するサブルーチン
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
subroutine fstr_set_current_config_to_mesh(hecmesh, fstrsolid, coord)
Definition: m_fstr.f90:1057
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
integer(kind=kint), parameter ilog
FILE HANDLER.
Definition: m_fstr.f90:91
integer(kind=kint), parameter ista
Definition: m_fstr.f90:92
subroutine fstr_recover_initial_config_to_mesh(hecmesh, fstrsolid, coord)
Definition: m_fstr.f90:1070
logical paracontactflag
PARALLEL CONTACT FLAG.
Definition: m_fstr.f90:84
This module provides functions to solve sparse system of \linear equitions in the case of contact ana...
subroutine, public solve_lineq_contact_init(hecmesh, hecmat, fstrmat, is_sym)
This subroutine.
subroutine, public solve_lineq_contact(hecmesh, hecmat, fstrmat, istat, rf, conmat)
This subroutine.
This module provides functions to calcualte contact stiff matrix.
Definition: fstr_contact.f90:6
logical function fstr_is_matrixstructure_changed(infoctchange)
subroutine fstr_scan_contact_state(cstep, sub_step, cont_step, dt, ctalgo, hecmesh, fstrsolid, infoctchange, b)
Scanning contact state.
logical function fstr_is_contact_conv(ctalgo, infoctchange, hecmesh)
logical function fstr_is_contact_active()
subroutine initialize_contact_output_vectors(fstrsolid, hecmat)
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