FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
m_step.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!-------------------------------------------------------------------------------
6module m_step
7 use hecmw
8 implicit none
9
10 include 'fstr_ctrl_util_f.inc'
11
12 integer, parameter :: stepstatic = 1
13 integer, parameter :: stepvisco = 2
14 integer, parameter :: stepfixedinc = 1
15 integer, parameter :: stepautoinc = 2
16
17 ! statistics of newton iteration
18 integer(kind=kint),parameter :: knstmaxit = 1 ! maximum number of newton iteration
19 integer(kind=kint),parameter :: knstsumit = 2 ! total number of newton iteration
20 integer(kind=kint),parameter :: knstciter = 3 ! number of contact iteration
21 integer(kind=kint),parameter :: knstdresn = 4 ! reason of not to converged
22
25 integer :: solution
26 integer :: inc_type
27 character( len=80 ) :: control
28 character( len=80 ) :: convcontrol
30 real(kind=kreal) :: converg
31 real(kind=kreal) :: maxres
32
33 integer :: num_substep
34 integer :: max_iter
35 integer :: max_contiter
36 integer :: amp_id
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
47 end type
48
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
59 end type
60
61contains
62
64 subroutine init_stepinfo( stepinfo )
65 type( step_info ), intent(out) :: stepinfo
66 stepinfo%solution = stepstatic
67 stepinfo%inc_type = stepfixedinc
68 stepinfo%num_substep = 1
69 stepinfo%max_iter = 50
70 stepinfo%max_contiter = 10
71 stepinfo%amp_id = -1
72 stepinfo%initdt = 1.d0
73 stepinfo%mindt = 1.d-4
74 stepinfo%maxdt = 1.d0
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
81 end subroutine
82
83 subroutine setup_stepinfo_starttime( stepinfos )
84 type( step_info ), pointer, intent(inout) :: stepinfos(:)
85
86 integer :: i
87
88 stepinfos(1)%starttime = 0.d0
89 do i=1,size(stepinfos)-1
90 stepinfos(i+1)%starttime = stepinfos(i)%starttime + stepinfos(i)%elapsetime
91 end do
92 end subroutine
93
95 logical function isboundaryactive( bnd, stepinfo )
96 integer, intent(in) :: bnd
97 type( step_info ), intent(in) :: stepinfo
98 isboundaryactive = .false.
99 if( .not. associated( stepinfo%Boundary ) ) return
100 if( any( stepinfo%Boundary== bnd ) ) isboundaryactive = .true.
101 end function
102
104 logical function isloadactive( bnd, stepinfo )
105 integer, intent(in) :: bnd
106 type( step_info ), intent(in) :: stepinfo
107 isloadactive = .false.
108 if( .not. associated( stepinfo%Load ) ) return
109 if( any( stepinfo%Load == bnd ) ) isloadactive = .true.
110 end function
111
113 logical function iscontactactive( bnd, stepinfo )
114 integer, intent(in) :: bnd
115 type( step_info ), intent(in) :: stepinfo
116 iscontactactive = .false.
117 if( .not. associated( stepinfo%Contact ) ) return
118 if( any( stepinfo%Contact== bnd ) ) iscontactactive = .true.
119 end function
120
122 subroutine free_stepinfo( step )
123 type(step_info), intent(inout) :: step
124 if( associated( step%Boundary ) ) deallocate( step%Boundary )
125 if( associated( step%Load ) ) deallocate( step%Load )
126 if( associated( step%Contact ) ) deallocate( step%Contact )
127 end subroutine
128
130 subroutine fstr_print_steps( nfile, steps )
131 integer, intent(in) :: nfile
132 type(step_info), intent(in) :: steps(:)
133 integer :: i, j, nstep, nbc
134 nstep = size(steps)
135
136 write( nfile, * ) "-----Information of steps:",nstep
137
138 do i=1,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 )
146 endif
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 )
151 endif
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 )
156 endif
157 enddo
158 end subroutine
159
161 subroutine init_aincparam( aincparam )
162 type( tparamautoinc ), intent(out) :: aincparam
163
164 aincparam%name = ''
165 aincparam%ainc_Rs = 0.25d0
166 aincparam%ainc_Rl = 1.25d0
167 aincparam%NRbound_s = 0
168 aincparam%NRbound_s(knstmaxit) = 10
169 aincparam%NRbound_s(knstsumit) = 50
170 aincparam%NRbound_s(knstciter) = 10
171 aincparam%NRbound_l = 0
172 aincparam%NRbound_l(knstmaxit) = 1
173 aincparam%NRbound_l(knstsumit) = 1
174 aincparam%NRbound_l(knstciter) = 1
175 aincparam%NRtimes_s = 1
176 aincparam%NRtimes_l = 2
177 aincparam%ainc_Rc = 0.25d0
178 aincparam%CBbound = 5
179 end subroutine
180
181end module
Definition: hecmw.f90:6
This module manages step infomation.
Definition: m_step.f90:6
logical function iscontactactive(bnd, stepinfo)
Is contact condition in this step active.
Definition: m_step.f90:114
subroutine fstr_print_steps(nfile, steps)
Print out step control.
Definition: m_step.f90:131
subroutine free_stepinfo(step)
Finalizer.
Definition: m_step.f90:123
subroutine init_stepinfo(stepinfo)
Initializer.
Definition: m_step.f90:65
integer, parameter stepvisco
Definition: m_step.f90:13
integer, parameter stepfixedinc
Definition: m_step.f90:14
subroutine init_aincparam(aincparam)
Initializer.
Definition: m_step.f90:162
logical function isboundaryactive(bnd, stepinfo)
Is boundary condition in this step active.
Definition: m_step.f90:96
integer(kind=kint), parameter knstmaxit
Definition: m_step.f90:18
integer, parameter stepautoinc
Definition: m_step.f90:15
integer(kind=kint), parameter knstdresn
Definition: m_step.f90:21
integer(kind=kint), parameter knstsumit
Definition: m_step.f90:19
integer(kind=kint), parameter knstciter
Definition: m_step.f90:20
subroutine setup_stepinfo_starttime(stepinfos)
Definition: m_step.f90:84
logical function isloadactive(bnd, stepinfo)
Is external load in this step active.
Definition: m_step.f90:105
integer, parameter stepstatic
Definition: m_step.f90:12
Step control such as active boundary condition, convergent condition etc.
Definition: m_step.f90:24