FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
heat_mat_ass_bc_RADIATE.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!-------------------------------------------------------------------------------
8contains
9 !C
10 !C***
11 !C*** MAT_ASS_RADIATE
12 !C***
13 !C
14 subroutine heat_mat_ass_bc_radiate( hecMESH, hecMAT, fstrHEAT, CTIME, DTIME, beta )
15
16 use m_fstr
19
20 implicit none
21 integer(kind=kint) :: k, icel, isuf, iam1, iam2, ic_type, isect, nn, is, j, mm, m, ic, ip
22 integer(kind=kint) :: inod, jp, jnod, isU, ieU, ik, isL, ieL
23 real(kind=kreal) :: ctime, dtime, tzero, qq, rr, sink, thick, beta
24 type(fstr_heat) :: fstrheat
25 type(hecmwst_matrix) :: hecMAT
26 type(hecmwst_local_mesh) :: hecMESH
27
28 real(kind=kreal) :: xx(20), yy(20), zz(20), tt(20)
29 real(kind=kreal) :: term1(64), term2(20), stiff(8,8)
30 integer(kind=kint) :: nodLocal(20), nsuf(8), nodSurf(8)
31
32 tzero = hecmesh%zero_temp
33 !C
34 do k = 1, fstrheat%R_SUF_tot
35 icel = fstrheat%R_SUF_elem(k)
36 isuf = fstrheat%R_SUF_surf(k)
37 iam1 = fstrheat%R_SUF_ampl(k,1)
38 iam2 = fstrheat%R_SUF_ampl(k,2)
39 call heat_get_amplitude ( fstrheat, iam1, ctime, qq )
40 rr = fstrheat%R_SUF_val (k,1) * qq
41 call heat_get_amplitude ( fstrheat, iam2, ctime, qq )
42 sink = fstrheat%R_SUF_val (k,2) * qq
43 ic_type = hecmesh%elem_type(icel)
44 isect = hecmesh%section_ID(icel)
45 !C**
46 nn = hecmw_get_max_node(ic_type)
47 !C**
48 is = hecmesh%elem_node_index(icel-1)
49 do j = 1, nn
50 nodlocal(j) = hecmesh%elem_node_item(is+j)
51 xx(j) = hecmesh%node( 3*nodlocal(j)-2 )
52 yy(j) = hecmesh%node( 3*nodlocal(j)-1 )
53 zz(j) = hecmesh%node( 3*nodlocal(j) )
54 tt(j) = fstrheat%TEMP( nodlocal(j) )
55 enddo
56 !C**
57 if ( ic_type.eq.231 ) then
58 is = hecmesh%section%sect_R_index(isect)
59 thick = hecmesh%section%sect_R_item(is)
60 mm=2
61 call heat_radiate_231(nn,xx,yy,zz,thick,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
62 elseif( ic_type.eq.232 ) then
63 is = hecmesh%section%sect_R_index(isect)
64 thick = hecmesh%section%sect_R_item(is)
65 mm=3
66 call heat_radiate_232(nn,xx,yy,zz,thick,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
67 elseif( ic_type.eq.241 ) then
68 is = hecmesh%section%sect_R_index(isect)
69 thick = hecmesh%section%sect_R_item(is)
70 mm=2
71 call heat_radiate_241(nn,xx,yy,zz,thick,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
72 elseif( ic_type.eq.242 ) then
73 is = hecmesh%section%sect_R_index(isect)
74 thick = hecmesh%section%sect_R_item(is)
75 mm=3
76 call heat_radiate_242(nn,xx,yy,zz,thick,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
77 elseif( ic_type.eq.341 ) then
78 mm=3
79 call heat_radiate_341(nn,xx,yy,zz,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
80 elseif( ic_type.eq.342 ) then
81 mm=6
82 call heat_radiate_342(nn,xx,yy,zz,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
83 elseif( ic_type.eq.351 ) then
84 mm=4
85 if( isuf.eq.1 .or. isuf.eq.2 ) mm=3
86 call heat_radiate_351(nn,xx,yy,zz,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
87 elseif( ic_type.eq.352 ) then
88 mm=8
89 if( isuf.eq.1 .or. isuf.eq.2 ) mm=6
90 call heat_radiate_352(nn,xx,yy,zz,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
91 elseif( ic_type.eq.361 ) then
92 mm=4
93 call heat_radiate_361(nn,xx,yy,zz,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
94 elseif( ic_type.eq.362 ) then
95 mm=8
96 call heat_radiate_362(nn,xx,yy,zz,tt,isuf,rr,sink,tzero,mm,term1,term2,nsuf)
97 elseif( ic_type.eq.731 ) then
98 call heat_radiate_731(nn,xx,yy,zz,tt,isuf,rr,sink,tzero,term1,term2)
99 mm=3
100 do m = 1, mm
101 nsuf(m) = m
102 enddo
103 elseif( ic_type.eq.741 ) then
104 call heat_radiate_741(nn,xx,yy,zz,tt,isuf,rr,sink,tzero,term1,term2)
105 mm=4
106 do m = 1, mm
107 nsuf(m) = m
108 enddo
109 endif
110 !C
111 do ip = 1, mm
112 nodsurf(ip) = nodlocal(nsuf(ip))
113 end do
114 !C
115 ic = 0
116 stiff = 0.d0
117 do ip = 1, mm
118 do jp = 1, mm
119 ic = ic + 1
120 stiff(ip,jp) = -term1(ic)
121 enddo
122 enddo
123 !C
124 call hecmw_mat_ass_elem(hecmat, mm, nodsurf, stiff)
125 !C
126 do ip = 1, mm
127 !$omp atomic
128 hecmat%B(nodsurf(ip)) = hecmat%B(nodsurf(ip)) - term2(ip)
129 end do
130 !C
131 enddo
132
133 end subroutine heat_mat_ass_bc_radiate
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
This moudle provide a function to get amplitude definition.
subroutine heat_get_amplitude(fstrheat, id, tt, qq, outofrange)
This module provides subroutines to generate heat radiate boundary.
subroutine heat_radiate_342(nn, xx, yy, zz, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_361(nn, xx, yy, zz, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_242(nn, xx, yy, zz, thick, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_352(nn, xx, yy, zz, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_231(nn, xx, yy, zz, thick, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_351(nn, xx, yy, zz, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_731(nn, xx, yy, zz, temp, ltype, rr, sink, tzero, term1, term2)
subroutine heat_radiate_341(nn, xx, yy, zz, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_741(nn, xx, yy, zz, temp, ltype, rr, sink, tzero, term1, term2)
subroutine heat_radiate_241(nn, xx, yy, zz, thick, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_362(nn, xx, yy, zz, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
subroutine heat_radiate_232(nn, xx, yy, zz, thick, temp, ltype, rr, sink, tzero, mm, term1, term2, nod)
This module provides a subroutine for setting heat radiate boundary conditions.
subroutine heat_mat_ass_bc_radiate(hecmesh, hecmat, fstrheat, ctime, dtime, beta)
Data for HEAT ANSLYSIS (fstrHEAT)
Definition: m_fstr.f90:394