FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
heat_init.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!-------------------------------------------------------------------------------
7contains
8
9 subroutine heat_init(hecMESH, fstrHEAT)
10 use m_fstr
11 implicit none
12 type(fstr_heat) :: fstrHEAT
13 type(hecmwst_local_mesh) :: hecMESH
14 integer(kind=kint) :: i, j
15
16 allocate(fstrheat%TEMP0(hecmesh%n_node))
17 allocate(fstrheat%TEMPC(hecmesh%n_node))
18 allocate(fstrheat%TEMP (hecmesh%n_node))
19 fstrheat%TEMP0 = 0.0d0
20 fstrheat%TEMPC = 0.0d0
21 fstrheat%TEMP = 0.0d0
22
23 if(hecmesh%hecmw_flag_initcon == 1)then
24 do i = 1, hecmesh%n_node
25 j = hecmesh%node_init_val_index(i)
26 fstrheat%TEMP0(i) = hecmesh%node_init_val_item(j)
27 fstrheat%TEMPC(i) = fstrheat%TEMP0(i)
28 fstrheat%TEMP (i) = fstrheat%TEMP0(i)
29 enddo
30 write(ilog,*) ' Initial condition of temperatures: OK'
31 endif
32 if( associated(g_initialcnd) ) then
33 do j=1,size(g_initialcnd)
34 if( g_initialcnd(j)%cond_name=="temperature" ) then
35 do i= 1, hecmesh%n_node
36 fstrheat%TEMP0(i)= g_initialcnd(j)%realval(i)
37 fstrheat%TEMPC(i)= fstrheat%TEMP0(i)
38 fstrheat%TEMP (i)= fstrheat%TEMP0(i)
39 enddo
40 exit
41 endif
42 enddo
43 write(ilog,*) ' Initial condition of temperatures: OK'
44 endif
45 end subroutine heat_init
46
47 subroutine heat_init_log(hecMESH)
48 use m_fstr
49 implicit none
50 type(hecmwst_local_mesh) :: hecMESH
51
52 if(hecmesh%my_rank == 0)then
53 write(imsg,*) '============================='
54 write(imsg,*) ' H E A T T R A N S F E R '
55 write(imsg,*) '============================='
56 write(ista,*)
57 write(ista,*)' ISTEP INCR ITER RESIDUAL IITER '
58 write(ista,*)'-------------------------------------------------'
59 endif
60 end subroutine heat_init_log
61
62 subroutine heat_finalize(fstrHEAT)
63 use m_fstr
64 implicit none
65 type(fstr_heat) :: fstrHEAT
66 deallocate(fstrheat%TEMP0)
67 deallocate(fstrheat%TEMPC)
68 deallocate(fstrheat%TEMP )
69 end subroutine heat_finalize
70
71 !C***
72 !C*** INIT_AMPLITUDE
73 !C***
74 subroutine heat_init_amplitude (hecMESH, fstrHEAT)
75
76 use m_fstr
77
78 implicit none
79 integer(kind=kint) :: namax, i, nn, is, iE, icou, j, k
80 real(kind=kreal) :: x1, y1, x2, y2
81 type(fstr_heat) :: fstrheat
82 type(hecmwst_local_mesh) :: hecMESH
83 !C
84 !C===
85 namax = 0
86 do i = 1, hecmesh%amp%n_amp
87 nn = hecmesh%amp%amp_index(i) - hecmesh%amp%amp_index(i-1)
88 namax = max(nn,namax)
89 enddo
90
91 fstrheat%AMPLITUDEtot= hecmesh%amp%n_amp
92
93 allocate (fstrheat%AMPLtab (fstrheat%AMPLITUDEtot) )
94 allocate (fstrheat%AMPL (fstrheat%AMPLITUDEtot,namax), &
95 fstrheat%AMPLtime(fstrheat%AMPLITUDEtot,namax) )
96
97 fstrheat%AMPLtab = 0
98 fstrheat%AMPL = 0.d0
99 fstrheat%AMPLtime = 0.d0
100
101 do i = 1, fstrheat%AMPLITUDEtot
102 is = hecmesh%amp%amp_index(i-1) + 1
103 ie = hecmesh%amp%amp_index(i)
104 nn = ie - is + 1
105
106 fstrheat%AMPLtab(i) = nn
107
108 icou = 0
109 do j = is, ie
110 icou = icou + 1
111 fstrheat%AMPL (i,icou) = hecmesh%amp%amp_val (j)
112 fstrheat%AMPLtime(i,icou) = hecmesh%amp%amp_table(j)
113 enddo
114 enddo
115 !C===
116
117 !C
118 !C +-----------+
119 !C | AMP-TABLE |
120 !C +-----------+
121 !C===
122 allocate ( fstrheat%AMPLfuncA( fstrheat%AMPLITUDEtot,namax+1 ) )
123 allocate ( fstrheat%AMPLfuncB( fstrheat%AMPLITUDEtot,namax+1 ) )
124
125 fstrheat%AMPLfuncA = 0.d0
126 fstrheat%AMPLfuncB = 0.d0
127 !C
128 !C--
129 do i = 1, fstrheat%AMPLITUDEtot
130
131 fstrheat%AMPLfuncA(i,1) = 0.d0
132 fstrheat%AMPLfuncB(i,1) = fstrheat%AMPL(i,1)
133
134 nn = fstrheat%AMPLtab(i)
135
136 do k = 2, nn
137
138 x1 = fstrheat%AMPLtime(i,k-1)
139 y1 = fstrheat%AMPL (i,k-1)
140 x2 = fstrheat%AMPLtime(i,k)
141 y2 = fstrheat%AMPL (i,k)
142
143 fstrheat%AMPLfuncA(i,k) = (y2-y1)/(x2-x1)
144 fstrheat%AMPLfuncB(i,k) = -(y2-y1)/(x2-x1)*x1 + y1
145
146 enddo
147
148 fstrheat%AMPLfuncA(i,nn+1) = 0.d0
149 fstrheat%AMPLfuncB(i,nn+1) = fstrheat%AMPL(i,nn)
150
151 enddo
152
153 !C===
154 end subroutine heat_init_amplitude
155 !C***
156 !C*** INIT_MATERIAL
157 !C***
158 subroutine heat_init_material (hecMESH, fstrHEAT)
159
160 use m_fstr
161
162 implicit none
163 integer(kind=kint) :: m1max, m2max, m3max, icou, im, jm, nn, ic, jS, jE, kc, km, k
164 real(kind=kreal) :: aa, bb
165 type(fstr_heat) :: fstrheat
166 type(hecmwst_local_mesh) :: hecMESH
167
168 !C
169 !C +----------+
170 !C | MATERIAL |
171 !C +----------+
172 !C===
173 fstrheat%MATERIALtot= hecmesh%material%n_mat
174
175 m1max = 0
176 m2max = 0
177 m3max = 0
178 icou = 0
179 do im = 1, hecmesh%material%n_mat
180 do jm = 1, 3
181 icou = icou + 1
182 nn = hecmesh%material%mat_TABLE_index(icou) - hecmesh%material%mat_TABLE_index(icou-1)
183 if( jm.eq.1 ) m1max = max(nn,m1max)
184 if( jm.eq.2 ) m2max = max(nn,m2max)
185 if( jm.eq.3 ) m3max = max(nn,m3max)
186 enddo
187 enddo
188
189 allocate (fstrheat%RHOtab (fstrheat%MATERIALtot), &
190 fstrheat%CPtab (fstrheat%MATERIALtot), &
191 fstrheat%CONDtab (fstrheat%MATERIALtot))
192 allocate (fstrheat%RHO (fstrheat%MATERIALtot,m1max), &
193 fstrheat%RHOtemp (fstrheat%MATERIALtot,m1max))
194 allocate (fstrheat%CP (fstrheat%MATERIALtot,m2max), &
195 fstrheat%CPtemp (fstrheat%MATERIALtot,m2max))
196 allocate (fstrheat%COND (fstrheat%MATERIALtot,m3max), &
197 fstrheat%CONDtemp(fstrheat%MATERIALtot,m3max))
198
199 fstrheat%RHO = 0.d0
200 fstrheat%CP = 0.d0
201 fstrheat%COND = 0.d0
202
203 fstrheat%RHOtemp = 0.d0
204 fstrheat%CPtemp = 0.d0
205 fstrheat%CONDtemp = 0.d0
206 fstrheat%RHOtab = 0
207 fstrheat%CPtab = 0
208 fstrheat%CONDtab = 0
209
210 ic = 0
211 do im = 1, fstrheat%MATERIALtot
212 do jm = 1, 3
213 ic = ic + 1
214 js = hecmesh%material%mat_TABLE_index(ic-1) + 1
215 je = hecmesh%material%mat_TABLE_index(ic )
216 nn = je - js + 1
217 if( jm.eq.1 ) fstrheat%RHOtab (im) = nn
218 if( jm.eq.2 ) fstrheat%CPtab (im) = nn
219 if( jm.eq.3 ) fstrheat%CONDtab(im) = nn
220
221 kc = 0
222 do km = js, je
223 kc = kc + 1
224 if( jm.eq.1 ) then
225 fstrheat%RHO (im,kc) = hecmesh%material%mat_VAL (km)
226 fstrheat%RHOtemp (im,kc) = hecmesh%material%mat_TEMP(km)
227 endif
228 if( jm.eq.2 ) then
229 fstrheat%CP (im,kc) = hecmesh%material%mat_VAL (km)
230 fstrheat%CPtemp (im,kc) = hecmesh%material%mat_TEMP(km)
231 endif
232 if( jm.eq.3 ) then
233 fstrheat%COND (im,kc) = hecmesh%material%mat_VAL (km)
234 fstrheat%CONDtemp(im,kc) = hecmesh%material%mat_TEMP(km)
235 endif
236 enddo
237 enddo
238 enddo
239 !C===
240
241 !C
242 !C +-----------+
243 !C | MAT-TABLE |
244 !C +-----------+
245 !C===
246 allocate (fstrheat%RHOfuncA (fstrheat%MATERIALtot, m1max+1) &
247 ,fstrheat%RHOfuncB (fstrheat%MATERIALtot, m1max+1))
248 allocate (fstrheat%CPfuncA (fstrheat%MATERIALtot, m2max+1) &
249 ,fstrheat%CPfuncB (fstrheat%MATERIALtot, m2max+1))
250 allocate (fstrheat%CONDfuncA(fstrheat%MATERIALtot, m3max+1) &
251 ,fstrheat%CONDfuncB(fstrheat%MATERIALtot, m3max+1))
252
253 fstrheat%RHOfuncA = 0.d0
254 fstrheat%RHOfuncB = 0.d0
255 fstrheat%CPfuncA = 0.d0
256 fstrheat%CPfuncB = 0.d0
257 fstrheat%CONDfuncA = 0.d0
258 fstrheat%CONDfuncB = 0.d0
259 !C
260 !C--RHO
261 do im = 1, fstrheat%MATERIALtot
262 fstrheat%RHOfuncB(im,1) = fstrheat%RHO(im,1)
263 do k = 2, fstrheat%RHOtab(im)
264 bb= fstrheat%RHO (im,k) - fstrheat%RHO (im,k-1)
265 aa= fstrheat%RHOtemp(im,k) - fstrheat%RHOtemp(im,k-1)
266 fstrheat%RHOfuncA(im,k) = bb/aa
267 fstrheat%RHOfuncB(im,k) = -(bb/aa)*fstrheat%RHOtemp(im,k-1) + fstrheat%RHO(im,k-1)
268 enddo
269 fstrheat%RHOfuncB(im,fstrheat%RHOtab(im)+1) = fstrheat%RHO(im,fstrheat%RHOtab(im))
270 enddo
271 !C
272 !C-- CP
273 do im = 1, fstrheat%MATERIALtot
274 fstrheat%CPfuncB(im,1) = fstrheat%CP(im,1)
275 do k = 2, fstrheat%CPtab(im)
276 bb= fstrheat%CP (im,k) - fstrheat%CP (im,k-1)
277 aa= fstrheat%CPtemp(im,k) - fstrheat%CPtemp(im,k-1)
278 fstrheat%CPfuncA(im,k) = bb/aa
279 fstrheat%CPfuncB(im,k) = -(bb/aa)*fstrheat%CPtemp(im,k-1) + fstrheat%CP(im,k-1)
280 enddo
281 fstrheat%CPfuncB(im,fstrheat%CPtab(im)+1) = fstrheat%CP(im,fstrheat%CPtab(im))
282 enddo
283 !C
284 !C-- COND.
285 do im = 1, fstrheat%MATERIALtot
286 fstrheat%CONDfuncB(im,1)= fstrheat%COND(im,1)
287 do k = 2, fstrheat%CONDtab(im)
288 bb = fstrheat%COND (im,k) - fstrheat%COND (im,k-1)
289 aa = fstrheat%CONDtemp(im,k) - fstrheat%CONDtemp(im,k-1)
290 fstrheat%CONDfuncA(im,k) = bb/aa
291 fstrheat%CONDfuncB(im,k) = -(bb/aa)*fstrheat%CONDtemp(im,k-1) + fstrheat%COND(im,k-1)
292 enddo
293 fstrheat%CONDfuncB(im,fstrheat%CONDtab(im)+1) = fstrheat%COND(im,fstrheat%CONDtab(im))
294 enddo
295 !C===
296 end subroutine heat_init_material
297
298end module m_heat_init
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
integer(kind=kint), parameter ilog
FILE HANDLER.
Definition: m_fstr.f90:91
integer(kind=kint), parameter ista
Definition: m_fstr.f90:92
type(tinitialcondition), dimension(:), pointer, save g_initialcnd
Definition: m_fstr.f90:135
This module provides functions to initialize heat analysis.
Definition: heat_init.f90:6
subroutine heat_init_log(hecmesh)
Definition: heat_init.f90:48
subroutine heat_finalize(fstrheat)
Definition: heat_init.f90:63
subroutine heat_init_amplitude(hecmesh, fstrheat)
Definition: heat_init.f90:75
subroutine heat_init(hecmesh, fstrheat)
Definition: heat_init.f90:10
subroutine heat_init_material(hecmesh, fstrheat)
Definition: heat_init.f90:159
Data for HEAT ANSLYSIS (fstrHEAT)
Definition: m_fstr.f90:394