FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_ctrl_heat.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 use m_fstr
8 use hecmw
9
10 include 'fstr_ctrl_util_f.inc'
11
12 private :: pc_strupr
13
14contains
15
16 subroutine pc_strupr( s )
17 implicit none
18 character(*) :: s
19 integer :: i, n, a, da
20
21 n = len_trim(s)
22 da = iachar('a') - iachar('A')
23 do i = 1, n
24 a = iachar(s(i:i))
25 if( a > iachar('Z')) then
26 a = a - da
27 s(i:i) = achar(a)
28 end if
29 end do
30 end subroutine pc_strupr
31
32
33
35 function fstr_ctrl_get_heat( ctrl, dt, etime, dtmin, deltmx, itmax, eps, tpname )
36 implicit none
37 integer(kind=kint) :: ctrl
38 real(kind=kreal),pointer :: dt(:)
39 real(kind=kreal),pointer :: etime(:)
40 real(kind=kreal),pointer :: dtmin(:)
41 real(kind=kreal),pointer :: deltmx(:)
42 integer(kind=kint),pointer :: itmax(:)
43 real(kind=kreal),pointer :: eps(:)
44 character(len=*), intent(out) :: tpname
45 integer(kind=kint) :: fstr_ctrl_get_heat
46
47 integer(kind=kint) :: result
48
50
51 tpname=""
52 if( fstr_ctrl_get_param_ex( ctrl, 'TIMEPOINTS ', '# ', 0, 'S', tpname )/= 0) return
53
54 ! JP-7
55 if( fstr_ctrl_get_data_array_ex( ctrl, 'rrrrir ', dt, etime, dtmin, deltmx, itmax, eps )/= 0) return
56
58 end function fstr_ctrl_get_heat
59
61 function fstr_ctrl_get_fixtemp( ctrl, amp, node_grp_name, node_grp_name_len, value )
62 implicit none
63 integer(kind=kint) :: ctrl
64 character(len=HECMW_NAME_LEN) :: amp
65 character(len=HECMW_NAME_LEN),target :: node_grp_name(:)
66 character(len=HECMW_NAME_LEN),pointer:: node_grp_name_p
67 integer(kind=kint) :: node_grp_name_len
68 real(kind=kreal), pointer :: value(:)
69 integer(kind=kint) :: fstr_ctrl_get_fixtemp
70
71 character(len=HECMW_NAME_LEN) :: data_fmt,ss
72
74
75 ! JP-8
76 if( fstr_ctrl_get_param_ex( ctrl, 'AMP ', '# ', 0, 'S', amp )/= 0) return
77
78 write(ss,*) node_grp_name_len
79 write(data_fmt,'(a,a,a)') 'S',trim(adjustl(ss)),'r '
80
81 node_grp_name_p => node_grp_name(1)
82 fstr_ctrl_get_fixtemp = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, node_grp_name_p, value )
83
84 end function fstr_ctrl_get_fixtemp
85
86
88 function fstr_ctrl_get_cflux( ctrl, amp, node_grp_name, node_grp_name_len, value )
89 implicit none
90 integer(kind=kint) :: ctrl
91 character(len=HECMW_NAME_LEN) :: amp
92 character(len=HECMW_NAME_LEN),target :: node_grp_name(:)
93 character(len=HECMW_NAME_LEN),pointer:: node_grp_name_p
94 integer(kind=kint) :: node_grp_name_len
95 real(kind=kreal),pointer :: value(:)
96 integer(kind=kint) :: fstr_ctrl_get_cflux
97
98 character(len=HECMW_NAME_LEN) :: data_fmt,ss
99
101
102 ! JP-9
103 if( fstr_ctrl_get_param_ex( ctrl, 'AMP ', '# ', 0, 'S', amp )/= 0) return
104
105 write(ss,*) node_grp_name_len
106 write(data_fmt,'(a,a,a)') 'S',trim(adjustl(ss)),'r '
107
108 node_grp_name_p => node_grp_name(1)
109 fstr_ctrl_get_cflux = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, node_grp_name_p, value )
110
111 end function fstr_ctrl_get_cflux
112
114 function fstr_ctrl_get_dflux( ctrl, amp, elem_grp_name, elem_grp_name_len, load_type, value )
115 implicit none
116 integer(kind=kint) :: ctrl
117 character(len=HECMW_NAME_LEN) :: amp
118 character(len=HECMW_NAME_LEN),target :: elem_grp_name(:)
119 character(len=HECMW_NAME_LEN),pointer:: elem_grp_name_p
120 integer(kind=kint) :: elem_grp_name_len
121 integer(kind=kint),pointer :: load_type(:)
122 real(kind=kreal),pointer :: value(:)
123 integer(kind=kint) :: fstr_ctrl_get_dflux
124
125 integer(kind=kint), parameter :: type_name_size = 5
126 integer(kind=kint) :: i, n
127 character(len=HECMW_NAME_LEN) :: data_fmt,s1,s2
128 character(len=type_name_size),pointer :: type_name_list(:)
129 character(len=type_name_size),pointer :: type_name_list_p
130 integer(kind=kint) :: rcode
131 integer(kind=kint) :: lid = -1
132
134
135 ! JP-10
136 if( fstr_ctrl_get_param_ex( ctrl, 'AMP ', '# ', 0, 'S', amp )/= 0) return
137
138 write(s1,*) elem_grp_name_len
139 write(s2,*) type_name_size
140 write(data_fmt,'(a,a,a,a,a)') 'S',trim(adjustl(s1)),'S',trim(adjustl(s2)),'r '
141
143 allocate( type_name_list(n) )
144
145 elem_grp_name_p => elem_grp_name(1)
146 type_name_list_p => type_name_list(1)
147 rcode = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, elem_grp_name_p, type_name_list_p, value)
148
149 if( rcode /= 0 ) then
150 deallocate( type_name_list )
151 return
152 end if
153
154 do i = 1, n
155 lid = -1;
156 call pc_strupr( type_name_list(i) )
157 if( type_name_list(i)(1:2) == 'BF' ) then; lid = 0
158 else if( type_name_list(i)(1:2) == 'S0' ) then; lid = 1
159 else if( type_name_list(i)(1:2) == 'S1' ) then; lid = 1
160 else if( type_name_list(i)(1:2) == 'S2' ) then; lid = 2
161 else if( type_name_list(i)(1:2) == 'S3' ) then; lid = 3
162 else if( type_name_list(i)(1:2) == 'S4' ) then; lid = 4
163 else if( type_name_list(i)(1:2) == 'S5' ) then; lid = 5
164 else if( type_name_list(i)(1:2) == 'S6' ) then; lid = 6
165 end if
166 if( lid < 0 ) then
167 write(ilog,*) 'Error : !DFLUX : Load type ',type_name_list(i),' is unknown'
168 deallocate( type_name_list )
169 return
170 end if
171 load_type(i) = lid
172 end do
173
174 deallocate( type_name_list )
176 end function fstr_ctrl_get_dflux
177
178 !* ----------------------------------------------------------------------------------------------- *!
180 !* ----------------------------------------------------------------------------------------------- *!
181
182 function fstr_ctrl_get_sflux( ctrl, amp, surface_grp_name, surface_grp_name_len, value )
183 implicit none
184 integer(kind=kint) :: ctrl
185 character(len=HECMW_NAME_LEN) :: amp
186 character(len=HECMW_NAME_LEN),target :: surface_grp_name(:)
187 character(len=HECMW_NAME_LEN),pointer:: surface_grp_name_p
188 integer(kind=kint) :: surface_grp_name_len
189 real(kind=kreal),pointer :: value(:)
190 integer(kind=kint) :: fstr_ctrl_get_sflux
191
192 character(len=HECMW_NAME_LEN) :: data_fmt,ss
193
195
196 ! JP-11
197 if( fstr_ctrl_get_param_ex( ctrl, 'AMP ', '# ', 0, 'S', amp )/= 0) return
198
199 write(ss,*) surface_grp_name_len
200 write(data_fmt,'(a,a,a)') 'S',trim(adjustl(ss)),'r '
201
202 surface_grp_name_p => surface_grp_name(1)
203 fstr_ctrl_get_sflux = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, surface_grp_name_p, value )
204 end function fstr_ctrl_get_sflux
205
206
207 !* ----------------------------------------------------------------------------------------------- *!
209 !* ----------------------------------------------------------------------------------------------- *!
210
211 function fstr_ctrl_get_film( ctrl, amp1, amp2, elem_grp_name, elem_grp_name_len, load_type, value, sink)
212 implicit none
213 integer(kind=kint) :: ctrl
214 character(len=HECMW_NAME_LEN) :: amp1
215 character(len=HECMW_NAME_LEN) :: amp2
216 character(len=HECMW_NAME_LEN),target :: elem_grp_name(:)
217 character(len=HECMW_NAME_LEN),pointer:: elem_grp_name_p
218 integer(kind=kint) :: elem_grp_name_len
219 integer(kind=kint),pointer :: load_type(:)
220 real(kind=kreal),pointer :: value(:)
221 real(kind=kreal),pointer :: sink(:)
222 integer(kind=kint) :: fstr_ctrl_get_film
223
224 integer(kind=kint),parameter :: type_name_size = 5
225 integer(kind=kint) :: i, n
226 character(len=HECMW_NAME_LEN) :: data_fmt,s1,s2
227 character(len=type_name_size),pointer :: type_name_list(:)
228 character(len=type_name_size),pointer :: type_name_list_p
229 integer(kind=kint) :: lid
230 integer(kind=kint) :: rcode
231
233
234 ! JP-12
235 if( fstr_ctrl_get_param_ex( ctrl, 'AMP1 ', '# ', 0, 'S', amp1 )/= 0) return
236 if( fstr_ctrl_get_param_ex( ctrl, 'AMP2 ', '# ', 0, 'S', amp2 )/= 0) return
237
238 write(s1,*) elem_grp_name_len
239 write(s2,*) type_name_size
240 write(data_fmt,'(a,a,a,a,a)') 'S',trim(adjustl(s1)),'S',trim(adjustl(s2)),'Rr '
241
243 allocate( type_name_list(n) )
244
245 elem_grp_name_p => elem_grp_name(1)
246 type_name_list_p => type_name_list(1)
247 rcode = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, elem_grp_name_p, type_name_list_p, value, sink)
248
249 if( rcode /= 0 ) then
250 deallocate( type_name_list )
251 return
252 end if
253
254 do i = 1, n
255 lid = -1;
256 call pc_strupr( type_name_list(i) )
257 if( type_name_list(i)(1:2) == 'F0' ) then; lid = 1
258 else if( type_name_list(i)(1:2) == 'F1' ) then; lid = 1
259 else if( type_name_list(i)(1:2) == 'F2' ) then; lid = 2
260 else if( type_name_list(i)(1:2) == 'F3' ) then; lid = 3
261 else if( type_name_list(i)(1:2) == 'F4' ) then; lid = 4
262 else if( type_name_list(i)(1:2) == 'F5' ) then; lid = 5
263 else if( type_name_list(i)(1:2) == 'F6' ) then; lid = 6
264 end if
265 if( lid < 0 ) then
266 write(ilog,*) 'Error : !FILM : Load type ',type_name_list(i),' is unknown'
267 deallocate( type_name_list )
268 return
269 end if
270 load_type(i) = lid
271 end do
272
273 deallocate( type_name_list )
274
276 end function fstr_ctrl_get_film
277
278 !* ----------------------------------------------------------------------------------------------- *!
280 !* ----------------------------------------------------------------------------------------------- *!
281
282 function fstr_ctrl_get_sfilm( ctrl, amp1, amp2, surface_grp_name, surface_grp_name_len, value, sink)
283 implicit none
284 integer(kind=kint) :: ctrl
285 character(len=HECMW_NAME_LEN) :: amp1
286 character(len=HECMW_NAME_LEN) :: amp2
287 character(len=HECMW_NAME_LEN),target :: surface_grp_name(:)
288 character(len=HECMW_NAME_LEN),pointer:: surface_grp_name_p
289 integer(kind=kint) :: surface_grp_name_len
290 real(kind=kreal),pointer :: value(:)
291 real(kind=kreal),pointer :: sink(:)
292 integer(kind=kint) :: fstr_ctrl_get_sfilm
293
294 character(len=HECMW_NAME_LEN) :: data_fmt,ss
295
297
298 ! JP-13
299 if( fstr_ctrl_get_param_ex( ctrl, 'AMP1 ', '# ', 0, 'S', amp1 )/= 0) return
300 if( fstr_ctrl_get_param_ex( ctrl, 'AMP2 ', '# ', 0, 'S', amp2 )/= 0) return
301
302 write(ss,*) surface_grp_name_len
303 write(data_fmt,'(a,a,a)') 'S',trim(adjustl(ss)),'Rr '
304
305 surface_grp_name_p => surface_grp_name(1)
306 fstr_ctrl_get_sfilm = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, surface_grp_name_p, value, sink )
307 end function fstr_ctrl_get_sfilm
308
309
310 !* ----------------------------------------------------------------------------------------------- *!
312 !* ----------------------------------------------------------------------------------------------- *!
313
314 function fstr_ctrl_get_radiate( ctrl, amp1, amp2, elem_grp_name, elem_grp_name_len, load_type, value, sink)
315 implicit none
316 integer(kind=kint) :: ctrl
317 character(len=HECMW_NAME_LEN) :: amp1
318 character(len=HECMW_NAME_LEN) :: amp2
319 character(len=HECMW_NAME_LEN),target :: elem_grp_name(:)
320 character(len=HECMW_NAME_LEN),pointer:: elem_grp_name_p
321 integer(kind=kint) :: elem_grp_name_len
322 integer(kind=kint),pointer :: load_type(:)
323 real(kind=kreal),pointer :: value(:)
324 real(kind=kreal),pointer :: sink(:)
325 integer(kind=kint) :: fstr_ctrl_get_radiate
326
327 integer(kind=kint),parameter :: type_name_size = 5
328 integer(kind=kint) :: i, n
329 character(len=HECMW_NAME_LEN) :: data_fmt,s1,s2
330 character(len=type_name_size),pointer :: type_name_list(:)
331 character(len=type_name_size),pointer :: type_name_list_p
332 integer(kind=kint) :: lid
333 integer(kind=kint) :: rcode
334
336
337 ! JP-14
338 if( fstr_ctrl_get_param_ex( ctrl, 'AMP1 ', '# ', 0, 'S', amp1 )/= 0) return
339 if( fstr_ctrl_get_param_ex( ctrl, 'AMP2 ', '# ', 0, 'S', amp2 )/= 0) return
340
341 write(s1,*) elem_grp_name_len
342 write(s2,*) type_name_size
343 write(data_fmt,'(a,a,a,a,a)') 'S',trim(adjustl(s1)),'S',trim(adjustl(s2)),'Rr '
344
346 allocate( type_name_list(n) )
347
348 elem_grp_name_p => elem_grp_name(1)
349 type_name_list_p => type_name_list(1)
350 rcode = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, elem_grp_name_p ,type_name_list_p, value, sink)
351
352 if( rcode /= 0 ) then
353 deallocate( type_name_list )
354 return
355 end if
356
357 do i = 1, n
358 lid = -1;
359 call pc_strupr( type_name_list(i) )
360 if( type_name_list(i)(1:2) == 'R0' ) then; lid = 1
361 else if( type_name_list(i)(1:2) == 'R1' ) then; lid = 1
362 else if( type_name_list(i)(1:2) == 'R2' ) then; lid = 2
363 else if( type_name_list(i)(1:2) == 'R3' ) then; lid = 3
364 else if( type_name_list(i)(1:2) == 'R4' ) then; lid = 4
365 else if( type_name_list(i)(1:2) == 'R5' ) then; lid = 5
366 else if( type_name_list(i)(1:2) == 'R6' ) then; lid = 6
367 end if
368 if( lid < 0 ) then
369 write(ilog,*) 'Error : !RADIATE : Load type ',type_name_list(i),' is unknown'
370 deallocate( type_name_list )
371 return
372 end if
373 load_type(i) = lid
374 end do
375
376 deallocate( type_name_list )
378 end function fstr_ctrl_get_radiate
379
380
381 !* ----------------------------------------------------------------------------------------------- *!
383 !* ----------------------------------------------------------------------------------------------- *!
384
385 function fstr_ctrl_get_sradiate( ctrl, amp1, amp2, surface_grp_name, surface_grp_name_len, value, sink)
386 implicit none
387 integer(kind=kint) :: ctrl
388 character(len=HECMW_NAME_LEN) :: amp1
389 character(len=HECMW_NAME_LEN) :: amp2
390 character(len=HECMW_NAME_LEN),target :: surface_grp_name(:)
391 character(len=HECMW_NAME_LEN),pointer:: surface_grp_name_p
392 integer(kind=kint) :: surface_grp_name_len
393 real(kind=kreal),pointer :: value(:)
394 real(kind=kreal),pointer :: sink(:)
395 integer(kind=kint) :: fstr_ctrl_get_sradiate
396
397 character(len=HECMW_NAME_LEN) :: data_fmt
398 character(len=HECMW_NAME_LEN) :: s1
399
401
402 ! JP-15
403 if( fstr_ctrl_get_param_ex( ctrl, 'AMP1 ', '# ', 0, 'S', amp1 )/= 0) return
404 if( fstr_ctrl_get_param_ex( ctrl, 'AMP2 ', '# ', 0, 'S', amp2 )/= 0) return
405
406 write(s1,*) surface_grp_name_len;
407 write(data_fmt,'(a,a,a)') 'S', trim(adjustl(s1)), 'Rr '
408
409 surface_grp_name_p => surface_grp_name(1)
410 fstr_ctrl_get_sradiate = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, surface_grp_name_p, value, sink )
411
412 end function fstr_ctrl_get_sradiate
413
414 !* ----------------------------------------------------------------------------------------------- *!
416 !* ----------------------------------------------------------------------------------------------- *!
417
418 function fstr_ctrl_get_weldline( ctrl, hecMESH, grp_name_len, weldline )
420 implicit none
421 integer(kind=kint), intent(in) :: ctrl
422 type(hecmwst_local_mesh), intent(in) :: hecmesh
423 integer(kind=kint), intent(in) :: grp_name_len
424 type(tweldline), intent(inout) :: weldline
425 integer(kind=kint) :: fstr_ctrl_get_weldline
426
427 character(len=HECMW_NAME_LEN) :: data_fmt
428 character(len=HECMW_NAME_LEN) :: s1, grp_id_name(1)
429 integer :: grp_id(1)
430
432 if( fstr_ctrl_get_data_ex( ctrl, 1, 'RRRR ', weldline%I, weldline%U, weldline%coe, weldline%v )/=0 ) return
433 write(s1,*) grp_name_len
434 write(data_fmt,'(a,a,a)') 'S', trim(adjustl(s1)), 'IRRRR '
435 if( fstr_ctrl_get_data_ex( ctrl, 2, data_fmt, grp_id_name, weldline%xyz, weldline%n1, &
436 weldline%n2, weldline%distol, weldline%tstart )/=0 ) return
437 call elem_grp_name_to_id( hecmesh, 'WELD_LINE ', 1, grp_id_name, grp_id )
438 weldline%egrpid = grp_id(1)
439
441 end function fstr_ctrl_get_weldline
442
443 !* ----------------------------------------------------------------------------------------------- *!
444end module fstr_ctrl_heat
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,...)
This module contains control file data obtaining functions for heat conductive analysis.
integer(kind=kint) function fstr_ctrl_get_heat(ctrl, dt, etime, dtmin, deltmx, itmax, eps, tpname)
Read in !HEAT.
integer(kind=kint) function fstr_ctrl_get_dflux(ctrl, amp, elem_grp_name, elem_grp_name_len, load_type, value)
Read in !DFLUX (heat)
integer(kind=kint) function fstr_ctrl_get_sflux(ctrl, amp, surface_grp_name, surface_grp_name_len, value)
Read in !SFLUX (heat)
integer(kind=kint) function fstr_ctrl_get_weldline(ctrl, hecmesh, grp_name_len, weldline)
Read in !WELD_LINE (heat)
integer(kind=kint) function fstr_ctrl_get_film(ctrl, amp1, amp2, elem_grp_name, elem_grp_name_len, load_type, value, sink)
Read in !FILM (heat)
integer(kind=kint) function fstr_ctrl_get_radiate(ctrl, amp1, amp2, elem_grp_name, elem_grp_name_len, load_type, value, sink)
Read in !RADIATE (heat)
integer(kind=kint) function fstr_ctrl_get_cflux(ctrl, amp, node_grp_name, node_grp_name_len, value)
Read in !CFLUX (heat)
integer(kind=kint) function fstr_ctrl_get_fixtemp(ctrl, amp, node_grp_name, node_grp_name_len, value)
Read in !FIXTEMP.
integer(kind=kint) function fstr_ctrl_get_sfilm(ctrl, amp1, amp2, surface_grp_name, surface_grp_name_len, value, sink)
Read in !SFILM (heat)
integer(kind=kint) function fstr_ctrl_get_sradiate(ctrl, amp1, amp2, surface_grp_name, surface_grp_name_len, value, sink)
Read in !SRADIATE (heat)
This module contains auxiliary functions in calculation setup.
subroutine elem_grp_name_to_id(hecmesh, header_name, n, grp_id_name, grp_id)
Definition: hecmw.f90:6
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
real(kind=kreal) eps
Definition: m_fstr.f90:126
real(kind=kreal) etime
Definition: m_fstr.f90:124
integer(kind=kint) itmax
Definition: m_fstr.f90:125
integer(kind=kint), parameter ilog
FILE HANDLER.
Definition: m_fstr.f90:91
real(kind=kreal) dt
ANALYSIS CONTROL for NLGEOM and HEAT.
Definition: m_fstr.f90:123
Data for weld line.
Definition: m_fstr.f90:597