29 subroutine fstr_newton( cstep, hecMESH, hecMAT, fstrSOLID, fstrPARAM, &
30 restrt_step_num, sub_step, ctime, dtime )
32 integer,
intent(in) :: cstep
33 type (hecmwST_local_mesh) :: hecMESH
34 type (hecmwST_matrix) :: hecMAT
35 type (fstr_solid) :: fstrSOLID
36 integer,
intent(in) :: sub_step
37 real(kind=kreal),
intent(in) :: ctime
38 real(kind=kreal),
intent(in) :: dtime
39 type (fstr_param) :: fstrPARAM
40 type (fstrST_matrix_contact_lagrange) :: fstrMAT
42 type (hecmwST_local_mesh),
pointer :: hecMESHmpc
43 type (hecmwST_matrix),
pointer :: hecMATmpc
44 integer(kind=kint) :: ndof
45 integer(kind=kint) :: i, iter
46 integer(kind=kint) :: stepcnt
47 integer(kind=kint) :: restrt_step_num
48 real(kind=kreal) :: tt0, tt, res, qnrm, rres, tincr, xnrm, dunrm, rxnrm
49 real(kind=kreal),
allocatable :: coord(:), p(:)
50 logical :: isLinear = .false.
52 call hecmw_mpc_mat_init(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
54 if(.not.
fstrpr%nlgeom)
then
58 hecmat%NDOF = hecmesh%n_dof
61 allocate(p(hecmesh%n_node*ndof))
62 allocate(coord(hecmesh%n_node*ndof))
65 if( fstrsolid%step_ctrl(cstep)%solution == stepstatic ) tincr = 0.d0
69 fstrsolid%dunode(:) = 0.0d0
70 fstrsolid%NRstat_i(:) = 0
72 call fstr_ass_load(cstep, ctime+dtime, hecmesh, hecmat, fstrsolid, fstrparam)
75 do iter=1,fstrsolid%step_ctrl(cstep)%max_iter
82 call hecmw_mpc_mat_ass(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
83 call hecmw_mpc_trans_rhs(hecmesh, hecmat, hecmatmpc)
84 call fstr_addbc(cstep, hecmesh, hecmatmpc, fstrsolid, fstrparam, fstrmat, stepcnt)
87 if( sub_step == restrt_step_num .and. iter == 1 ) hecmatmpc%Iarray(98) = 1
89 hecmatmpc%Iarray(97) = 2
91 hecmatmpc%Iarray(97) = 1
95 call solve_lineq(hecmeshmpc,hecmatmpc)
97 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
102 do i = 1, hecmesh%n_node*ndof
103 fstrsolid%dunode(i) = fstrsolid%dunode(i) + hecmat%X(i)
110 if( fstrsolid%DLOAD_follow /= 0 .or. fstrsolid%CLOAD_ngrp_rot /= 0 ) &
111 &
call fstr_ass_load(cstep, ctime+dtime, hecmesh, hecmat, fstrsolid, fstrparam )
118 call hecmw_innerproduct_r(hecmesh, ndof, hecmat%B, hecmat%B, res)
120 call hecmw_innerproduct_r(hecmesh, ndof, hecmat%X, hecmat%X, xnrm)
122 call hecmw_innerproduct_r(hecmesh, ndof, fstrsolid%QFORCE, fstrsolid%QFORCE, qnrm)
124 if (qnrm < 1.0d-8) qnrm = 1.0d0
128 call hecmw_innerproduct_r(hecmesh, ndof, fstrsolid%dunode, fstrsolid%dunode, dunrm)
133 if( hecmesh%my_rank == 0 )
then
134 if (qnrm == 1.0d0)
then
135 write(*,
"(a,i8,a,1pe11.4,a,1pe11.4)")
" iter:", iter,
", residual(abs):", rres,
", disp.corr.:", rxnrm
137 write(*,
"(a,i8,a,1pe11.4,a,1pe11.4)")
" iter:", iter,
", residual:", rres,
", disp.corr.:", rxnrm
140 if( hecmw_mat_get_flag_diverged(hecmat) ==
kno )
then
141 if( rres < fstrsolid%step_ctrl(cstep)%converg )
exit
142 if( rxnrm < fstrsolid%step_ctrl(cstep)%converg )
exit
146 if( iter == fstrsolid%step_ctrl(cstep)%max_iter .or. rres > fstrsolid%step_ctrl(cstep)%maxres .or. rres /= rres )
then
147 if( hecmesh%my_rank == 0)
then
148 write(
ilog,
'(a,i5,a,i5)')
'### Fail to Converge : at total_step=', cstep,
' sub_step=', sub_step
149 write( *,
'(a,i5,a,i5)')
' ### Fail to Converge : at total_step=', cstep,
' sub_step=', sub_step
151 fstrsolid%NRstat_i(knstmaxit) = max(fstrsolid%NRstat_i(knstmaxit),iter)
152 fstrsolid%NRstat_i(knstsumit) = fstrsolid%NRstat_i(knstsumit) + iter
153 fstrsolid%CutBack_stat = fstrsolid%CutBack_stat + 1
154 if( iter == fstrsolid%step_ctrl(cstep)%max_iter ) fstrsolid%NRstat_i(knstdresn) = 1
155 if( rres > fstrsolid%step_ctrl(cstep)%maxres .or. rres /= rres ) fstrsolid%NRstat_i(knstdresn) = 2
161 fstrsolid%NRstat_i(knstmaxit) = max(fstrsolid%NRstat_i(knstmaxit),iter)
162 fstrsolid%NRstat_i(knstsumit) = fstrsolid%NRstat_i(knstsumit) + iter
166 do i=1,hecmesh%n_node*ndof
167 fstrsolid%unode(i) = fstrsolid%unode(i) + fstrsolid%dunode(i)
172 fstrsolid%CutBack_stat = 0
175 call hecmw_mpc_mat_finalize(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
183 restart_step_num, restart_substep_num, sub_step, ctime, dtime, infoCTChange )
186 integer,
intent(in) :: cstep
187 type (hecmwST_local_mesh) :: hecMESH
188 type (hecmwST_matrix) :: hecMAT
189 type (fstr_solid) :: fstrSOLID
190 integer,
intent(in) :: sub_step
191 real(kind=kreal),
intent(in) :: ctime
192 real(kind=kreal),
intent(in) :: dtime
193 type (fstr_param) :: fstrPARAM
194 type (fstr_info_contactChange) :: infoCTChange
195 type (fstrST_matrix_contact_lagrange) :: fstrMAT
197 type (hecmwST_local_mesh),
pointer :: hecMESHmpc
198 type (hecmwST_matrix),
pointer :: hecMATmpc
199 integer(kind=kint) :: ndof
200 integer(kind=kint) :: ctAlgo
201 integer(kind=kint) :: i, iter
202 integer(kind=kint) :: al_step, n_al_step, stepcnt
203 real(kind=kreal) :: tt0, tt, res, res0, res1, maxv, relres, tincr
204 integer(kind=kint) :: restart_step_num, restart_substep_num
205 logical :: convg, ctchange
206 integer(kind=kint) :: n_node_global
207 real(kind=kreal),
allocatable :: coord(:)
209 call hecmw_mpc_mat_init(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
212 n_node_global = hecmesh%nn_internal
213 call hecmw_allreduce_i1(hecmesh,n_node_global,hecmw_sum)
215 ctalgo = fstrparam%contact_algo
217 hecmat%NDOF = hecmesh%n_dof
220 fstrsolid%NRstat_i(:) = 0
222 allocate(coord(hecmesh%n_node*ndof))
225 if( fstrsolid%step_ctrl(cstep)%solution == stepstatic ) tincr = 0.0d0
227 if( cstep == 1 .and. sub_step == restart_substep_num )
then
228 if(hecmesh%my_rank==0)
write(*,*)
"---Scanning initial contact state---"
230 if(hecmesh%my_rank==0)
write(*,*)
237 call fstr_ass_load(cstep, ctime+dtime, hecmesh, hecmat, fstrsolid, fstrparam)
240 n_al_step = fstrsolid%step_ctrl(cstep)%max_contiter
243 do al_step = 1, n_al_step
245 if( hecmesh%my_rank == 0)
then
246 if( n_al_step > 1 )
then
247 print *,
"Augmentation step:", al_step
248 write(
imsg, *)
"Augmentation step:", al_step
252 fstrsolid%dunode(:) = 0.0d0
259 do iter = 1,fstrsolid%step_ctrl(cstep)%max_iter
267 maxv = hecmw_mat_diag_max( hecmat, hecmesh )
275 call hecmw_mpc_mat_ass(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
276 call hecmw_mpc_trans_rhs(hecmesh, hecmat, hecmatmpc)
277 call fstr_addbc(cstep, hecmesh,hecmatmpc,fstrsolid,fstrparam,fstrmat,stepcnt)
280 if( sub_step == restart_step_num .and. iter == 1 ) hecmatmpc%Iarray(98) = 1
282 hecmatmpc%Iarray(97) = 2
284 hecmatmpc%Iarray(97) = 1
288 call solve_lineq(hecmeshmpc,hecmatmpc)
290 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
292 if( hecmesh%n_dof == 3 )
then
293 call hecmw_update_3_r (hecmesh, hecmat%X, hecmat%NP)
294 if( hecmesh%my_rank == 0 )
then
295 write(
imsg, *)
'hecmw_update_3_R: OK'
297 else if( hecmesh%n_dof == 2 )
then
298 call hecmw_update_2_r (hecmesh, hecmat%X, hecmat%NP)
299 if( hecmesh%my_rank == 0 )
then
300 write(
imsg, *)
'hecmw_update_2_R: OK'
307 do i = 1, hecmesh%n_node*ndof
308 fstrsolid%dunode(i) = fstrsolid%dunode(i)+hecmat%X(i)
315 if( fstrsolid%DLOAD_follow /= 0 .or. fstrsolid%CLOAD_ngrp_rot /= 0 ) &
316 call fstr_ass_load(cstep, ctime+dtime, hecmesh, hecmat, fstrsolid, fstrparam)
327 res = sqrt(res)/n_node_global
328 if( iter == 1 ) res0 = res
329 if( res0 == 0.0d0 )
then
332 relres = dabs( res1-res )/res0
335 if( hecmesh%my_rank == 0 )
then
336 write(*,
'(a,i3,a,2e15.7)')
' - Residual(',iter,
') =', res, relres
340 if( res < fstrsolid%step_ctrl(cstep)%converg .or. &
341 relres < fstrsolid%step_ctrl(cstep)%converg )
exit
345 if( iter == fstrsolid%step_ctrl(cstep)%max_iter .or. res > fstrsolid%step_ctrl(cstep)%maxres .or. res /= res )
then
346 if( hecmesh%my_rank == 0)
then
347 write( *,
'(a,i5,a,i5)')
' ### Fail to Converge : at total_step=', cstep,
' sub_step=', sub_step
349 fstrsolid%NRstat_i(knstmaxit) = max(fstrsolid%NRstat_i(knstmaxit),iter)
350 fstrsolid%NRstat_i(knstsumit) = fstrsolid%NRstat_i(knstsumit) + iter
351 fstrsolid%NRstat_i(knstciter) = al_step
352 fstrsolid%CutBack_stat = fstrsolid%CutBack_stat + 1
353 if( iter == fstrsolid%step_ctrl(cstep)%max_iter ) fstrsolid%NRstat_i(knstdresn) = 1
354 if( res > fstrsolid%step_ctrl(cstep)%maxres .or. res /= res ) fstrsolid%NRstat_i(knstdresn) = 2
361 fstrsolid%NRstat_i(knstmaxit) = max(fstrsolid%NRstat_i(knstmaxit),iter)
362 fstrsolid%NRstat_i(knstsumit) = fstrsolid%NRstat_i(knstsumit) + iter
367 if(
associated(fstrsolid%contacts) )
then
369 call fstr_scan_contact_state( cstep, sub_step, al_step, dtime, ctalgo, hecmesh, fstrsolid, infoctchange, hecmat%B )
370 if( infoctchange%contact2free+infoctchange%contact2neighbor+infoctchange%free2contact > 0 ) &
380 do i=1,hecmesh%n_node*ndof
381 fstrsolid%unode(i) = fstrsolid%unode(i)+fstrsolid%dunode(i)
384 if( convg .and. (.not.ctchange) )
exit
387 if( al_step >= fstrsolid%step_ctrl(cstep)%max_contiter )
then
388 if( hecmesh%my_rank == 0)
then
389 write( *,
'(a,i5,a,i5)')
' ### Contact failed to Converge : at total_step=', cstep,
' sub_step=', sub_step
391 fstrsolid%NRstat_i(knstciter) = al_step
392 fstrsolid%CutBack_stat = fstrsolid%CutBack_stat + 1
393 fstrsolid%NRstat_i(knstdresn) = 3
400 fstrsolid%NRstat_i(knstciter) = al_step-1
405 fstrsolid%CutBack_stat = 0
406 call hecmw_mpc_mat_finalize(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
413 restart_step_num, restart_substep_num, sub_step, ctime, dtime, infoCTChange, conMAT )
419 integer,
intent(in) :: cstep
420 type (hecmwST_local_mesh) :: hecMESH
421 type (hecmwST_matrix) :: hecMAT
422 type (fstr_solid) :: fstrSOLID
423 integer,
intent(in) :: sub_step
424 real(kind=kreal),
intent(in) :: ctime
425 real(kind=kreal),
intent(in) :: dtime
426 type (fstr_param) :: fstrPARAM
427 type (fstr_info_contactChange) :: infoCTChange
428 type (fstrST_matrix_contact_lagrange) :: fstrMAT
429 type (hecmwST_matrix),
optional :: conMAT
431 type (hecmwST_local_mesh),
pointer :: hecMESHmpc
432 type (hecmwST_matrix),
pointer :: hecMATmpc
433 integer(kind=kint) :: ndof
434 integer(kind=kint) :: ctAlgo
435 integer(kind=kint) :: i, iter, max_iter_contact
436 integer(kind=kint) :: stepcnt, count_step
437 real(kind=kreal) :: tt0, tt, res, res0, res1, relres, tincr, resx
438 integer(kind=kint) :: restart_step_num, restart_substep_num
439 logical :: is_mat_symmetric
440 integer(kind=kint) :: n_node_global
441 integer(kind=kint) :: contact_changed_global
442 integer(kint) :: nndof
443 real(kreal) :: q_residual,x_residual
444 real(kind=kreal),
allocatable :: coord(:)
445 integer(kind=kint) :: istat
447 call hecmw_mpc_mat_init(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
450 n_node_global = hecmesh%nn_internal
451 call hecmw_allreduce_i1(hecmesh,n_node_global,hecmw_sum)
454 write(*, *)
' This type of direct solver is not yet available in such case ! '
455 write(*, *)
' Please use intel MKL direct solver !'
456 call hecmw_abort( hecmw_comm_get_comm() )
459 ctalgo = fstrparam%contact_algo
461 hecmat%NDOF = hecmesh%n_dof
464 fstrsolid%NRstat_i(:) = 0
466 allocate(coord(hecmesh%n_node*ndof))
469 if( fstrsolid%step_ctrl(cstep)%solution == stepstatic ) tincr = 0.0d0
471 fstrsolid%dunode(:) = 0.0d0
473 if( cstep==1 .and. sub_step==restart_substep_num )
then
477 call hecmw_mat_copy_profile( hecmat, conmat )
486 elseif( hecmat%Iarray(99)==4 )
then
487 write(*, *)
' This type of direct solver is not yet available in such case ! '
488 write(*, *)
' Please change the solver type to intel MKL direct solver !'
489 call hecmw_abort(hecmw_comm_get_comm())
497 call fstr_ass_load(cstep, ctime+dtime, hecmesh, hecmat, fstrsolid, fstrparam)
500 call hecmw_mat_clear_b(conmat)
512 fstrsolid%dunode(:) = 0.0d0
516 loopforcontactanalysis:
do while( .true. )
517 count_step = count_step+1
524 do iter = 1, fstrsolid%step_ctrl(cstep)%max_iter
525 call hecmw_barrier(hecmesh)
526 if(
myrank == 0 ) print *,
'-------------------------------------------------'
527 call hecmw_barrier(hecmesh)
534 call hecmw_mat_clear( conmat )
548 call hecmw_mpc_mat_ass(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
549 call hecmw_mpc_trans_rhs(hecmesh, hecmat, hecmatmpc)
551 call fstr_addbc(cstep, hecmesh, hecmatmpc, fstrsolid, fstrparam, fstrmat, stepcnt, conmat)
553 call fstr_addbc(cstep, hecmesh, hecmatmpc, fstrsolid, fstrparam, fstrmat, stepcnt)
556 nndof = hecmat%N*hecmat%ndof
570 call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
572 if( istat /= 0 )
then
573 if( hecmesh%my_rank == 0)
then
574 write( *,
'(a,i5,a,i5)')
' ### Fail to Converge : at total_step=', cstep,
' sub_step=', sub_step
576 fstrsolid%NRstat_i(knstdresn) = 4
577 fstrsolid%CutBack_stat = fstrsolid%CutBack_stat + 1
583 call hecmw_innerproduct_r(hecmesh,ndof,hecmat%X,hecmat%X,resx)
584 resx = sqrt(resx)/n_node_global
586 if( hecmesh%my_rank==0 )
then
587 write(*,
'(a,i3,a,e15.7)')
' - ResidualX (',iter,
') =',resx
588 write(*,
'(a,i3,a,e15.7)')
' - ResidualX+LAG(',iter,
') =',sqrt(x_residual)/n_node_global
589 write(*,
'(a,i3,a,e15.7)')
' - ResidualQ (',iter,
') =',sqrt(q_residual)/n_node_global
593 do i = 1, hecmesh%n_node*ndof
594 fstrsolid%dunode(i) = fstrsolid%dunode(i) + hecmat%X(i)
599 do i = 1, fstrmat%num_lagrange
600 fstrmat%lagrange(i) = fstrmat%lagrange(i)+hecmat%X(hecmesh%n_node*ndof+i)
608 if( fstrsolid%DLOAD_follow /= 0 .or. fstrsolid%CLOAD_ngrp_rot /= 0 ) &
609 call fstr_ass_load(cstep, ctime+dtime, hecmesh, hecmat, fstrsolid, fstrparam)
619 call hecmw_mat_clear_b( conmat )
632 res = sqrt(res)/n_node_global
633 if( iter == 1 ) res0 = res
634 if( res0 == 0.0d0 )
then
637 relres = dabs( res1-res )/res0
639 if( hecmesh%my_rank == 0 )
then
640 write(*,
'(a,i3,a,2e15.7)')
' - Residual(',iter,
') =',res,relres
644 if( res < fstrsolid%step_ctrl(cstep)%converg .or. &
645 relres < fstrsolid%step_ctrl(cstep)%converg )
exit
649 if( iter == fstrsolid%step_ctrl(cstep)%max_iter .or. res > fstrsolid%step_ctrl(cstep)%maxres .or. res /= res )
then
650 if( hecmesh%my_rank == 0)
then
651 write( *,
'(a,i5,a,i5)')
' ### Fail to Converge : at total_step=', cstep,
' sub_step=', sub_step
653 fstrsolid%NRstat_i(knstmaxit) = max(fstrsolid%NRstat_i(knstmaxit),iter)
654 fstrsolid%NRstat_i(knstsumit) = fstrsolid%NRstat_i(knstsumit) + iter
655 fstrsolid%NRstat_i(knstciter) = count_step
656 fstrsolid%CutBack_stat = fstrsolid%CutBack_stat + 1
657 if( iter == fstrsolid%step_ctrl(cstep)%max_iter ) fstrsolid%NRstat_i(knstdresn) = 1
658 if( res > fstrsolid%step_ctrl(cstep)%maxres .or. res /= res ) fstrsolid%NRstat_i(knstdresn) = 2
665 fstrsolid%NRstat_i(knstmaxit) = max(fstrsolid%NRstat_i(knstmaxit),iter)
666 fstrsolid%NRstat_i(knstsumit) = fstrsolid%NRstat_i(knstsumit) + iter
668 call fstr_scan_contact_state( cstep, sub_step, count_step, dtime, ctalgo, hecmesh, fstrsolid, infoctchange, hecmat%B )
671 write(*, *)
' This type of direct solver is not yet available in such case ! '
672 write(*, *)
' Please use intel MKL direct solver !'
673 call hecmw_abort( hecmw_comm_get_comm() )
677 contact_changed_global = 0
685 contact_changed_global = 1
690 call hecmw_allreduce_i1(hecmesh, contact_changed_global, hecmw_max)
692 if (contact_changed_global > 0)
then
693 call hecmw_mat_clear_b( hecmat )
695 call hecmw_mat_clear_b( conmat )
701 if( count_step >= fstrsolid%step_ctrl(cstep)%max_contiter )
then
702 if( hecmesh%my_rank == 0)
then
703 write( *,
'(a,i5,a,i5)')
' ### Contact failed to Converge : at total_step=', cstep,
' sub_step=', sub_step
705 fstrsolid%NRstat_i(knstciter) = count_step
706 fstrsolid%CutBack_stat = fstrsolid%CutBack_stat + 1
707 fstrsolid%NRstat_i(knstdresn) = 3
720 call hecmw_mat_clear_b( conmat )
727 enddo loopforcontactanalysis
729 fstrsolid%NRstat_i(knstciter) = count_step
733 do i = 1, hecmesh%n_node*ndof
734 fstrsolid%unode(i) = fstrsolid%unode(i)+fstrsolid%dunode(i)
741 fstrsolid%CutBack_stat = 0
742 call hecmw_mpc_mat_finalize(hecmesh, hecmat, hecmeshmpc, hecmatmpc)
This module provides a function to deal with prescribed displacement.
subroutine fstr_addbc(cstep, hecmesh, hecmat, fstrsolid, fstrparam, fstrmat, iter, conmat)
Add Essential Boundary Conditions.
This module provides functions to take into acount external load.
subroutine fstr_ass_load(cstep, ctime, hecmesh, hecmat, fstrsolid, fstrparam)
This subroutine assmble following external force into fstrSOLIDGL and hecMATB afterwards.
This module provides functions on nonlinear analysis.
subroutine fstr_newton_contactslag(cstep, hecmesh, hecmat, fstrsolid, fstrparam, fstrmat, restart_step_num, restart_substep_num, sub_step, ctime, dtime, infoctchange, conmat)
This subroutine solve nonlinear solid mechanics problems by Newton-Raphson method....
subroutine fstr_newton(cstep, hecmesh, hecmat, fstrsolid, fstrparam, restrt_step_num, sub_step, ctime, dtime)
This subroutine solve nonlinear solid mechanics problems by Newton-Raphson method.
subroutine fstr_newton_contactalag(cstep, hecmesh, hecmat, fstrsolid, fstrparam, restart_step_num, restart_substep_num, sub_step, ctime, dtime, infoctchange)
This subroutine solve nonlinear solid mechanics problems by Newton-Raphson method combined with Neste...
This module provides function to calcualte residual of nodal force.
real(kind=kreal) function, public fstr_get_x_norm_contact(hecmat, fstrmat, hecmesh)
real(kind=kreal) function, public fstr_get_norm_contact(flag, hecmesh, hecmat, fstrsolid, fstrmat)
Calculate square norm.
real(kind=kreal) function, public fstr_get_residual(force, hecmesh)
Calculate magnitude of a real vector.
real(kind=kreal) function, public fstr_get_norm_para_contact(hecmat, fstrmat, conmat, hecmesh)
subroutine, public fstr_update_ndforce(cstep, hecmesh, hecmat, fstrsolid, conmat)
This module provides functions to read in and write out restart fiels.
This module provides functions to deal with spring force.
subroutine fstr_addspring(cstep, hecmesh, hecmat, fstrsolid, fstrparam)
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.
integer(kind=kint) myrank
PARALLEL EXECUTION.
subroutine fstr_set_current_config_to_mesh(hecmesh, fstrsolid, coord)
integer(kind=kint), parameter imsg
integer(kind=kint), parameter ilog
FILE HANDLER.
type(fstr_param), target fstrpr
GLOBAL VARIABLE INITIALIZED IN FSTR_SETUP.
subroutine fstr_recover_initial_config_to_mesh(hecmesh, fstrsolid, coord)
integer(kind=kint), parameter kno
logical paracontactflag
PARALLEL CONTACT FLAG.
This modules just summarizes all modules used in static analysis.
This module provides functions to output result.