FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_ctrl_freq.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!-------------------------------------------------------------------------------
7
8!-----------------------------------------------------------------------------!
10!-----------------------------------------------------------------------------!
11 subroutine fstr_setup_fload( ctrl, counter, P )
12 !---- args
13 integer(kind=kint) :: ctrl
14 integer(kind=kint) :: counter
15 type(fstr_param_pack) :: P
16 !---- vals
17 integer(kind=kint) :: rcode
18 character(HECMW_NAME_LEN) :: amp
19 integer(kind=kint) :: amp_id
20 character(HECMW_NAME_LEN), pointer :: grp_id_name(:)
21 real(kind=kreal), pointer :: val_ptr(:)
22 integer(kind=kint), pointer :: id_ptr(:), type_ptr(:)
23 integer(kind=kint) :: i, n, old_size, new_size
24 integer(kind=kint) :: gid, loadcase
25 !---- body
26
27 if( p%SOLID%file_type /= kbcffstr) return
28
29 !read grpid
30 gid = 1
31 rcode = fstr_ctrl_get_param_ex( ctrl, 'GRPID ', '# ', 0, 'I', gid )
32 !read loadcase (real=1:default, img=2)
33 loadcase = kfloadcase_re
34 rcode = fstr_ctrl_get_param_ex( ctrl, 'LOAD CASE ', '# ', 0, 'I', loadcase)
35 !write(*,*) "loadcase=", loadcase
36 !pause
37
38 !read the num of dataline
40 if( n == 0 ) return
41 old_size = p%FREQ%FLOAD_ngrp_tot
42 new_size = old_size + n
43
44 !expand data array
45 p%FREQ%FLOAD_ngrp_tot = new_size
46 call fstr_expand_integer_array( p%FREQ%FLOAD_ngrp_GRPID, old_size, new_size )
47 call fstr_expand_integer_array( p%FREQ%FLOAD_ngrp_ID, old_size, new_size )
48 call fstr_expand_integer_array( p%FREQ%FLOAD_ngrp_TYPE, old_size, new_size )
49 call fstr_expand_integer_array( p%FREQ%FLOAD_ngrp_DOF, old_size, new_size )
50 call fstr_expand_real_array ( p%FREQ%FLOAD_ngrp_valre, old_size, new_size )
51 call fstr_expand_real_array ( p%FREQ%FLOAD_ngrp_valim, old_size, new_size )
52
53 !fill bc data
54 allocate( grp_id_name(n) )
55 if(loadcase == kfloadcase_re) then
56 val_ptr => p%FREQ%FLOAD_ngrp_valre(old_size+1:)
57 else if(loadcase == kfloadcase_im) then
58 val_ptr => p%FREQ%FLOAD_ngrp_valim(old_size+1:)
59 else
60 !error
61 write(*,*) "Error this load set is not defined!"
62 write(ilog,*) "Error this load set is not defined!"
63 stop
64 end if
65 id_ptr => p%FREQ%FLOAD_ngrp_DOF(old_size+1:)
66 type_ptr => p%FREQ%FLOAD_ngrp_TYPE(old_size+1:)
67 val_ptr = 0.0d0
68 rcode = fstr_ctrl_get_fload( ctrl, grp_id_name, hecmw_name_len, id_ptr, val_ptr)
69 if( rcode /= 0 ) call fstr_ctrl_err_stop
70 p%FREQ%FLOAD_ngrp_GRPID(old_size+1:new_size) = gid
71 call nodesurf_grp_name_to_id_ex( p%MESH, '!FLOAD', n, grp_id_name, &
72 p%FREQ%FLOAD_ngrp_ID(old_size+1:), p%FREQ%FLOAD_ngrp_TYPE(old_size+1:))
73
74 deallocate( grp_id_name )
75 return
76
77 contains
78
79 function fstr_ctrl_get_fload(ctrl, node_id, node_id_len, dof_id, value)
80! include 'fstr_ctrl_util_f.inc' !Fortran->C interface for fstr ctrl API
81 integer(kind=kint) :: ctrl
82 character(len=HECMW_NAME_LEN),target :: node_id(:) !Node group name
83 integer(kind=kint), pointer :: dof_id(:)
84 integer(kind=kint) :: node_id_len
85 real(kind=kreal), pointer :: value(:)
86 integer(kind=kint) :: fstr_ctrl_get_FLOAD !return value
87 character(len=HECMW_NAME_LEN) :: data_fmt, ss
88 character(len=HECMW_NAME_LEN),pointer :: node_id_p
89
90 write(ss,*) node_id_len
91 write(data_fmt, '(a,a,a)') 'S', trim(adjustl(ss)), 'IR '
92
93 node_id_p => node_id(1)
94 fstr_ctrl_get_fload = fstr_ctrl_get_data_array_ex(ctrl, data_fmt, node_id_p, dof_id, value)
95 end function
96
97 end subroutine
98
99!-----------------------------------------------------------------------------!
101!-----------------------------------------------------------------------------!
102 subroutine fstr_setup_eigenread( ctrl, counter, P )
103 !---- args
104 integer(kind=kint) :: ctrl
105 integer(kind=kint) :: counter
106 type(fstr_param_pack) :: P
107 !---- vals
108 integer(kind=kint) :: filename_len
109 character(len=HECMW_NAME_LEN) :: datafmt, ss
110 !---- body
111
112 filename_len = hecmw_filename_len
113 write(ss,*) filename_len
114 write(datafmt, '(a,a,a)') 'S', trim(adjustl(ss)), ' '
115
116 if( fstr_ctrl_get_data_ex( ctrl, 1, datafmt, p%FREQ%eigenlog_filename ) /= 0) return
117 if( fstr_ctrl_get_data_ex( ctrl, 2, 'ii ', p%FREQ%start_mode, p%FREQ%end_mode ) /= 0) return
118
119 return
120
121 end subroutine
122
subroutine fstr_setup_eigenread(ctrl, counter, p)
Read in !EIGENREAD !
subroutine fstr_setup_fload(ctrl, counter, p)
This source file contains subroutine for reading control data for harmonic response analysis (this im...
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_line_n(int *ctrl)
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,...)