FrontISTR  5.7.1
Large-scale structural analysis program with finit element method
fstr_ctrl_dynamic.f90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
2 ! Copyright (c) 2019 FrontISTR Commos
3 ! This software is released under the MIT License, see LICENSE.txt
4 !-------------------------------------------------------------------------------
6 
8  use m_fstr
9  use hecmw
11  private :: fstr_ctrl_get_nval
12 contains
13 
14  !* ----------------------------------------------------------------------------------------------- *!
15  !* fstr_ctrl_get_nval *!
16  !* ----------------------------------------------------------------------------------------------- *!
17 
18  function fstr_ctrl_get_nval( ctrl, iType, amp, node_id, node_id_len, dof_ids, dof_ide, value )
19  implicit none
20  integer(kind=kint) :: ctrl
21  integer(kind=kint) :: iType
22  character(len=HECMW_NAME_LEN) :: amp
23  character(len=HECMW_NAME_LEN) :: node_id(:)
24  integer(kind=kint) :: node_id_len
25  integer(kind=kint),pointer :: dof_ids (:)
26  integer(kind=kint),pointer :: dof_ide (:)
27  real(kind=kreal),pointer :: value(:)
28  integer(kind=kint) :: fstr_ctrl_get_nval
29  integer(kind=kint) :: rcode
30 
31  character(len=HECMW_NAME_LEN) :: data_fmt,ss
32  write(ss,*) node_id_len
33  write(data_fmt,'(a,a,a)') 'S',trim(adjustl(ss)),'IIr '
34 
35  itype=2
36  ss='INITIAL,TRANSIT '
37  fstr_ctrl_get_nval = -1
38  rcode=fstr_ctrl_get_param_ex( ctrl, 'TYPE ', ss, 0, 'P', itype )
39 
40  if( fstr_ctrl_get_param_ex( ctrl, 'AMP ', '# ', 0, 'S', amp )/= 0) return
41  fstr_ctrl_get_nval = &
42  fstr_ctrl_get_data_array_ex( ctrl, data_fmt, node_id, dof_ids, dof_ide, value )
43 
44  end function fstr_ctrl_get_nval
45 
47  function fstr_ctrl_get_velocity( ctrl, vType, amp, node_id, node_id_len, dof_ids, dof_ide, value )
48  implicit none
49  integer(kind=kint) :: ctrl
50  integer(kind=kint) vtype
51  character(len=HECMW_NAME_LEN) :: amp
52  character(len=HECMW_NAME_LEN) :: node_id(:)
53  integer(kind=kint) :: node_id_len
54  integer(kind=kint),pointer :: dof_ids (:)
55  integer(kind=kint),pointer :: dof_ide (:)
56  real(kind=kreal),pointer :: value(:)
57  integer(kind=kint) :: fstr_ctrl_get_velocity
58 
60  fstr_ctrl_get_nval( ctrl, vtype, amp, node_id, node_id_len, dof_ids, dof_ide, value )
61 
62  end function fstr_ctrl_get_velocity
63 
65  function fstr_ctrl_get_acceleration( ctrl, aType, amp, node_id, node_id_len, dof_ids, dof_ide, value )
66  implicit none
67  integer(kind=kint) :: ctrl
68  integer(kind=kint) :: atype
69  character(len=HECMW_NAME_LEN) :: amp
70  character(len=HECMW_NAME_LEN) :: node_id(:)
71  integer(kind=kint) :: node_id_len
72  integer(kind=kint),pointer :: dof_ids (:)
73  integer(kind=kint),pointer :: dof_ide (:)
74  real(kind=kreal),pointer :: value(:)
75  integer(kind=kint) :: fstr_ctrl_get_acceleration
76 
78  fstr_ctrl_get_nval( ctrl, atype, amp, node_id, node_id_len, dof_ids, dof_ide, value )
79 
80  end function fstr_ctrl_get_acceleration
81 
83  function fstr_ctrl_get_dynamic( ctrl, nlgeom, &
84  idx_eqa, idx_resp, n_step, t_start, t_end, t_delta, &
85  gamma, beta, idx_mas, idx_dmp, ray_m, ray_k, &
86  nout, node_id, node_id_len, nout_monit, iout_list )
87  implicit none
88  integer(kind=kint) :: ctrl
89 
90  ! ANALYSIS TYPE CONTROL
91  logical :: nlgeom
92  integer(kind=kint) :: idx_eqa
93  integer(kind=kint) :: idx_resp
94 
95  ! TIME CONTROL
96  integer(kind=kint) :: n_step
97  real(kind=kreal) :: t_start
98  real(kind=kreal) :: t_end
99  real(kind=kreal) :: t_delta
100 
101  ! Newmark-beta parameter
102  real(kind=kreal) :: gamma
103  real(kind=kreal) :: beta
104 
105  ! mass matrix control
106  integer(kind=kint) :: idx_mas
107 
108  ! damping control
109  integer(kind=kint) :: idx_dmp
110  real(kind=kreal) :: ray_m
111  real(kind=kreal) :: ray_k
112 
113  ! OUTPUT CONTROL
114  integer(kind=kint) :: nout
115  character(len=HECMW_NAME_LEN) :: node_id
116  integer(kind=kint) :: node_id_len
117  integer(kind=kint) :: nout_monit
118  integer(kind=kint) :: iout_list(6)
119 
120  integer(kind=kint) :: rcode, nlflag
121  character(len=80) :: s
122  character(len=HECMW_NAME_LEN) :: data_fmt,ss
123 
124  integer(kind=kint) :: fstr_ctrl_get_dynamic
125 
127 
128 
129  s = 'LINEAR,NONLINEAR '
130  nlflag=0
131  rcode = fstr_ctrl_get_param_ex( ctrl, 'TYPE ', s, 1, 'P', nlflag )
132  if( nlflag/=0 ) nlgeom = (nlflag==2)
133 
134  if( fstr_ctrl_get_data_ex( ctrl, 1, 'ii ', idx_eqa, idx_resp )/=0 ) return
135  if( fstr_ctrl_get_data_ex( ctrl, 2, 'rrir ', t_start, t_end, n_step, t_delta )/=0 ) return
136  if( fstr_ctrl_get_data_ex( ctrl, 3, 'rr ', gamma, beta )/=0 ) return
137  if( fstr_ctrl_get_data_ex( ctrl, 4, 'iirr ', idx_mas, idx_dmp, ray_m, ray_k )/=0 ) return
138  write(ss,*) node_id_len
139  write(data_fmt,'(a,a,a)') 'iS',trim(adjustl(ss)),'i '
140  if( fstr_ctrl_get_data_ex( ctrl, 5, data_fmt, nout, node_id, nout_monit )/=0 ) return
141  if( fstr_ctrl_get_data_ex( ctrl, 6, 'iiiiii ', &
142  iout_list(1), iout_list(2), iout_list(3), iout_list(4), iout_list(5), iout_list(6) )/=0 ) return
144 
145  end function fstr_ctrl_get_dynamic
146 
147 end module fstr_ctrl_dynamic
int fstr_ctrl_get_param_ex(int *ctrl, const char *param_name, const char *value_list, int *necessity, char *type, void *val)
int fstr_ctrl_get_data_array_ex(int *ctrl, const char *format,...)
int fstr_ctrl_get_data_ex(int *ctrl, int *line_no, const char *format,...)
This module contains control file data obtaining functions for dynamic analysis.
integer(kind=kint) function fstr_ctrl_get_dynamic(ctrl, nlgeom, idx_eqa, idx_resp, n_step, t_start, t_end, t_delta, gamma, beta, idx_mas, idx_dmp, ray_m, ray_k, nout, node_id, node_id_len, nout_monit, iout_list)
Read in !DYNAMIC.
integer(kind=kint) function fstr_ctrl_get_velocity(ctrl, vType, amp, node_id, node_id_len, dof_ids, dof_ide, value)
Read in !VELOCITY.
integer(kind=kint) function fstr_ctrl_get_acceleration(ctrl, aType, amp, node_id, node_id_len, dof_ids, dof_ide, value)
Read in !ACCELERATION.
Definition: hecmw.f90:6
This module defines common data and basic structures for analysis.
Definition: m_fstr.f90:15