10 include
'fstr_ctrl_util_f.inc'
27 character( len=80 ) :: control
28 character( len=80 ) :: convcontrol
30 real(kind=kreal) :: converg
31 real(kind=kreal) :: maxres
33 integer :: num_substep
35 integer :: max_contiter
37 real(kind=kreal) :: initdt
38 real(kind=kreal) :: elapsetime
39 real(kind=kreal) :: mindt
40 real(kind=kreal) :: maxdt
41 real(kind=kreal) :: starttime
42 integer,
pointer :: boundary(:)=>null()
43 integer,
pointer :: load(:)=>null()
44 integer,
pointer :: contact(:)=>null()
45 integer :: timepoint_id
46 integer :: aincparam_id
50 character(HECMW_NAME_LEN) :: name
51 real(kind=kreal) :: ainc_rs
52 real(kind=kreal) :: ainc_rl
53 integer( kind=kint ) :: nrbound_s(10)
54 integer( kind=kint ) :: nrbound_l(10)
55 integer( kind=kint ) :: nrtimes_s
56 integer( kind=kint ) :: nrtimes_l
57 real(kind=kreal) :: ainc_rc
58 integer( kind=kint ) :: cbbound
65 type(
step_info ),
intent(out) :: stepinfo
68 stepinfo%num_substep = 1
69 stepinfo%max_iter = 50
70 stepinfo%max_contiter = 10
72 stepinfo%initdt = 1.d0
73 stepinfo%mindt = 1.d-4
75 stepinfo%elapsetime = 1.d0
76 stepinfo%starttime = 0.d0
77 stepinfo%converg = 1.d-3
78 stepinfo%maxres = 1.d+10
79 stepinfo%timepoint_id = 0
80 stepinfo%AincParam_id = 0
84 type(
step_info ),
pointer,
intent(inout) :: stepinfos(:)
88 stepinfos(1)%starttime = 0.d0
89 do i=1,
size(stepinfos)-1
90 stepinfos(i+1)%starttime = stepinfos(i)%starttime + stepinfos(i)%elapsetime
96 integer,
intent(in) :: bnd
99 if( .not.
associated( stepinfo%Boundary ) )
return
105 integer,
intent(in) :: bnd
106 type(
step_info ),
intent(in) :: stepinfo
108 if( .not.
associated( stepinfo%Load ) )
return
114 integer,
intent(in) :: bnd
115 type(
step_info ),
intent(in) :: stepinfo
117 if( .not.
associated( stepinfo%Contact ) )
return
124 if(
associated( step%Boundary ) )
deallocate( step%Boundary )
125 if(
associated( step%Load ) )
deallocate( step%Load )
126 if(
associated( step%Contact ) )
deallocate( step%Contact )
131 integer,
intent(in) :: nfile
133 integer :: i, j, nstep, nbc
136 write( nfile, * )
"-----Information of steps:",nstep
139 write( nfile, * )
" -----Step:",i
140 write(nfile,*) steps(i)%solution, steps(i)%elapsetime, steps(i)%converg, &
141 steps(i)%num_substep, steps(i)%max_iter
142 if(
associated( steps(i)%Boundary ) )
then
143 nbc =
size( steps(i)%Boundary )
144 write(nfile,*)
" Boundary conditions"
145 write(nfile,*) ( steps(i)%Boundary(j),j=1,nbc )
147 if(
associated( steps(i)%Load ) )
then
148 nbc =
size( steps(i)%Load )
149 write(nfile,*)
" External load conditions"
150 write(nfile,*) ( steps(i)%Load(j),j=1,nbc )
152 if(
associated( steps(i)%Contact ) )
then
153 nbc =
size( steps(i)%Contact )
154 write(nfile,*)
" Contact conditions"
155 write(nfile,*) ( steps(i)%Contact(j),j=1,nbc )
165 aincparam%ainc_Rs = 0.25d0
166 aincparam%ainc_Rl = 1.25d0
167 aincparam%NRbound_s = 0
171 aincparam%NRbound_l = 0
175 aincparam%NRtimes_s = 1
176 aincparam%NRtimes_l = 2
177 aincparam%ainc_Rc = 0.25d0
178 aincparam%CBbound = 5
This module manages step infomation.
logical function iscontactactive(bnd, stepinfo)
Is contact condition in this step active.
subroutine fstr_print_steps(nfile, steps)
Print out step control.
subroutine free_stepinfo(step)
Finalizer.
subroutine init_stepinfo(stepinfo)
Initializer.
integer, parameter stepvisco
integer, parameter stepfixedinc
subroutine init_aincparam(aincparam)
Initializer.
logical function isboundaryactive(bnd, stepinfo)
Is boundary condition in this step active.
integer(kind=kint), parameter knstmaxit
integer, parameter stepautoinc
integer(kind=kint), parameter knstdresn
integer(kind=kint), parameter knstsumit
integer(kind=kint), parameter knstciter
subroutine setup_stepinfo_starttime(stepinfos)
logical function isloadactive(bnd, stepinfo)
Is external load in this step active.
integer, parameter stepstatic
Step control such as active boundary condition, convergent condition etc.