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