FrontISTR  5.7.0
Large-scale structural analysis program with finit element method
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 !-------------------------------------------------------------------------------
6 module m_heat_init
7 contains
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 
298 end module m_heat_init
m_heat_init
This module provides functions to initialize heat analysis.
Definition: heat_init.f90:6
m_heat_init::heat_init_material
subroutine heat_init_material(hecMESH, fstrHEAT)
Definition: heat_init.f90:159
m_fstr::fstr_heat
Data for HEAT ANSLYSIS (fstrHEAT)
Definition: m_fstr.f90:425
m_fstr
This module defines common data and basic structures for analysis.
Definition: m_fstr.f90:15
m_fstr::g_initialcnd
type(tinitialcondition), dimension(:), pointer, save g_initialcnd
Definition: m_fstr.f90:151
m_fstr::ista
integer(kind=kint), parameter ista
Definition: m_fstr.f90:108
m_heat_init::heat_init
subroutine heat_init(hecMESH, fstrHEAT)
Definition: heat_init.f90:10
m_fstr::ilog
integer(kind=kint), parameter ilog
FILE HANDLER.
Definition: m_fstr.f90:107
m_heat_init::heat_init_log
subroutine heat_init_log(hecMESH)
Definition: heat_init.f90:48
m_heat_init::heat_finalize
subroutine heat_finalize(fstrHEAT)
Definition: heat_init.f90:63
m_heat_init::heat_init_amplitude
subroutine heat_init_amplitude(hecMESH, fstrHEAT)
Definition: heat_init.f90:75
m_fstr::imsg
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:110