30 fstrDYNAMIC,fstrRESULT,fstrPARAM,fstrCPL, restrt_step_num )
33 integer,
intent(in) :: cstep
34 type(hecmwst_local_mesh) :: hecMESH
35 type(hecmwst_matrix) :: hecMAT
38 type(hecmwst_result_data) :: fstrRESULT
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
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(:)
67 call hecmw_mpc_mat_init(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
71 n_node_global = hecmesh%nn_internal
72 call hecmw_allreduce_i1(hecmesh,n_node_global,hecmw_sum)
74 hecmat%NDOF=hecmesh%n_dof
79 allocate(coord(hecmesh%n_node*ndof))
82 time_1 = hecmw_wtime()
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'
89 call hecmw_abort( hecmw_comm_get_comm())
93 if(fstrdynamic%idx_mas == 1)
then
94 call setmass(fstrsolid,hecmesh,hecmat,fstreig)
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 !'
101 call hecmw_abort( hecmw_comm_get_comm())
104 hecmat%Iarray(98) = 1
105 hecmat%Iarray(97) = 1
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
118 if( restrt_step_num == 1 )
then
123 fstrdynamic%VEC3(:) =0.0d0
127 do i= restrt_step_num, fstrdynamic%n_step
128 if(ndof == 4 .and. hecmesh%my_rank==0)
write(*,
'(a,i5)')
"iter: ",i
130 fstrdynamic%i_step = i
131 fstrdynamic%t_curr = fstrdynamic%t_delta * i
133 if(hecmesh%my_rank==0)
then
135 write(
ista,
'('' time step='',i10,'' time='',1pe13.4e3)') i,fstrdynamic%t_curr
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)
146 fstrsolid%dunode(:) =0.d0
149 if( fstrparam%fg_couple == 1)
then
150 if( fstrparam%fg_couple_type==1 .or. &
151 fstrparam%fg_couple_type==3 .or. &
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 )
159 if (.not.
associated(hecmat0))
then
160 call fstr_stiffmatrix( hecmesh, hecmat, fstrsolid, fstrdynamic%t_curr, fstrdynamic%t_delta )
162 call hecmw_mat_init(hecmat0)
163 call hecmw_mat_copy_profile(hecmat, hecmat0)
164 call hecmw_mat_copy_val(hecmat, hecmat0)
166 call hecmw_mat_copy_val(hecmat0, hecmat)
170 if( fstrdynamic%ray_k/=0.d0 .or. fstrdynamic%ray_m/=0.d0 )
then
172 hecmat%X(j) = fstrdynamic%VEC2(j) - b3*fstrsolid%dunode(j)
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)
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)
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
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 )
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
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 )
219 do j = 1 ,nn*hecmat%NP
220 hecmat%D(j) = c1* hecmat%D(j)
222 do j = 1 ,nn*hecmat%NPU
223 hecmat%AU(j) = c1* hecmat%AU(j)
225 do j = 1 ,nn*hecmat%NPL
226 hecmat%AL(j) = c1*hecmat%AL(j)
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)
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)
245 call hecmw_innerproduct_r(hecmesh, ndof, hecmatmpc%B, hecmatmpc%B, bsize)
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
261 if( iexit .ne. 1 )
then
262 if( fstrparam%nlgeom )
then
264 hecmatmpc%Iarray(97) = 2
266 hecmatmpc%Iarray(97) = 1
270 call solve_lineq(hecmeshmpc,hecmatmpc)
271 if( fstrparam%nlgeom )
then
275 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
277 do j=1,hecmesh%n_node*ndof
278 fstrsolid%dunode(j) = fstrsolid%dunode(j)+hecmat%X(j)
282 & fstrdynamic%t_delta, iter, fstrdynamic%strainEnergy )
284 if(.not. fstrparam%nlgeom)
exit
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
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
305 kkk1 = fstrcpl%coupled_node(j)*3
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 )
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)
325 kkk1 = fstrcpl%coupled_node(j)*2
327 fstrcpl%disp (kkk0-1) = fstrsolid%unode(kkk1-1) + fstrsolid%dunode(kkk1-1)
328 fstrcpl%disp (kkk0 ) = fstrsolid%unode(kkk1 ) + fstrsolid%dunode(kkk1 )
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)
343 select case ( fstrparam%fg_couple_type )
348 if( revocap_flag==0 ) cycle
351 if( revocap_flag==0 )
then
364 fstrdynamic%kineticEnergy = 0.0d0
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)
373 fstrsolid%unode(j) = fstrsolid%unode(j)+fstrsolid%dunode(j)
374 fstrdynamic%DISP(j,2) = fstrsolid%unode(j)
376 fstrdynamic%kineticEnergy = fstrdynamic%kineticEnergy + &
377 0.5d0*fstreig%mass(j)*fstrdynamic%VEL(j,2)*fstrdynamic%VEL(j,2)
381 if( fstrdynamic%restart_nout > 0 .and. &
382 (mod(i,fstrdynamic%restart_nout).eq.0 .or. i.eq.fstrdynamic%n_step) )
then
395 time_2 = hecmw_wtime()
397 if( hecmesh%my_rank == 0 )
then
398 write(
ista,
'(a,f10.2,a)')
' solve (sec) :', time_2 - time_1,
's'
402 call hecmw_mpc_mat_finalize(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
403 if (
associated(hecmat0))
then
404 call hecmw_mat_finalize(hecmat0)
412 ,fstrDYNAMIC,fstrRESULT,fstrPARAM &
413 ,fstrCPL,fstrMAT,restrt_step_num,infoCTChange &
425 integer,
intent(in) :: cstep
426 type(hecmwst_local_mesh) :: hecMESH
427 type(hecmwst_matrix) :: hecMAT
430 type(hecmwst_result_data) :: fstrRESULT
435 type(fstr_info_contactchange) :: infoCTChange
436 type(hecmwst_matrix),
optional :: conMAT
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
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
454 integer(kind=kint) :: restrt_step_num
456 integer(kind=kint) :: ctAlgo
457 integer(kind=kint) :: max_iter_contact, count_step
458 integer(kind=kint) :: stepcnt
459 real(kind=kreal) :: maxdlag
461 logical :: is_mat_symmetric
462 integer(kind=kint) :: n_node_global
463 integer(kind=kint) :: contact_changed_global
466 integer(kind=kint) :: nndof,npdof
467 real(kind=kreal),
allocatable :: tmp_conb(:)
469 real(kind=kreal),
allocatable :: coord(:)
471 call hecmw_mpc_mat_init(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
474 n_node_global = hecmesh%nn_internal
475 call hecmw_allreduce_i1(hecmesh,n_node_global,hecmw_sum)
477 ctalgo = fstrparam%contact_algo
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())
485 hecmat%NDOF=hecmesh%n_dof
491 allocate(coord(hecmesh%n_node*ndof))
497 time_1 = hecmw_wtime()
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'
506 call hecmw_abort( hecmw_comm_get_comm())
512 if(fstrdynamic%idx_mas == 1)
then
514 call setmass(fstrsolid,hecmesh,hecmat,fstreig)
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 !'
521 call hecmw_abort( hecmw_comm_get_comm())
524 hecmat%Iarray(98) = 1
525 hecmat%Iarray(97) = 1
529 if( restrt_step_num == 1 .and. fstrdynamic%VarInitialize .and. fstrdynamic%ray_m /= 0.0d0 ) &
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
546 if( restrt_step_num == 1 )
then
551 fstrdynamic%VEC3(:) =0.d0
555 call fstr_scan_contact_state(cstep, restrt_step_num, 0, fstrdynamic%t_delta, ctalgo, hecmesh, fstrsolid, infoctchange, hecmat%B)
558 call hecmw_mat_copy_profile( hecmat, conmat )
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())
579 do i= restrt_step_num, fstrdynamic%n_step
581 fstrdynamic%i_step = i
582 fstrdynamic%t_curr = fstrdynamic%t_delta * i
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
590 fstrsolid%dunode(:) =0.d0
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)
601 loopforcontactanalysis:
do while( .true. )
602 count_step = count_step + 1
609 do iter = 1, fstrsolid%step_ctrl(cstep)%max_iter
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
614 hecmat%X(j) = fstrdynamic%VEC2(j) - b3*fstrsolid%dunode(j)
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)
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)
633 do j = 1 ,nn*hecmat%NP
634 hecmat%D(j) = c1* hecmat%D(j)
636 do j = 1 ,nn*hecmat%NPU
637 hecmat%AU(j) = c1* hecmat%AU(j)
639 do j = 1 ,nn*hecmat%NPL
640 hecmat%AL(j) = c1*hecmat%AL(j)
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)
652 call hecmw_mat_clear( conmat )
653 call hecmw_mat_clear_b( conmat )
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
679 call hecmw_mpc_mat_ass(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
680 call hecmw_mpc_trans_rhs(hecmesh, hecmat, hecmatmpc)
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)
686 call dynamic_mat_ass_bc (hecmesh, hecmatmpc, fstrsolid, fstrdynamic, fstrparam, fstrmat, stepcnt)
696 call hecmw_innerproduct_r(hecmesh,ndof,hecmatmpc%B,hecmatmpc%B,res)
698 res = sqrt(res)/n_node_global
700 if( iter==1 ) res0=res
701 if( res0==0.d0 )
then
704 relres = dabs(res1-res)/res0
710 elseif( maxdlag == 0.0d0)
then
713 call hecmw_allreduce_r1(hecmesh, maxdlag, hecmw_max)
715 if( hecmesh%my_rank==0 )
then
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
727 if( (res<fstrsolid%step_ctrl(cstep)%converg .or. &
728 relres<fstrsolid%step_ctrl(cstep)%converg) .and. maxdlag<1.0d-4 )
exit
731 if( iter>1 .and. res>res1 )rf=0.5d0*rf
743 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
746 call hecmw_update_3_r (hecmesh, hecmat%X, hecmat%NP)
749 do j=1,hecmesh%n_node*ndof
750 fstrsolid%dunode(j) = fstrsolid%dunode(j)+hecmat%X(j)
753 & fstrdynamic%t_delta,iter, fstrdynamic%strainEnergy )
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))
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
777 call fstr_scan_contact_state(cstep, i, count_step, fstrdynamic%t_delta, ctalgo, hecmesh, fstrsolid, infoctchange, hecmat%B)
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())
786 contact_changed_global=0
788 exit loopforcontactanalysis
796 contact_changed_global=1
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 )
805 if( count_step > max_iter_contact )
exit loopforcontactanalysis
808 enddo loopforcontactanalysis
812 fstrdynamic%kineticEnergy = 0.0d0
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)
821 fstrsolid%unode(j) = fstrsolid%unode(j)+fstrsolid%dunode(j)
822 fstrdynamic%DISP(j,2) = fstrsolid%unode(j)
824 fstrdynamic%kineticEnergy = fstrdynamic%kineticEnergy + &
825 0.5d0*fstreig%mass(j)*fstrdynamic%VEL(j,2)*fstrdynamic%VEL(j,2)
829 if( fstrdynamic%restart_nout > 0 .and. &
830 (mod(i,fstrdynamic%restart_nout).eq.0 .or. i.eq.fstrdynamic%n_step) )
then
832 infoctchange%contactNode_current)
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'
853 call hecmw_mpc_mat_finalize(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
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 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.
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.
subroutine fstr_updatenewton(hecmesh, hecmat, fstrsolid, time, tincr, iter, strainenergy)
変位/応力・ひずみ/内力のアップデート
subroutine fstr_updatestate(hecmesh, fstrsolid, tincr)
Update elastiplastic status.
This module defined coomon data and basic structures for analysis.
subroutine fstr_set_current_config_to_mesh(hecmesh, fstrsolid, coord)
integer(kind=kint), parameter imsg
integer(kind=kint), parameter ilog
FILE HANDLER.
integer(kind=kint), parameter ista
subroutine fstr_recover_initial_config_to_mesh(hecmesh, fstrsolid, coord)
logical paracontactflag
PARALLEL CONTACT FLAG.
Data for coupling analysis.
Data for DYNAMIC ANSLYSIS (fstrDYNAMIC)
Package of data used by Lanczos eigenvalue solver.
FSTR INNER CONTROL PARAMETERS (fstrPARAM)