OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
i18dst3.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!|| i18dst3_mod ../engine/source/interfaces/int18/i18dst3.F
25!||--- called by ------------------------------------------------------
26!|| i18main_kine_i ../engine/source/interfaces/int18/i18main_kine.F
27!|| i7mainf ../engine/source/interfaces/int07/i7mainf.F
28!||====================================================================
29 MODULE i18dst3_mod
30 CONTAINS
31!||====================================================================
32!|| i18dst3 ../engine/source/interfaces/int18/i18dst3.F
33!||--- called by ------------------------------------------------------
34!|| i18main_kine_i ../engine/source/interfaces/int18/i18main_kine.F
35!|| i7mainf ../engine/source/interfaces/int07/i7mainf.F
36!||--- uses -----------------------------------------------------
37!|| ale_connectivity_mod ../common_source/modules/ale/ale_connectivity_mod.F
38!|| multi_fvm_mod ../common_source/modules/ale/multi_fvm_mod.F90
39!|| tri7box ../engine/share/modules/tri7box.F
40!||====================================================================
41 SUBROUTINE i18dst3(
42 1 JLT ,CAND_N ,CAND_E ,CN_LOC ,CE_LOC ,
43 2 X1 ,X2 ,X3 ,X4 ,Y1 ,
44 3 Y2 ,Y3 ,Y4 ,Z1 ,Z2 ,
45 4 Z3 ,Z4 ,XI ,YI ,ZI ,
46 5 NX1 ,NX2 ,NX3 ,NX4 ,NY1 ,
47 6 NY2 ,NY3 ,NY4 ,NZ1 ,NZ2 ,
48 7 NZ3 ,NZ4 ,LB1 ,LB2 ,LB3 ,
49 8 LB4 ,LC1 ,LC2 ,LC3 ,LC4 ,
50 9 P1 ,P2 ,P3 ,P4 ,IX1 ,
51 A IX2 ,IX3 ,IX4 ,NSVG ,STIF ,
52 B JLT_NEW ,GAPV ,CAND_P ,ALE_NE_CONNECT ,
53 C INDEX ,VXI ,VYI ,ITAB ,XCELL ,
54 D VZI ,MSI ,KINI ,
55 E IGAP ,MULTI_FVM ,S_XCELL_REMOTE ,XCELL_REMOTE)
56C-----------------------------------------------
57C D e s c r i p t i o n
58C-----------------------------------------------
59C This subroutine is computing GAP value (if not constant)
60C and also penetration into the gap.
61C-----------------------------------------------
62C M o d u l e s
63C-----------------------------------------------
64 USE ale_connectivity_mod
65 USE multi_fvm_mod
66 USE tri7box
67C-----------------------------------------------
68C I m p l i c i t T y p e s
69C-----------------------------------------------
70#include "implicit_f.inc"
71#include "comlock.inc"
72#include "com04_c.inc"
73#include "tabsiz_c.inc"
74C-----------------------------------------------
75C G l o b a l P a r a m e t e r s
76C-----------------------------------------------
77#include "mvsiz_p.inc"
78C-----------------------------------------------
79C D u m m y A r g u m e n t s
80C-----------------------------------------------
81 INTEGER JLT, JLT_NEW, CAND_N(*),CN_LOC(MVSIZ),
82 . CAND_E(*),CE_LOC(MVSIZ), NSVG(MVSIZ), KINI(*),ICURV
83 INTEGER IX1(MVSIZ), IX2(MVSIZ), IX3(MVSIZ), IX4(MVSIZ),
84 . INDEX(MVSIZ)
85 INTEGER, INTENT(IN) :: ITAB(NUMNOD),IGAP
86 INTEGER, INTENT(in) :: S_XCELL_REMOTE
87 my_real
88 . NX1(MVSIZ), NX2(MVSIZ), NX3(MVSIZ), NX4(MVSIZ),
89 . NY1(MVSIZ), NY2(MVSIZ), NY3(MVSIZ), NY4(MVSIZ),
90 . NZ1(MVSIZ), NZ2(MVSIZ), NZ3(MVSIZ), NZ4(MVSIZ),
91 . LB1(MVSIZ), LB2(MVSIZ), LB3(MVSIZ), LB4(MVSIZ),
92 . LC1(MVSIZ), LC2(MVSIZ), LC3(MVSIZ), LC4(MVSIZ),
93 . KB1(MVSIZ), KB2(MVSIZ), KB3(MVSIZ), KB4(MVSIZ),
94 . KC1(MVSIZ), KC2(MVSIZ), KC3(MVSIZ), KC4(MVSIZ),
95 . x1(mvsiz), x2(mvsiz), x3(mvsiz), x4(mvsiz),
96 . y1(mvsiz), y2(mvsiz), y3(mvsiz), y4(mvsiz),
97 . z1(mvsiz), z2(mvsiz), z3(mvsiz), z4(mvsiz),
98 . xi(mvsiz), yi(mvsiz), zi(mvsiz), stif(mvsiz),
99 . p1(mvsiz), p2(mvsiz), p3(mvsiz), p4(mvsiz),
100 . gapv(mvsiz), cand_p(*),
101 . vxi(mvsiz), vyi(mvsiz), vzi(mvsiz), msi(mvsiz),
102 . xcell(3,sxcell)
103 TYPE(multi_fvm_struct), INTENT(IN) :: MULTI_FVM
104 TYPE(t_connectivity), INTENT(IN) :: ALE_NE_CONNECT
105 my_real, DIMENSION(S_XCELL_REMOTE), INTENT(in) :: XCELL_REMOTE
106C-----------------------------------------------
107C L o c a l V a r i a b l e s
108C-----------------------------------------------
109 INTEGER I, IG, J, IEL, IAD1, IAD2, ELEM_ID
110 my_real
111 . X0(MVSIZ), Y0(MVSIZ), Z0(MVSIZ),
112 . AL1(MVSIZ), AL2(MVSIZ), AL3(MVSIZ), AL4(MVSIZ),
113 . X01(MVSIZ), X02(MVSIZ), X03(MVSIZ), X04(MVSIZ),
114 . Y01(MVSIZ), Y02(MVSIZ), Y03(MVSIZ), Y04(MVSIZ),
115 . Z01(MVSIZ), Z02(MVSIZ), Z03(MVSIZ), Z04(MVSIZ),
116 . XI1(MVSIZ), XI2(MVSIZ), XI3(MVSIZ), XI4(MVSIZ),
117 . yi1(mvsiz), yi2(mvsiz), yi3(mvsiz), yi4(mvsiz),
118 . zi1(mvsiz), zi2(mvsiz), zi3(mvsiz), zi4(mvsiz),
119 . pene2(mvsiz),
120 . hlb1(mvsiz), hlc1(mvsiz), hlb2(mvsiz),hlc2(mvsiz),
121 . hlb3(mvsiz), hlc3(mvsiz), hlb4(mvsiz),hlc4(mvsiz)
122 my_real
123 . s2,d1,d2,d3,d4,dl,
124 . x12,x23,x34,x41,xi0,sx1,sx2,sx3,sx4,sx0,
125 . y12,y23,y34,y41,yi0,sy1,sy2,sy3,sy4,sy0,
126 . z12,z23,z34,z41,zi0,sz1,sz2,sz3,sz4,sz0,
127 . gap2(mvsiz), ds2,t1,t2,t3,
128 . al1num,al2num,al3num,al4num,al1den,al2den,al3den,al4den,
129 . x23d,y23d,z23d,x34d,y34d,z34d,x41d,y41d,z41d,
130 . x12d,y12d,z12d,gap2d,xi0d,yi0d,zi0d,s2d, la, hla, aaa,
131 . xi0v(mvsiz), yi0v(mvsiz), zi0v(mvsiz)
132 my_real
133 . nnx1(mvsiz), nnx2(mvsiz), nnx3(mvsiz), nnx4(mvsiz),
134 . nny1(mvsiz), nny2(mvsiz), nny3(mvsiz), nny4(mvsiz),
135 . nnz1(mvsiz), nnz2(mvsiz), nnz3(mvsiz), nnz4(mvsiz)
136 LOGICAL IS_GAP_CONSTANT, IS_GAP_VARIABLE
137 INTEGER :: NODE_ID
138C-----------------------------------------------
139C S o u r c e L i n e s
140C-----------------------------------------------
141
142 IS_GAP_CONSTANT = .true.
143 is_gap_variable = .false.
144 IF(igap == 1)THEN
145 is_gap_constant = .false.
146 is_gap_variable = .true.
147 ENDIF
148
149 !Reference point for triangular decomposition
150 DO i=1,jlt
151 IF(ix3(i) /= ix4(i))THEN
152 !mean point is the top of each 4 triangles composing the 4-nodes face
153 x0(i) = fourth*(x1(i)+x2(i)+x3(i)+x4(i))
154 y0(i) = fourth*(y1(i)+y2(i)+y3(i)+y4(i))
155 z0(i) = fourth*(z1(i)+z2(i)+z3(i)+z4(i))
156 ELSE
157 !already a triangle
158 x0(i) = x3(i)
159 y0(i) = y3(i)
160 z0(i) = z3(i)
161 ENDIF
162 ENDDO
163
164 !COMPUTING VECTORS using each top of the face (from salve point I and from reference point 0 on the face)
165 DO i=1,jlt
166 ! Vector 01 : from centroid to 1st point of the MAINface
167 x01(i) = x1(i) - x0(i)
168 y01(i) = y1(i) - y0(i)
169 z01(i) = z1(i) - z0(i)
170 ! Vector 02 : from centroid to 2nd point of the MAIN face
171 x02(i) = x2(i) - x0(i)
172 y02(i) = y2(i) - y0(i)
173 z02(i) = z2(i) - z0(i)
174 ! Vector 03 : from centroid to 3rd point of the MAIN face
175 x03(i) = x3(i) - x0(i)
176 y03(i) = y3(i) - y0(i)
177 z03(i) = z3(i) - z0(i)
178 ! Vector 04 : from centroid to 4th point of the MAIN face
179 x04(i) = x4(i) - x0(i)
180 y04(i) = y4(i) - y0(i)
181 z04(i) = z4(i) - z0(i)
182 ! Vector IOV : from SECONDARY node to centroid of the MAIN face
183 xi0v(i) = x0(i) - xi(i)
184 yi0v(i) = y0(i) - yi(i)
185 zi0v(i) = z0(i) - zi(i)
186 ! Vector 1I : from salve node to 1st point of the face
187 xi1(i) = x1(i) - xi(i)
188 yi1(i) = y1(i) - yi(i)
189 zi1(i) = z1(i) - zi(i)
190 ! Vector 2I : from salve node to 2nd point of the face
191 xi2(i) = x2(i) - xi(i)
192 yi2(i) = y2(i) - yi(i)
193 zi2(i) = z2(i) - zi(i)
194 ! Vector 3I : from salve node to 3rd point of the face
195 xi3(i) = x3(i) - xi(i)
196 yi3(i) = y3(i) - yi(i)
197 zi3(i) = z3(i) - zi(i)
198 ! Vector 4I : from salve node to 4th point of the face
199 xi4(i) = x4(i) - xi(i)
200 yi4(i) = y4(i) - yi(i)
201 zi4(i) = z4(i) - zi(i)
202 ENDDO
203
204 DO i=1,jlt
205 ! IO (x) I1
206 sx1 = yi0v(i)*zi1(i) - zi0v(i)*yi1(i)
207 sy1 = zi0v(i)*xi1(i) - xi0v(i)*zi1(i)
208 sz1 = xi0v(i)*yi1(i) - yi0v(i)*xi1(i)
209 ! IO (x) I2
210 sx2 = yi0v(i)*zi2(i) - zi0v(i)*yi2(i)
211 sy2 = zi0v(i)*xi2(i) - xi0v(i)*zi2(i)
212 sz2 = xi0v(i)*yi2(i) - yi0v(i)*xi2(i)
213 ! IO (x) I3
214 sx3 = yi0v(i)*zi3(i) - zi0v(i)*yi3(i)
215 sy3 = zi0v(i)*xi3(i) - xi0v(i)*zi3(i)
216 sz3 = xi0v(i)*yi3(i) - yi0v(i)*xi3(i)
217 ! IO (x) I4
218 sx4 = yi0v(i)*zi4(i) - zi0v(i)*yi4(i)
219 sy4 = zi0v(i)*xi4(i) - xi0v(i)*zi4(i)
220 sz4 = xi0v(i)*yi4(i) - yi0v(i)*xi4(i)
221
222 !---TRIANGLE_1-----------------------
223 ! normal vector (2Sn, where ||n||=1 )
224 sx0 = y01(i)*z02(i) - z01(i)*y02(i)
225 sy0 = z01(i)*x02(i) - x01(i)*z02(i)
226 sz0 = x01(i)*y02(i) - y01(i)*x02(i)
227 s2 = one/max(em30,sx0*sx0 + sy0*sy0 + sz0*sz0)
228 nnx1(i) = sx0
229 nny1(i) = sy0
230 nnz1(i) = sz0
231 lb1(i) = -(sx0*sx2 + sy0*sy2 + sz0*sz2) * s2
232 lc1(i) = (sx0*sx1 + sy0*sy1 + sz0*sz1) * s2
233
234 !---TRIANGLE_2-----------------------
235 ! normal vector (2Sn, where ||n||=1 )
236 sx0 = y02(i)*z03(i) - z02(i)*y03(i)
237 sy0 = z02(i)*x03(i) - x02(i)*z03(i)
238 sz0 = x02(i)*y03(i) - y02(i)*x03(i)
239 s2 = 1./max(em30,sx0*sx0 + sy0*sy0 + sz0*sz0)
240 nnx2(i) = sx0
241 nny2(i) = sy0
242 nnz2(i) = sz0
243 lb2(i) = -(sx0*sx3 + sy0*sy3 + sz0*sz3) * s2
244 lc2(i) = (sx0*sx2 + sy0*sy2 + sz0*sz2) * s2
245
246 !---TRIANGLE_3-----------------------
247 ! normal vector (2Sn, where ||n||=1 )
248 sx0 = y03(i)*z04(i) - z03(i)*y04(i)
249 sy0 = z03(i)*x04(i) - x03(i)*z04(i)
250 sz0 = x03(i)*y04(i) - y03(i)*x04(i)
251 s2 = one/max(em30,sx0*sx0 + sy0*sy0 + sz0*sz0)
252 nnx3(i) = sx0
253 nny3(i) = sy0
254 nnz3(i) = sz0
255 lb3(i) = -(sx0*sx4 + sy0*sy4 + sz0*sz4) * s2
256 lc3(i) = (sx0*sx3 + sy0*sy3 + sz0*sz3) * s2
257
258 !---TRIANGLE_3-----------------------
259 ! normal vector (2Sn, where ||n||=1 )
260 sx0 = y04(i)*z01(i) - z04(i)*y01(i)
261 sy0 = z04(i)*x01(i) - x04(i)*z01(i)
262 sz0 = x04(i)*y01(i) - y04(i)*x01(i)
263 s2 = one/max(em30,sx0*sx0 + sy0*sy0 + sz0*sz0)
264 nnx4(i) = sx0
265 nny4(i) = sy0
266 nnz4(i) = sz0
267 lb4(i) = -(sx0*sx1 + sy0*sy1 + sz0*sz1) * s2
268 lc4(i) = (sx0*sx4 + sy0*sy4 + sz0*sz4) * s2
269 ENDDO
270 IF(igap == 1)THEN !---VARIABLE GAP, CALCULATED HERE
271 IF (.NOT.multi_fvm%IS_USED) THEN
272 DO i=1,jlt
273 node_id = nsvg(i)
274 ! ----------
275 ! local node
276 IF(node_id>0) THEN
277 iad1 = ale_ne_connect%IAD_CONNECT(node_id)
278 iad2 = ale_ne_connect%IAD_CONNECT(node_id + 1) - 1
279 dl = zero
280 DO iel = iad1, iad2
281 elem_id = ale_ne_connect%CONNECTED(iel)
282 dl=max(dl, xcell(1,elem_id))
283 ENDDO
284 ELSE
285 node_id = abs(node_id)
286 dl = xcell_remote(node_id)
287 ENDIF
288 ! ----------
289 gapv(i) = three_half*dl
290 enddo!next I
291 ELSEIF(multi_fvm%IS_USED) THEN
292 DO i=1,jlt
293 IF(nsvg(i)>0) THEN
294 elem_id = nsvg(i) - numnod
295 dl=xcell(1,elem_id)
296 ELSE
297 node_id = abs(nsvg(i))
298 dl = xcell_remote(node_id)
299 ENDIF
300 !gap value
301 gapv(i) = three_half*dl
302 enddo!next I
303 ENDIF
304 ELSE !---CONSTANT GAP INPUT BY USER (IGAP=0=
305 gapv(1:jlt)=gapv(1)
306 ENDIF
307
308 IF(is_gap_constant)THEN
309 !Gap is constant for interface 18 : !GAPV(1:JLT)=GAP
310 aaa =gapv(1)*gapv(1)
311 DO i=1,jlt
312 gap2(i)=aaa
313 ENDDO
314 ELSE
315 DO i=1,jlt
316 gap2(i)=gapv(i)*gapv(i)
317 ENDDO
318 ENDIF
319
320 DO i=1,jlt
321 aaa = one/max(em30,x01(i)*x01(i)+y01(i)*y01(i)+z01(i)*z01(i))
322 hlc1(i)= lc1(i)*abs(lc1(i))*aaa
323 hlb4(i)= lb4(i)*abs(lb4(i))*aaa
324 al1(i) = -(xi0v(i)*x01(i)+yi0v(i)*y01(i)+zi0v(i)*z01(i))*aaa
325 al1(i) = max(zero,min(one,al1(i)))
326 aaa = one/max(em30,x02(i)*x02(i)+y02(i)*y02(i)+z02(i)*z02(i))
327 hlc2(i)= lc2(i)*abs(lc2(i))*aaa
328 hlb1(i)= lb1(i)*abs(lb1(i))*aaa
329 al2(i) = -(xi0v(i)*x02(i)+yi0v(i)*y02(i)+zi0v(i)*z02(i))*aaa
330 al2(i) = max(zero,min(one,al2(i)))
331 aaa = one/max(em30,x03(i)*x03(i)+y03(i)*y03(i)+z03(i)*z03(i))
332 hlc3(i)= lc3(i)*abs(lc3(i))*aaa
333 hlb2(i)= lb2(i)*abs(lb2(i))*aaa
334 al3(i) = -(xi0v(i)*x03(i)+yi0v(i)*y03(i)+zi0v(i)*z03(i))*aaa
335 al3(i) = max(zero,min(one,al3(i)))
336 aaa = one/max(em30,x04(i)*x04(i)+y04(i)*y04(i)+z04(i)*z04(i))
337 hlc4(i)= lc4(i)*abs(lc4(i))*aaa
338 hlb3(i)= lb3(i)*abs(lb3(i))*aaa
339 al4(i) = -(xi0v(i)*x04(i)+yi0v(i)*y04(i)+zi0v(i)*z04(i))*aaa
340 al4(i) = max(zero,min(one,al4(i)))
341 ENDDO
342
343 !computes local coordinates (LB1,LC1) on triangle_1 which minimize the distance from SECONDARY point I.
344 ! (XI,YI,ZI) is SECONDARY point.
345 ! let (LB1,LV1) be a point of triangle_1 : 0 < LB1 < 1 and 0 < LC1 < 1
346 ! find (LB1,LC1) which minimize the distance from SECONDARY point I
347 DO i=1,jlt
348 x12 = x2(i) - x1(i)
349 y12 = y2(i) - y1(i)
350 z12 = z2(i) - z1(i)
351 la = one - lb1(i) - lc1(i)
352 aaa = one / max(em20,x12*x12+y12*y12+z12*z12)
353 hla= la*abs(la) * aaa
354 IF(la < zero.AND.hla <= hlb1(i).AND.hla <= hlc1(i))THEN
355 lb1(i) = (xi2(i)*x12+yi2(i)*y12+zi2(i)*z12) * aaa
356 lb1(i) = max(zero,min(one,lb1(i)))
357 lc1(i) = one - lb1(i)
358 ELSEIF(lb1(i) < zero.AND.hlb1(i) <= hlc1(i).AND.hlb1(i) <= hla)THEN
359 lb1(i) = zero
360 lc1(i) = al2(i)
361 ELSEIF(lc1(i) < zero.AND.hlc1(i) <= hla.AND.hlc1(i) <= hlb1(i))THEN
362 lc1(i) = zero
363 lb1(i) = al1(i)
364 ENDIF
365 ENDDO
366
367 !computes local coordinates (LB2,LC2) on triangle_2 which minimize the distance from SECONDARY point I.
368 ! (XI,YI,ZI) is SECONDARY point.
369 ! let (LB2,LV2) be a point of triangle_2 : 0 < LB2 < 1 and 0 < LC2 < 1
370 ! find (LB2,LC2) which minimize the distance from SECONDARY point I
371 DO i=1,jlt
372 x23 = x3(i) - x2(i)
373 y23 = y3(i) - y2(i)
374 z23 = z3(i) - z2(i)
375 la = one - lb2(i) - lc2(i)
376 aaa = one / max(em20,x23*x23+y23*y23+z23*z23)
377 hla= la*abs(la) * aaa
378 IF(la < zero.AND.hla <= hlb2(i).AND.hla <= hlc2(i))THEN
379 lb2(i) = (xi3(i)*x23+yi3(i)*y23+zi3(i)*z23)*aaa
380 lb2(i) = max(zero,min(one,lb2(i)))
381 lc2(i) = one - lb2(i)
382 ELSEIF(lb2(i) < zero.AND.hlb2(i) <= hlc2(i).AND.hlb2(i) <= hla)THEN
383 lb2(i) = zero
384 lc2(i) = al3(i)
385 ELSEIF(lc2(i) < zero.AND.hlc2(i) <= hla.AND.hlc2(i) <= hlb2(i))THEN
386 lc2(i) = zero
387 lb2(i) = al2(i)
388 ENDIF
389 ENDDO
390
391 !computes local coordinates (LB3,LC3) on triangle_3 which minimize the distance from SECONDARY point I.
392 ! (XI,YI,ZI) is SECONDARY point.
393 ! let (LB3,LV3) be a point of triangle_3 : 0 < LB3 < 1 and 0 < LC3 < 1
394 ! find (LB3,LC3) which minimize the distance from SECONDARY point I
395 DO i=1,jlt
396 x34 = x4(i) - x3(i)
397 y34 = y4(i) - y3(i)
398 z34 = z4(i) - z3(i)
399 la = one - lb3(i) - lc3(i)
400 aaa = one / max(em20,x34*x34+y34*y34+z34*z34)
401 hla= la*abs(la) * aaa
402 IF(la < zero.AND.hla <= hlb3(i).AND.hla <= hlc3(i))THEN
403 lb3(i) = (xi4(i)*x34+yi4(i)*y34+zi4(i)*z34)*aaa
404 lb3(i) = max(zero,min(one,lb3(i)))
405 lc3(i) = one - lb3(i)
406 ELSEIF(lb3(i) < zero.AND.hlb3(i) <= hlc3(i).AND.hlb3(i) <= hla)THEN
407 lb3(i) = zero
408 lc3(i) = al4(i)
409 ELSEIF(lc3(i) < zero.AND.hlc3(i) <= hla.AND.hlc3(i) <= hlb3(i))THEN
410 lc3(i) = zero
411 lb3(i) = al3(i)
412 ENDIF
413 ENDDO
414
415 !computes local coordinates (LB4,LC4) on triangle_4 which minimize the distance from SECONDARY point I.
416 ! (XI,YI,ZI) is SECONDARY point.
417 ! let (LB4,LV4) be a point of triangle_4 : 0 < LB4 < 1 and 0 < LC4 < 1
418 ! find (LB4,LC4) which minimize the distance from SECONDARY point I
419 DO i=1,jlt
420 x41 = x1(i) - x4(i)
421 y41 = y1(i) - y4(i)
422 z41 = z1(i) - z4(i)
423 la = one - lb4(i) - lc4(i)
424 aaa = one / max(em20,x41*x41+y41*y41+z41*z41)
425 hla= la*abs(la) * aaa
426 IF(la < zero.AND.hla <= hlb4(i).AND.hla <= hlc4(i))THEN
427 lb4(i) = (xi1(i)*x41+yi1(i)*y41+zi1(i)*z41)*aaa
428 lb4(i) = max(zero,min(one,lb4(i)))
429 lc4(i) = one - lb4(i)
430 ELSEIF(lb4(i) < zero.AND.hlb4(i) <= hlc4(i).AND.hlb4(i) <= hla)THEN
431 lb4(i) = zero
432 lc4(i) = al1(i)
433 ELSEIF(lc4(i) < zero.AND.hlc4(i) <= hla.AND.hlc4(i) <= hlb4(i))THEN
434 lc4(i) = zero
435 lb4(i) = al4(i)
436 ENDIF
437 ENDDO
438
439 !Use minimal distance from each triangle to compute penetration
440 DO i=1,jlt
441 !----MINIMAL DISTANCE FROM TRIANGLE_1
442 nx1(i) = xi(i)-(x0(i) + lb1(i)*x01(i) + lc1(i)*x02(i))
443 ny1(i) = yi(i)-(y0(i) + lb1(i)*y01(i) + lc1(i)*y02(i))
444 nz1(i) = zi(i)-(z0(i) + lb1(i)*z01(i) + lc1(i)*z02(i))
445 p1(i) = nx1(i)*nx1(i) + ny1(i)*ny1(i) +nz1(i)*nz1(i)
446 !----minimal distance from triangle_2
447 nx2(i) = xi(i)-(x0(i) + lb2(i)*x02(i) + lc2(i)*x03(i))
448 ny2(i) = yi(i)-(y0(i) + lb2(i)*y02(i) + lc2(i)*y03(i))
449 nz2(i) = zi(i)-(z0(i) + lb2(i)*z02(i) + lc2(i)*z03(i))
450 p2(i) = nx2(i)*nx2(i) + ny2(i)*ny2(i) +nz2(i)*nz2(i)
451 !----MINIMAL DISTANCE FROM TRIANGLE_3
452 nx3(i) = xi(i)-(x0(i) + lb3(i)*x03(i) + lc3(i)*x04(i))
453 ny3(i) = yi(i)-(y0(i) + lb3(i)*y03(i) + lc3(i)*y04(i))
454 nz3(i) = zi(i)-(z0(i) + lb3(i)*z03(i) + lc3(i)*z04(i))
455 p3(i) = nx3(i)*nx3(i) + ny3(i)*ny3(i) +nz3(i)*nz3(i)
456 !----MINIMAL DISTANCE FROM TRIANGLE_4
457 nx4(i) = xi(i)-(x0(i) + lb4(i)*x04(i) + lc4(i)*x01(i))
458 ny4(i) = yi(i)-(y0(i) + lb4(i)*y04(i) + lc4(i)*y01(i))
459 nz4(i) = zi(i)-(z0(i) + lb4(i)*z04(i) + lc4(i)*z01(i))
460 p4(i) = nx4(i)*nx4(i) + ny4(i)*ny4(i) +nz4(i)*nz4(i)
461 !---------estimation calculated to skip node out of gap
462 d1 = max(zero, gap2(i) - p1(i))
463 d2 = max(zero, gap2(i) - p2(i))
464 d3 = max(zero, gap2(i) - p3(i))
465 d4 = max(zero, gap2(i) - p4(i))
466 pene2(i) = max(d1,d2,d3,d4) !PENE2 = GAP^2 - DIST^2 : PENE**2 =0 => PENE=0
467
468 ENDDO
469
470
471 !Store local coordinates on each triangles for retaines nodes (PENE**2 >0)
472 DO i=1,jlt
473 IF(pene2(i) == zero.OR.stif(i) == zero)THEN
474 cand_p(index(i))=0
475 ENDIF
476 ENDDO
477 DO i=1,jlt
478 IF(pene2(i) /= zero.AND.stif(i) /= zero)THEN
479 jlt_new = jlt_new + 1
480 cn_loc(jlt_new) = cand_n(i)
481 ce_loc(jlt_new) = cand_e(i)
482 ix1(jlt_new) = ix1(i)
483 ix2(jlt_new) = ix2(i)
484 ix3(jlt_new) = ix3(i)
485 ix4(jlt_new) = ix4(i)
486 nsvg(jlt_new) = nsvg(i)
487 nx1(jlt_new) = nnx1(i)
488 nx2(jlt_new) = nnx2(i)
489 nx3(jlt_new) = nnx3(i)
490 nx4(jlt_new) = nnx4(i)
491 ny1(jlt_new) = nny1(i)
492 ny2(jlt_new) = nny2(i)
493 ny3(jlt_new) = nny3(i)
494 ny4(jlt_new) = nny4(i)
495 nz1(jlt_new) = nnz1(i)
496 nz2(jlt_new) = nnz2(i)
497 nz3(jlt_new) = nnz3(i)
498 nz4(jlt_new) = nnz4(i)
499 p1(jlt_new) = p1(i)
500 p2(jlt_new) = p2(i)
501 p3(jlt_new) = p3(i)
502 p4(jlt_new) = p4(i)
503 lb1(jlt_new) = lb1(i)
504 lb2(jlt_new) = lb2(i)
505 lb3(jlt_new) = lb3(i)
506 lb4(jlt_new) = lb4(i)
507 lc1(jlt_new) = lc1(i)
508 lc2(jlt_new) = lc2(i)
509 lc3(jlt_new) = lc3(i)
510 lc4(jlt_new) = lc4(i)
511 stif(jlt_new) = stif(i)
512 gapv(jlt_new) = gapv(i)
513 index(jlt_new) = index(i)
514 kini(jlt_new) = kini(i)
515 vxi(jlt_new) = vxi(i)
516 vyi(jlt_new) = vyi(i)
517 vzi(jlt_new) = vzi(i)
518 msi(jlt_new) = msi(i)
519 ENDIF
520 ENDDO
521C---------------------
522 RETURN
523 END SUBROUTINE i18dst3
524C===============================================================================
525 END MODULE i18dst3_mod
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
ale_connectivity_mod
Definition ale_connectivity_mod.F:223
i18dst3_mod
Definition i18dst3.F:29
i18dst3_mod::i18dst3
subroutine i18dst3(jlt, cand_n, cand_e, cn_loc, ce_loc, x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4, xi, yi, zi, nx1, nx2, nx3, nx4, ny1, ny2, ny3, ny4, nz1, nz2, nz3, nz4, lb1, lb2, lb3, lb4, lc1, lc2, lc3, lc4, p1, p2, p3, p4, ix1, ix2, ix3, ix4, nsvg, stif, jlt_new, gapv, cand_p, ale_ne_connect, index, vxi, vyi, itab, xcell, vzi, msi, kini, igap, multi_fvm, s_xcell_remote, xcell_remote)
Definition i18dst3.F:56