OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
phase_detection.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!|| phase_detection ../starter/source/initial_conditions/inivol/phase_detection.F
25!||--- called by ------------------------------------------------------
26!|| getphase ../starter/source/initial_conditions/inivol/getphase.F
27!||--- calls -----------------------------------------------------
28!|| find_closest_node ../starter/source/initial_conditions/inivol/find_closest_node.F
29!|| in_out_side ../starter/source/initial_conditions/inivol/in_out_side.F
30!|| phase_propagation ../starter/source/initial_conditions/inivol/phase_propagation.F
31!||--- uses -----------------------------------------------------
32!||====================================================================
33 SUBROUTINE phase_detection(
34 . NPARG,NGROUP,NUMELS,NUMELQ,NUMELTG,NUMNOD,NSURF,N2D,
35 . LEADING_DIMENSION,NB_CELL_X,NB_CELL_Y,NB_CELL_Z,NB_BOX_LIMIT,
36 . IPARG,IXS,IXQ,IXTG,X,ID_SURFACE,
37 . CELL,CELL_POSITION,NODAL_PHASE,CLOSEST_NODE_ID,
38 . NNOD2SURF,KNOD2SURF,INOD2SURF,
39 . NOD_NORMAL,NSEG_USED,SEGTOSURF,NSEG,SURF_ELTYP,SURFACE_NODES,SWIFTSURF)
40!$COMMENT
41! PHASE_DETECTION description
42! PHASE_DETECTION find the pseudo distance node - surface
43! pseudo distance = 1 (node in the phase 1), -1 (node in the phase -1), or a real distance
44! there are 3 cases :
45! * a node is far from the surface and the pseudo distance of neighbour nodes is known
46! --> apply the pseudo distance of neighbour nodes to the node
47! * a node is far from the surface and the pseudo distance of neighbour nodes is unknown
48! --> compute the pseudo distance and apply it to the node
49! * a node is close to a surface --> need to compute the real distance to the surface
50!
51! PHASE_DETECTION organization :
52! - loop over the ALE element group
53! - for each ALE element group, loop over the NEL elements
54! - loop over the nodes of each element
55! * if all the nodes are far from the surface
56! * if the pseudo distance of at least 1 node is known --> apply it to the other nodes
57! * if the pseudo distance of all nodes is unknown --> need to compute the pseudo distance &
58! apply it to the other nodes
59! * if at least one node is near a surface --> need to compute the real distance
60!$ENDCOMMENT
61 use element_mod , only :nixs,nixq,nixtg
62C-----------------------------------------------
63C I m p l i c i t T y p e s
64C-----------------------------------------------
65#include "implicit_f.inc"
66C-----------------------------------------------
67C D u m m y A r g u m e n t s
68C-----------------------------------------------
69 INTEGER, INTENT(IN) :: NPARG,NGROUP,NUMELS,NUMELQ,NUMELTG,NUMNOD,NSURF,N2D
70 INTEGER, INTENT(IN) :: LEADING_DIMENSION
71 INTEGER, INTENT(IN) :: NB_BOX_LIMIT ! maximum number of cell
72 INTEGER, INTENT(IN) :: NB_CELL_X,NB_CELL_Y,NB_CELL_Z
73 INTEGER, DIMENSION(NPARG,NGROUP), INTENT(IN) :: IPARG ! group data
74 INTEGER, DIMENSION(NIXS,NUMELS),INTENT(IN), TARGET :: IXS ! solid data
75 INTEGER, DIMENSION(NIXQ,NUMELQ),INTENT(IN), TARGET :: IXQ ! quad data
76 INTEGER, DIMENSION(NIXTG,NUMELTG),INTENT(IN), TARGET :: IXTG ! triangle data
77 INTEGER, DIMENSION(NUMNOD), INTENT(INOUT) :: NODAL_PHASE ! phase of nodes (in / out / near the surface)
78 INTEGER, DIMENSION(NUMNOD), INTENT(INOUT) :: CLOSEST_NODE_ID ! list of closest node id
79 INTEGER, DIMENSION(3,NUMNOD), INTENT(IN) :: CELL_POSITION ! position of node/cell
80 INTEGER, DIMENSION(NB_CELL_X,NB_CELL_Y,NB_CELL_Z), INTENT(INOUT) :: CELL ! phase of the voxcell
81 my_real, DIMENSION(3,NUMNOD), INTENT(IN) :: x ! position
82 INTEGER, INTENT(IN) :: ID_SURFACE ! id of the surface
83 INTEGER, DIMENSION(NSEG,4), INTENT(IN) :: SURFACE_NODES ! list of nodes for each segment of the surface
84 INTEGER, INTENT(IN) :: NNOD2SURF,NSEG_USED ! size of SEGTOSURF & INOD2SURF arrays
85 INTEGER, DIMENSION(NUMNOD+1), INTENT(IN) :: KNOD2SURF
86 INTEGER, DIMENSION(NNOD2SURF,NUMNOD), INTENT(IN) :: INOD2SURF
87 my_real, DIMENSION(3,NUMNOD), INTENT(IN) :: nod_normal
88 INTEGER, DIMENSION(NSEG_USED), INTENT(IN) :: SEGTOSURF
89 INTEGER, INTENT(IN) :: NSEG ! number of segment of the current surface
90 INTEGER, DIMENSION(NSEG), INTENT(IN) :: SURF_ELTYP ! type of surface (shell or triangle)
91 INTEGER, DIMENSION(NSURF), INTENT(IN) :: SWIFTSURF
92C-----------------------------------------------
93C L o c a l V a r i a b l e s
94C-----------------------------------------------
95 LOGICAL :: BOOL,CONDITION
96 INTEGER :: I,J,NG
97 INTEGER :: MTN,NEL,NFT,ITY,ISOLNOD,INIVOL
98 INTEGER :: INOD,NODE_NUMBER,FIRST,SURF_NODE_NUMBER,NODE_ID,CLOSEST_NODE
99 INTEGER, DIMENSION(:,:), POINTER :: IXP
100 INTEGER, DIMENSION(:), ALLOCATABLE :: TAG_ARRAY,SURF_NODE_LIST
101
102 INTEGER :: IX,IY,IZ,NEXT_NODE
103 INTEGER :: MY_PHASE,OLD_PHASE
104 INTEGER :: UNKNOWN_CELL_INDEX,LAST_UNKNOWN_CELL,CURRENT_UNKNOWN_CELL
105 INTEGER, DIMENSION(8,3) :: POSITION_SAVE
106 my_real :: dist,eps
107
108 INTEGER :: info
109
110 INTEGER, DIMENSION(:), ALLOCATABLE :: KEY1,KEY2,ID_LIST
111 my_real, DIMENSION(:), ALLOCATABLE :: x_tmp
112 my_real, DIMENSION(3) :: xn
113 my_real, DIMENSION(:,:), ALLOCATABLE :: local_x
114 INTEGER, DIMENSION(:), ALLOCATABLE:: LIST_NODE ! list of node close to the surface
115 INTEGER, DIMENSION(32) :: ID_NODE_SAVE
116 integer, target :: nothing(1,1) !< dummy for indirection
117C-----------------------------------------------
118 ALLOCATE( tag_array(numnod) )
119 ALLOCATE( surf_node_list(numnod) )
120 ALLOCATE( list_node(numnod) )
121 tag_array(1:numnod) = 0
122 surf_node_number = 0
123
124 ! ----------------------
125 ! sort the node of the surface according a direction
126 ! this sorting is useful to reduce the elapsed time of the
127 ! NNS algo
128 DO j=1,4
129 DO i=1,nseg
130 node_id = surface_nodes(i,j)
131 IF(tag_array(node_id) == 0) THEN
132 tag_array(node_id) = 1
133 surf_node_number = surf_node_number + 1
134 surf_node_list(surf_node_number) = node_id
135 ENDIF
136 ENDDO
137 ENDDO
138
139 ALLOCATE( x_tmp(surf_node_number) )
140 ALLOCATE( key2(surf_node_number) )
141 DO i=1,surf_node_number
142 x_tmp(i) = x(leading_dimension,surf_node_list(i))
143 key2(i) = i
144 ENDDO
145
146 ! Sort according to chosen direction
147 CALL myqsort(surf_node_number,x_tmp,key2,info)
148 DEALLOCATE( x_tmp )
149 ! ----------------------
150
151 tag_array(1:numnod) = 0
152 next_node = 0
153 ALLOCATE( local_x(3,8) )
154
155 ! -----------------------
156 ! loop over the solid / quad / triangle elements with 51/151 material
157 DO ng=1,ngroup
158 mtn = iparg(1,ng) ! material law
159 nel = iparg(2,ng) ! number of element
160 nft = iparg(3,ng) ! adress of first element
161 ity = iparg(5,ng) ! type of element
162 isolnod = iparg(28,ng)
163 inivol = iparg(53,ng)
164 IF(inivol <= 0) cycle
165 IF(mtn /= 51 .AND. mtn /= 151) cycle
166 IF(n2d == 0 .AND. ity /= 1)THEN
167 cycle
168 ELSEIF(n2d > 0 .AND. ity /= 7 .AND. ity /= 2)THEN
169 cycle
170 ENDIF
171 ! ---------------
172 ! depending on the king of element
173 IF(ity == 1) THEN
174 first = 1
175 node_number = 8
176 ixp => ixs(1:node_number+1,nft+1:nft+nel)
177 ELSEIF(ity == 2) THEN
178 first = 2
179 node_number = 4
180 ixp => ixq(1:node_number+1,nft+1:nft+nel)
181 ELSEIF(ity == 7) THEN
182 first = 2
183 node_number = 3
184 ixp => ixtg(1:node_number+1,nft+1:nft+nel)
185 ELSE
186 first = -huge(first)
187 node_number = -huge(node_number)
188 ixp => nothing
189 ENDIF
190 ! ---------------
191
192 ! ---------------
193 ! loop over the elements of the group
194 DO j=1,nel
195 old_phase = 0
196 my_phase = 0
197 bool = .false.
198 condition = .true.
199 i = 1
200 unknown_cell_index = 0
201 position_save(1:node_number,1:3) = 0
202 last_unknown_cell = 0
203
204 ! ---------------
205 ! loop over the elements of the group
206 DO WHILE (condition)
207 node_id = ixp(1+i,j) ! node id
208 ix = cell_position(1,node_id) ! position in the grid
209 iy = cell_position(2,node_id) ! position in the grid
210 iz = cell_position(3,node_id) ! position in the grid
211 ! ---------------
212 ! the cell (ix,iy,iz) is crossed by a surface --> need to compute the real distance
213 IF(cell(ix,iy,iz)==2) THEN
214 bool = .true.
215 condition = .false.
216 ! ---------------
217 ! the pseudo distance of the cell (ix,iy,iz) is known --> apply it to the cell
218 ELSEIF(cell(ix,iy,iz) == 1 .OR. cell(ix,iy,iz) == -1) THEN
219 old_phase = my_phase
220 my_phase = cell(ix,iy,iz)
221 ! ---------------
222 ! the pseudo distance of the cell (ix,iy,iz) is unknown --> need to find the pseudo distance
223 ELSEIF(cell(ix,iy,iz) == 0) THEN
224 current_unknown_cell = ix + 1000 * iy + 1000**2 * iz
225 IF(last_unknown_cell /= current_unknown_cell) THEN
226 unknown_cell_index = unknown_cell_index + 1
227 position_save(unknown_cell_index,1) = ix
228 position_save(unknown_cell_index,2) = iy
229 position_save(unknown_cell_index,3) = iz
230 last_unknown_cell = current_unknown_cell
231 id_node_save(unknown_cell_index) = node_id
232 ENDIF
233 ENDIF
234 ! ---------------
235 i = i + 1
236 IF( i > node_number ) condition = .false.
237 ENDDO
238 ! -------------
239
240
241
242 IF(bool) THEN
243 ! -------------
244 ! current element is near a surface, need to compute the distance to the surface
245 DO i=1,node_number
246 node_id = ixp(1+i,j)
247 IF(tag_array(node_id) == 0) THEN
248 tag_array(node_id) = 1
249 next_node = next_node + 1
250 list_node(next_node) = node_id
251 ENDIF
252 ENDDO
253 ELSE
254 ! -------------
255 ! current element is far from a surface, 2 cases :
256 ! * nodes of element are in a non tagged cell --> need to find the phase of the cells
257 ! * at least 1 node is in a tagged cell --> apply the phase to the element & the adjacent cells
258
259 ! -------------
260 ! i found a phase, apply it to the nodes
261 IF(my_phase /= 0) THEN
262 DO i=1,node_number
263 node_id = ixp(1+i,j)
264 nodal_phase(node_id) = my_phase
265 ENDDO
266 DO i=1,unknown_cell_index
267 ix = position_save(i,1)
268 iy = position_save(i,2)
269 iz = position_save(i,3)
270 cell(ix,iy,iz) = my_phase
271 ENDDO
272 ! -------------
273 ! i need to find the phase of the current cells and extend it to the empty cells
274 ELSE
275 ! --------------------
276 ! find the nearest node
277 ALLOCATE( id_list(unknown_cell_index) )
278 ALLOCATE( key1(unknown_cell_index) )
279 DO i=1,unknown_cell_index
280 ix = position_save(i,1)
281 iy = position_save(i,2)
282 iz = position_save(i,3)
283 node_id = id_node_save(i)
284 local_x(1,i) = x(1,node_id)
285 local_x(2,i) = x(2,node_id)
286 local_x(3,i) = x(3,node_id)
287 key1(i) = i
288 ENDDO
289 eps = 1d-6
290 CALL find_closest_node(leading_dimension,unknown_cell_index,surf_node_number,numnod,
291 . local_x,x,eps,
292 . key1,key2,surf_node_list,id_list)
293 ! --------------------
294
295 ! --------------------
296 ! compute the pseudo distance
297 DO i=1,unknown_cell_index
298 inod = id_list(i)
299 xn(1:3) = local_x(1:3,i)
300 dist = zero
301 CALL in_out_side( inod,inod2surf,knod2surf,nnod2surf,x,
302 . xn,dist,nseg,surf_eltyp,nod_normal,
303 . surface_nodes,swiftsurf,id_surface,segtosurf )
304 ix = position_save(i,1)
305 iy = position_save(i,2)
306 iz = position_save(i,3)
307 my_phase = int(dist)
308 cell(ix,iy,iz) = my_phase
309 ENDDO
310
311 ! --------------------
312 ! save the pseudo distance
313 DO i=1,node_number
314 node_id = ixp(1+i,j)
315 ix = cell_position(1,node_id)
316 iy = cell_position(2,node_id)
317 iz = cell_position(3,node_id)
318 my_phase = cell(ix,iy,iz)
319 nodal_phase(node_id) = my_phase
320 ENDDO
321 ! --------------------
322 ! extend the phase to the others cells
323 DO i=1,unknown_cell_index
324 ix = position_save(i,1)
325 iy = position_save(i,2)
326 iz = position_save(i,3)
327 CALL phase_propagation(ix,iy,iz,nb_cell_x,nb_cell_y,nb_cell_z,nb_box_limit,cell)
328 ENDDO
329 ! --------------------
330 DEALLOCATE( id_list )
331 DEALLOCATE( key1 )
332 ENDIF
333 ! -------------
334 ENDIF
335 ENDDO
336 ! ---------------
337 ENDDO
338 ! -----------------------
339
340 DEALLOCATE( local_x )
341
342
343
344 ! -----------------------
345 ! find the nearest node
346 ALLOCATE( local_x(3,next_node) )
347 ALLOCATE( id_list(next_node) )
348 ALLOCATE( key1(next_node) )
349
350 DO i=1,next_node
351 node_id = list_node(i)
352 local_x(1:3,i) = x(1:3,node_id)
353 key1(i) = i
354 ENDDO
355
356 eps = 1d-6
357
358 ! --------------------
359 CALL find_closest_node(leading_dimension,next_node,surf_node_number,numnod,
360 . local_x,x,eps,
361 . key1,key2,surf_node_list,id_list)
362 ! --------------------
363
364 ! --------------------
365 ! compute the pseudo distance & save the closest node id
366 DO i=1,next_node
367 closest_node = id_list(i)
368 node_id = list_node(i)
369 xn(1:3) = local_x(1:3,i)
370 dist = zero
371 CALL in_out_side( closest_node,inod2surf,knod2surf,nnod2surf,x,
372 . xn,dist,nseg,surf_eltyp,nod_normal,
373 . surface_nodes,swiftsurf,id_surface,segtosurf )
374 my_phase = int(dist)
375 nodal_phase(node_id) = my_phase
376 closest_node_id(node_id) = closest_node
377 ENDDO
378 ! --------------------
379
380 ! -----------------------
381
382 DEALLOCATE( key2 )
383 DEALLOCATE( tag_array )
384 DEALLOCATE( surf_node_list )
385 DEALLOCATE( list_node )
386 DEALLOCATE( local_x )
387
388 RETURN
389 END SUBROUTINE phase_detection
my_real
#define my_real
Definition cppsort.cpp:32
find_closest_node
subroutine find_closest_node(leading_direction, n1, n2, n3, x1, x2, eps, key1, key2, id_x2, id_list)
Definition find_closest_node.F:31
in_out_side
subroutine in_out_side(inod, inod2surf, knod2surf, nnod2surf, x, xn, dist, nseg, surf_eltyp, nod_normal, surf_nodes, swiftsurf, idsurf, segtosurf)
Definition in_out_side.F:35
myqsort
subroutine myqsort(n, a, perm, error)
Definition myqsort.F:51
phase_detection
subroutine phase_detection(nparg, ngroup, numels, numelq, numeltg, numnod, nsurf, n2d, leading_dimension, nb_cell_x, nb_cell_y, nb_cell_z, nb_box_limit, iparg, ixs, ixq, ixtg, x, id_surface, cell, cell_position, nodal_phase, closest_node_id, nnod2surf, knod2surf, inod2surf, nod_normal, nseg_used, segtosurf, nseg, surf_eltyp, surface_nodes, swiftsurf)
Definition phase_detection.F:40
phase_propagation
subroutine phase_propagation(ix, iy, iz, nb_cell_x, nb_cell_y, nb_cell_z, nb_box_limit, cell)
Definition phase_propagation.F:29