OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
radiation.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| radiation ../engine/source/constraints/thermic/radiation.f
25!||--- called by ------------------------------------------------------
26!|| resol ../engine/source/engine/resol.F
27!||--- calls -----------------------------------------------------
28!|| finter ../engine/source/tools/curve/finter.F
29!|| omp_get_thread_num ../engine/source/engine/openmp_stub.F90
30!||--- uses -----------------------------------------------------
31!|| glob_therm_mod ../common_source/modules/mat_elem/glob_therm_mod.F90
32!|| python_funct_mod ../common_source/modules/python_mod.F90
33!|| sensor_mod ../common_source/modules/sensor_mod.F90
34!||====================================================================
35 SUBROUTINE radiation(IBCR ,FRADIA ,NPC ,TF , X ,
36 1 TEMP ,NSENSOR ,SENSOR_TAB,FTHE ,IAD ,
37 2 FTHESKY, PYTHON ,GLOB_THERM)
38C-----------------------------------------------
39C M o d u l e s
40C-----------------------------------------------
41 USE python_funct_mod
42 USE sensor_mod
43 use glob_therm_mod
44C-----------------------------------------------
45C I m p l i c i t T y p e s
46C-----------------------------------------------
47#include "implicit_f.inc"
48#include "param_c.inc"
49C-----------------------------------------------
50C C o m m o n B l o c k s
51C-----------------------------------------------
52#include "com04_c.inc"
53#include "com08_c.inc"
54#include "parit_c.inc"
55C-----------------------------------------------,
56C D u m m y A r g u m e n t s
57C-----------------------------------------------
58 type (glob_therm_) ,intent(inout) :: glob_therm
59 INTEGER ,INTENT(IN) :: NSENSOR
60 INTEGER NPC(*),IAD(4,*)
61 INTEGER IBCR(GLOB_THERM%NIRADIA,*)
62C REAL
63 my_real
64 . fradia(glob_therm%LFACTHER,*), tf(*), x(3,*),
65 . fthesky(lsky), temp(*), fthe(*)
66 TYPE (SENSOR_STR_) ,DIMENSION(NSENSOR) ,INTENT(IN) :: SENSOR_TAB
67 TYPE(PYTHON_) :: PYTHON
68C-----------------------------------------------
69C L o c a l V a r i a b l e s
70C-----------------------------------------------
71 INTEGER NL, N1, N2, N3, N4, ISENS,
72 . IFUNC_OLD,IFUNC,IAD1,IAD2,IAD3,IAD4
73 INTEGER IFLAG
74 my_real :: NX, NY, NZ, DYDX, TS, FLUX,
75 . ts_old, fcx, fcy, t_inf, te, area,
76 . startt, stopt, fcy_old, emisig, offg
77 my_real finter
78 my_real :: heat_radia_l ! thread-local
79 EXTERNAL finter
80 INTEGER :: OMP_GET_THREAD_NUM,ITSK
81 EXTERNAL omp_get_thread_num
82 INTEGER S1
83 INTEGER :: ISMOOTH
84
85C=======================================================================
86 ismooth = 0
87 s1 = numnod
88 ifunc_old = 0
89 ts_old = zero
90 fcy_old= zero
91 heat_radia_l = zero
92 t_inf = zero
93 n1 = 0
94 n2 = 0
95 n3 = 0
96 n4 = 0
97 IF(iparit==0) THEN
98 itsk = omp_get_thread_num()
99C-----------------------------------------------------------
100C CODE SPMD PARITH/OFF OU SMP NE PAS OUBLIER LE CODE P/ON !
101C-----------------------------------------------------------
102!$OMP DO SCHEDULE(GUIDED)
103 DO nl=1,glob_therm%NUMRADIA
104 offg = fradia(6,nl)
105 IF(offg <= zero) cycle
106 startt = fradia(4,nl)
107 stopt = fradia(5,nl)
108 isens = ibcr(6,nl)
109 IF(isens == 0)THEN
110 ts = tt*glob_therm%THEACCFACT - startt
111 ELSE
112 startt = startt + sensor_tab(isens)%TSTART
113 stopt = stopt + sensor_tab(isens)%TSTART
114 ts = tt*glob_therm%THEACCFACT -(startt + sensor_tab(isens)%TSTART)
115 ENDIF
116C
117 IF(tt*glob_therm%THEACCFACT < startt) cycle
118 IF(tt*glob_therm%THEACCFACT > stopt ) cycle
119 n1 = ibcr(1,nl)
120 n2 = ibcr(2,nl)
121 n3 = ibcr(3,nl)
122 n4 = ibcr(4,nl)
123 ifunc = ibcr(5,nl)
124 fcy = fradia(1,nl)
125 fcx = fradia(2,nl)
126 emisig= fradia(3,nl)
127C----------------------
128C RADIATION FLUX
129C----------------------
130 IF(ifunc_old /= ifunc .OR. ts_old /= ts .OR. fcy_old /= fcy ) THEN
131 ismooth = 0
132 IF (ifunc > 0) ismooth = npc(2*nfunct+ifunc+1)
133 IF(ismooth < 0) THEN
134 CALL python_call_funct1d(python, -ismooth,ts*fcx, t_inf)
135 t_inf = fcy*t_inf
136 ELSE
137 t_inf = fcy*finter(ifunc, ts*fcx,npc,tf,dydx)
138 ENDIF
139 ifunc_old = ifunc
140 ts_old = ts
141 fcy_old= fcy
142 ENDIF
143C ANALYSE 3D
144 IF(n4 > 0)THEN
145C
146 nx= (x(2,n3)-x(2,n1))*(x(3,n4)-x(3,n2))
147 + -(x(3,n3)-x(3,n1))*(x(2,n4)-x(2,n2))
148 ny= (x(3,n3)-x(3,n1))*(x(1,n4)-x(1,n2))
149 + -(x(1,n3)-x(1,n1))*(x(3,n4)-x(3,n2))
150 nz= (x(1,n3)-x(1,n1))*(x(2,n4)-x(2,n2))
151 + -(x(2,n3)-x(2,n1))*(x(1,n4)-x(1,n2))
152C
153 te = fourth*(temp(n1) + temp(n2) + temp(n3) + temp(n4))
154 area = half*sqrt( nx*nx + ny*ny + nz*nz)
155 flux = area*emisig*(t_inf**4 - te**4)*dt1*glob_therm%THEACCFACT
156 heat_radia_l = heat_radia_l + flux
157 flux = fourth*flux
158C
159 fthe(s1*itsk+n1) = fthe(s1*itsk+n1) + flux
160 fthe(s1*itsk+n2) = fthe(s1*itsk+n2) + flux
161 fthe(s1*itsk+n3)= fthe(s1*itsk+n3) + flux
162 fthe(s1*itsk+n4)= fthe(s1*itsk+n4) + flux
163
164C
165 ELSEIF(n3 > 0) THEN !TRUE TRIANGLES
166C
167 nx= (x(2,n3)-x(2,n1))*(x(3,n3)-x(3,n2))
168 + -(x(3,n3)-x(3,n1))*(x(2,n3)-x(2,n2))
169 ny= (x(3,n3)-x(3,n1))*(x(1,n3)-x(1,n2))
170 + -(x(1,n3)-x(1,n1))*(x(3,n3)-x(3,n2))
171 nz= (x(1,n3)-x(1,n1))*(x(2,n3)-x(2,n2))
172 + -(x(2,n3)-x(2,n1))*(x(1,n3)-x(1,n2))
173C
174 te = third*(temp(n1) + temp(n2) + temp(n3) )
175 area = half*sqrt( nx*nx + ny*ny + nz*nz)
176 flux = area*emisig*(t_inf**4 - te**4)*dt1*glob_therm%THEACCFACT
177 heat_radia_l = heat_radia_l + flux
178 flux = third*flux
179C
180 fthe(s1*itsk+n1) = fthe(s1*itsk+n1) + flux
181 fthe(s1*itsk+n2) = fthe(s1*itsk+n2) + flux
182 fthe(s1*itsk+n3)= fthe(s1*itsk+n3) + flux
183 ELSE
184C ANALYSE 2D
185 ny= -x(3,n2)+x(3,n1)
186 nz= x(2,n2)-x(2,n1)
187C
188 te = half*(temp(n1) + temp(n2) )
189 area = sqrt( ny*ny + nz*nz)
190 flux = area*emisig*(t_inf**4 - te**4)*dt1*glob_therm%THEACCFACT
191 glob_therm%HEAT_RADIA = glob_therm%HEAT_RADIA + flux
192 flux = half*flux
193C
194 fthe(s1*itsk+n1)=fthe(s1*itsk+n1) + flux
195 fthe(s1*itsk+n2)=fthe(s1*itsk+n2) + flux
196C
197 ENDIF
198 ENDDO
199!$OMP END DO
200
201!$OMP CRITICAL
202 glob_therm%HEAT_RADIA = glob_therm%HEAT_RADIA + heat_radia_l
203!$OMP END CRITICAL
204C
205 ELSE
206C-------------------------
207C CODE SPMD PARITH/ON
208C CODE NON VECTORIEL
209C-------------------------
210!$OMP DO SCHEDULE(GUIDED)
211 DO nl=1,glob_therm%NUMRADIA
212 startt = fradia(4,nl)
213 stopt = fradia(5,nl)
214 offg = fradia(6,nl)
215 isens = ibcr(6,nl)
216 IF(isens == 0)THEN
217 ts = tt*glob_therm%THEACCFACT - startt
218 ELSE
219 startt = startt + sensor_tab(isens)%TSTART
220 stopt = stopt + sensor_tab(isens)%TSTART
221 ts = tt*glob_therm%THEACCFACT -(startt + sensor_tab(isens)%TSTART)
222 ENDIF
223 iflag = 1
224 IF(tt*glob_therm%THEACCFACT < startt) iflag = 0
225 IF(tt*glob_therm%THEACCFACT > stopt ) iflag = 0
226 IF(offg <= zero) iflag = 0
227C---------------------
228C RADIATION FLUX
229C---------------------
230 IF(iflag==1) THEN
231 n1 =ibcr(1,nl)
232 n2 =ibcr(2,nl)
233 n3 =ibcr(3,nl)
234 n4 =ibcr(4,nl)
235 ifunc = ibcr(5,nl)
236 fcy = fradia(1,nl)
237 fcx = fradia(2,nl)
238 emisig= fradia(3,nl)
239 IF(ifunc_old /= ifunc .OR. ts_old /= ts) THEN
240 IF (ifunc > 0) ismooth = npc(2*nfunct+ifunc+1)
241 IF(ismooth < 0 .AND. ifunc > 0) THEN
242 CALL python_call_funct1d(python, -ismooth,ts*fcx, t_inf)
243 ELSE
244 t_inf = finter(ifunc,ts*fcx,npc,tf,dydx)
245 ENDIF
246 ifunc_old = ifunc
247 ts_old = ts
248 ENDIF
249 t_inf = fcy*t_inf
250C ANALYSE 3D
251 IF(n4 > 0)THEN
252 nx= (x(2,n3)-x(2,n1))*(x(3,n4)-x(3,n2))
253 + -(x(3,n3)-x(3,n1))*(x(2,n4)-x(2,n2))
254 ny= (x(3,n3)-x(3,n1))*(x(1,n4)-x(1,n2))
255 + -(x(1,n3)-x(1,n1))*(x(3,n4)-x(3,n2))
256 nz= (x(1,n3)-x(1,n1))*(x(2,n4)-x(2,n2))
257 + -(x(2,n3)-x(2,n1))*(x(1,n4)-x(1,n2))
258 te = fourth*(temp(n1) + temp(n2) + temp(n3) + temp(n4))
259 area = half*sqrt( nx*nx + ny*ny + nz*nz)
260 flux = area*emisig*(t_inf**4 - te**4)*dt1*glob_therm%THEACCFACT
261 heat_radia_l = heat_radia_l + flux
262 flux = fourth*flux
263C
264 iad1 = iad(1,nl)
265 fthesky(iad1) = flux
266 iad2 = iad(2,nl)
267 fthesky(iad2) = flux
268 iad3 = iad(3,nl)
269 fthesky(iad3) = flux
270 iad4 = iad(4,nl)
271 fthesky(iad4) = flux
272C
273 ELSEIF( n3 > 0) THEN
274C TRUE TRIANGLES
275 nx= (x(2,n3)-x(2,n1))*(x(3,n3)-x(3,n2))
276 + -(x(3,n3)-x(3,n1))*(x(2,n3)-x(2,n2))
277 ny= (x(3,n3)-x(3,n1))*(x(1,n3)-x(1,n2))
278 + -(x(1,n3)-x(1,n1))*(x(3,n3)-x(3,n2))
279 nz= (x(1,n3)-x(1,n1))*(x(2,n3)-x(2,n2))
280 + -(x(2,n3)-x(2,n1))*(x(1,n3)-x(1,n2))
281C
282 te = third*(temp(n1) + temp(n2) + temp(n3) )
283 area = half*sqrt( nx*nx + ny*ny + nz*nz)
284 flux = area*emisig*(t_inf**4 - te**4)*dt1*glob_therm%THEACCFACT
285 heat_radia_l = heat_radia_l + flux
286 flux = third*flux
287C
288 iad1 = iad(1,nl)
289 fthesky(iad1) = flux
290C
291 iad2 = iad(2,nl)
292 fthesky(iad2) = flux
293C
294 iad3 = iad(3,nl)
295 fthesky(iad3) = flux
296C
297 ELSE
298C ANALYSE 2D
299 ny= -x(3,n2)+x(3,n1)
300 nz= x(2,n2)-x(2,n1)
301
302C
303 te = half*(temp(n1) + temp(n2) )
304 area = sqrt( ny*ny + nz*nz)
305 flux = area*emisig*(t_inf**4 - te**4)*dt1*glob_therm%THEACCFACT
306 heat_radia_l = heat_radia_l + flux
307 flux = half*flux
308C
309 iad1 = iad(1,nl)
310 fthesky(iad1) = flux
311C
312 iad2 = iad(2,nl)
313 fthesky(iad2) = flux
314 ENDIF
315C
316 ELSE ! IFLAG=0
317 iad1 = iad(1,nl)
318 fthesky(iad1) = zero
319 iad2 = iad(2,nl)
320 fthesky(iad2) = zero
321 IF(n4 > 0)THEN
322 iad3 = iad(3,nl)
323 fthesky(iad3) = zero
324 iad4 = iad(4,nl)
325 fthesky(iad4) = zero
326 ELSEIF(n3 > 0)THEN
327 iad3 = iad(3,nl)
328 fthesky(iad3) = zero
329 ENDIF
330 ENDIF
331 ENDDO
332!$OMP END DO
333
334!$OMP CRITICAL
335 glob_therm%HEAT_RADIA = glob_therm%HEAT_RADIA + heat_radia_l
336!$OMP END CRITICAL
337
338C
339 ENDIF
340C
341 RETURN
342 END
my_real
#define my_real
Definition cppsort.cpp:32
area
subroutine area(d1, x, x2, y, y2, eint, stif0)
Definition hm_read_prop32.F:567
radiation
subroutine radiation(ibcr, fradia, npc, tf, x, temp, nsensor, sensor_tab, fthe, iad, fthesky, python, glob_therm)
Definition radiation.F:38