OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
i11pen3.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!|| i11pen3_vox ../engine/source/interfaces/intsort/i11pen3.F
25!||--- called by ------------------------------------------------------
26!|| i11sto_vox ../engine/source/interfaces/intsort/i11sto.F
27!||--- uses -----------------------------------------------------
28!|| tri7box ../engine/share/modules/tri7box.F
29!||====================================================================
30 SUBROUTINE i11pen3_vox(JLT ,CAND_S,CAND_M,GAPMIN,DRAD,
31 . MARGE ,GAP_S ,GAP_M ,GAP_S_L,GAP_M_L ,
32 . IGAP ,X ,IRECTS ,IRECTM ,PENE ,
33 . NRTS ,DGAPLOAD)
34
35C-----------------------------------------------
36C M o d u l e s
37C-----------------------------------------------
38 USE tri7box
39C-----------------------------------------------
40C I m p l i c i t T y p e s
41C-----------------------------------------------
42#include "implicit_f.inc"
43C-----------------------------------------------
44C G l o b a l P a r a m e t e r s
45C-----------------------------------------------
46#include "mvsiz_p.inc"
47C-----------------------------------------------
48C D u m m y A r g u m e n t s
49C-----------------------------------------------
50 INTEGER JLT, NRTS,IGAP
51 INTEGER IRECTS(2,*), IRECTM(2,*),CAND_S(*),CAND_M(*)
52 my_real
53 . gapmin,marge
54 my_real
55 . x(3,*), pene(mvsiz),
56 . gap_s(*), gap_m(*), gap_s_l(*), gap_m_l(*)
57 my_real , INTENT(IN) :: dgapload,drad
58C-----------------------------------------------
59C L o c a l V a r i a b l e s
60C-----------------------------------------------
61 INTEGER I, IG,N1,N2,M1,M2,NI,L,J
62 my_real
63 . XS12,YS12,ZS12,XM12,YM12,ZM12,XA,XB,
64 . xs2,xm2,xsm,xs2m2,ys2,ym2,ysm,ys2m2,zs2,zm2,zsm,zs2m2,
65 . xx,yy,zz,als,alm,det,
66 . gap2, x11, x12, x13, x21, x22, x23,
67 . xmax1,ymax1,zmax1,xmax2,ymax2,zmax2,
68 . xmin1,ymin1,zmin1,xmin2,ymin2,zmin2,dd,gapv(mvsiz)
69 my_real :: maxgapv
70C-----------------------------------------------
71 IF(igap==0)THEN
72 DO i=1,jlt
73 gapv(i)=max(drad,gapmin+dgapload)+marge
74 ENDDO
75 ELSE
76 DO i=1,jlt
77 l = cand_s(i)
78 IF( l <= nrts) THEN
79 gapv(i)=gap_s(l)+gap_m(cand_m(i))
80 IF(igap == 3)
81 . gapv(i)=min(gap_s_l(l)+gap_m_l(cand_m(i)),gapv(i))
82 gapv(i)=max(drad,max(gapmin,gapv(i))+dgapload)+marge
83 ELSE
84 gapv(i)=xrem(16,l-nrts )+gap_m(cand_m(i))
85 IF(igap == 3)
86 . gapv(i)=min(xrem(17,l-nrts)+gap_m_l(cand_m(i)),gapv(i))
87 gapv(i)=max(drad,max(gapmin,gapv(i))+dgapload)+marge
88 ENDIF ! L <= NRTS
89 ENDDO ! JLT
90 ENDIF ! IGAP ==0
91
92c MAXGAPV = GAPV(1)
93c DO I = 2,JLT
94c MAXGAPV = MAX(MAXGAPV,GAPV(I))
95c ENDDO
96c WRITE(6,*) __FILE__,__LINE__,"MAXGAPV=",MAXGAPV
97C--------------------------------------------------------
98C
99C--------------------------------------------------------
100C F = [A*X1+(1-A)*X2-B*X3-(1-B)*X4]^2 + [..Y..]^2 + [..Z..]^2
101C DF/DA = 0 = (X1-X2)(A(X1-X2)+X2-X4 +B(X4-X3))+...
102C DF/DA = 0 = A(X1-X2)^2 +X2-X4 + B(X1-X2)(X4-X3))+...
103C DF/DA = 0 = A[(X1-X2)^2 + (Y1-Y2)^2 + (Z1-Z2)^2]
104C + B[(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
105C + (X1-X2)(X2-X4) + (Y1-Y2)(Y2-Y4) + (Z1-Z2)(Z2-Z4)
106C DF/DB = 0 = (X4-X3)(A(X1-X2)+X2-X4 +B(X4-X3))+...
107C DF/DB = 0 = B[(X4-X3)^2 + (Y4-Y3)^2 + (Z4-Z3)^2]
108C + A[(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
109C + (X4-X3)(X2-X4) + (Y4-Y3)(Y2-Y4) + (Z4-Z3)(Z2-Z4)
110C XS2 = [(X1-X2)^2 + (Y1-Y2)^2 + (Z1-Z2)^2]
111C XM2 = [(X4-X3)^2 + (Y4-Y3)^2 + (Z4-Z3)^2]
112C XSM = [(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
113C XA = (X1-X2)(X2-X4) + (Y1-Y2)(Y2-Y4) + (Z1-Z2)(Z2-Z4)
114C XB = (X4-X3)(X2-X4) + (Y4-Y3)(Y2-Y4) + (Z4-Z3)(Z2-Z4)
115C A XS2 + B XSM + XA = 0
116C A XSM + B XM2 + XB = 0
117C
118C A = -(XA + B XSM)/XS2
119C -(XA + B XSM)*XSM + B XM2*XS2 + XB*XS2 = 0
120C -B XSM*XSM + B XM2*XS2 + XB*XS2-XA*XSM = 0
121C B*(XM2*XS2 - XSM*XSM) = -XB*XS2+XA*XSM
122C B = (XA*XSM-XB*XS2) / (XM2*XS2 - XSM*XSM)
123C A = (XB*XSM-XA*XM2) / (XM2*XS2 - XSM*XSM)
124C--------------------------------------------------------
125
126C
127 DO i=1,jlt
128 l = cand_s(i)
129 IF(l<=nrts) THEN
130 ni=0
131 n1=irects(1,cand_s(i))
132 n2=irects(2,cand_s(i))
133 x11 = x(1,n1)
134 x12 = x(2,n1)
135 x13 = x(3,n1)
136 x21 = x(1,n2)
137 x22 = x(2,n2)
138 x23 = x(3,n2)
139 ELSE
140 ni = l - nrts
141 x11 = xrem(1,ni)
142 x12 = xrem(2,ni)
143 x13 = xrem(3,ni)
144 x21 = xrem(8,ni)
145 x22 = xrem(9,ni)
146 x23 = xrem(10,ni)
147 END IF
148 m1=irectm(1,cand_m(i))
149 m2=irectm(2,cand_m(i))
150
151
152c calcul d'un minorant de la distance
153
154 xmax1 = max(x11,x21)
155 ymax1 = max(x12,x22)
156 zmax1 = max(x13,x23)
157 xmax2 = max(x(1,m1),x(1,m2))
158 ymax2 = max(x(2,m1),x(2,m2))
159 zmax2 = max(x(3,m1),x(3,m2))
160 xmin1 = min(x11,x21)
161 ymin1 = min(x12,x22)
162 zmin1 = min(x13,x23)
163 xmin2 = min(x(1,m1),x(1,m2))
164 ymin2 = min(x(2,m1),x(2,m2))
165 zmin2 = min(x(3,m1),x(3,m2))
166 dd = max(xmin1-xmax2,ymin1-ymax2,zmin1-zmax2,
167 . xmin2-xmax1,ymin2-ymax1,zmin2-zmax1)
168 IF(dd > gapv(i))THEN
169 pene(i) = zero
170 cycle
171 ENDIF
172
173c calcul de la distance^2
174
175 xs12 = x21-x11
176 ys12 = x22-x12
177 zs12 = x23-x13
178 xs2m2 = x(1,m2)-x21
179 ys2m2 = x(2,m2)-x22
180 zs2m2 = x(3,m2)-x23
181 xs2 = xs12*xs12 + ys12*ys12 + zs12*zs12
182 xm12 = x(1,m2)-x(1,m1)
183 ym12 = x(2,m2)-x(2,m1)
184 zm12 = x(3,m2)-x(3,m1)
185 xm2 = xm12*xm12 + ym12*ym12 + zm12*zm12
186 xsm = - (xs12*xm12 + ys12*ym12 + zs12*zm12)
187 xa = xs12*xs2m2 + ys12*ys2m2 + zs12*zs2m2
188 xb = -xm12*xs2m2 - ym12*ys2m2 - zm12*zs2m2
189 det = xm2*xs2 - xsm*xsm
190 det = max(em20,det)
191C
192 alm = (xa*xsm-xb*xs2) / det
193 xs2 = max(xs2,em20)
194 xm2 = max(xm2,em20)
195 alm=min(one,max(zero,alm))
196 als = -(xa + alm*xsm) / xs2
197 als=min(one,max(zero,als))
198 alm = -(xb + als*xsm) / xm2
199 alm=min(one,max(zero,alm))
200
201C PENE = GAP^2 - DIST^2 UTILISE POUR TESTER SI NON NUL
202
203 xx = als*x11 + (one-als)*x21
204 . - alm*x(1,m1) - (one-alm)*x(1,m2)
205 yy = als*x12 + (one-als)*x22
206 . - alm*x(2,m1) - (one-alm)*x(2,m2)
207 zz = als*x13 + (one-als)*x23
208 . - alm*x(3,m1) - (one-alm)*x(3,m2)
209 gap2=gapv(i)*gapv(i)
210 pene(i) = gap2- xx*xx - yy*yy - zz*zz
211C
212 END DO
213
214
215C
216 RETURN
217 END
218C===============================================================================
219!||====================================================================
220!|| i11pen3 ../engine/source/interfaces/intsort/i11pen3.F
221!||--- called by ------------------------------------------------------
222!|| i11sto ../engine/source/interfaces/intsort/i11sto.F
223!||--- uses -----------------------------------------------------
224!|| tri7box ../engine/share/modules/tri7box.F
225!||====================================================================
226 SUBROUTINE i11pen3(JLT ,CAND_S,CAND_M,GAP ,X ,
227 . IRECTS,IRECTM,PENE ,NRTS )
228C-----------------------------------------------
229C M o d u l e s
230C-----------------------------------------------
231 USE tri7box
232C-----------------------------------------------
233C I m p l i c i t T y p e s
234C-----------------------------------------------
235#include "implicit_f.inc"
236C-----------------------------------------------
237C G l o b a l P a r a m e t e r s
238C-----------------------------------------------
239#include "mvsiz_p.inc"
240C-----------------------------------------------
241C D u m m y A r g u m e n t s
242C-----------------------------------------------
243 INTEGER JLT, NRTS
244 INTEGER IRECTS(2,*), IRECTM(2,*),CAND_S(*),CAND_M(*)
245 my_real
246 . GAP
247 my_real
248 . x(3,*), pene(mvsiz)
249C-----------------------------------------------
250C L o c a l V a r i a b l e s
251C-----------------------------------------------
252 INTEGER I, IG,N1,N2,M1,M2,NI,L,J
253 my_real
254 . XS12,YS12,ZS12,XM12,YM12,ZM12,XA,XB,
255 . XS2,XM2,XSM,XS2M2,YS2,YM2,YSM,YS2M2,ZS2,ZM2,ZSM,ZS2M2,
256 . xx,yy,zz,als,alm,det,
257 . gap2, x11, x12, x13, x21, x22, x23,
258 . xmax1,ymax1,zmax1,xmax2,ymax2,zmax2,
259 . xmin1,ymin1,zmin1,xmin2,ymin2,zmin2,dd
260C-----------------------------------------------
261 gap2=gap*gap
262C--------------------------------------------------------
263C
264C--------------------------------------------------------
265C F = [A*X1+(1-A)*X2-B*X3-(1-B)*X4]^2 + [..Y..]^2 + [..Z..]^2
266C DF/DA = 0 = (X1-X2)(A(X1-X2)+X2-X4 +B(X4-X3))+...
267C DF/DA = 0 = A(X1-X2)^2 +X2-X4 + B(X1-X2)(X4-X3))+...
268C DF/DA = 0 = A[(X1-X2)^2 + (Y1-Y2)^2 + (Z1-Z2)^2]
269C + B[(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
270C + (X1-X2)(X2-X4) + (Y1-Y2)(Y2-Y4) + (Z1-Z2)(Z2-Z4)
271C DF/DB = 0 = (X4-X3)(A(X1-X2)+X2-X4 +B(X4-X3))+...
272C DF/DB = 0 = B[(X4-X3)^2 + (Y4-Y3)^2 + (Z4-Z3)^2]
273C + A[(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
274C + (X4-X3)(X2-X4) + (Y4-Y3)(Y2-Y4) + (Z4-Z3)(Z2-Z4)
275C XS2 = [(X1-X2)^2 + (Y1-Y2)^2 + (Z1-Z2)^2]
276C XM2 = [(X4-X3)^2 + (Y4-Y3)^2 + (Z4-Z3)^2]
277C XSM = [(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
278C XA = (X1-X2)(X2-X4) + (Y1-Y2)(Y2-Y4) + (Z1-Z2)(Z2-Z4)
279C XB = (X4-X3)(X2-X4) + (Y4-Y3)(Y2-Y4) + (Z4-Z3)(Z2-Z4)
280C A XS2 + B XSM + XA = 0
281C A XSM + B XM2 + XB = 0
282C
283C A = -(XA + B XSM)/XS2
284C -(XA + B XSM)*XSM + B XM2*XS2 + XB*XS2 = 0
285C -B XSM*XSM + B XM2*XS2 + XB*XS2-XA*XSM = 0
286C B*(XM2*XS2 - XSM*XSM) = -XB*XS2+XA*XSM
287C B = (XA*XSM-XB*XS2) / (XM2*XS2 - XSM*XSM)
288C A = (XB*XSM-XA*XM2) / (XM2*XS2 - XSM*XSM)
289C--------------------------------------------------------
290C
291 DO i=1,jlt
292 l = cand_s(i)
293 IF(l<=nrts) THEN
294 ni=0
295 n1=irects(1,cand_s(i))
296 n2=irects(2,cand_s(i))
297 x11 = x(1,n1)
298 x12 = x(2,n1)
299 x13 = x(3,n1)
300 x21 = x(1,n2)
301 x22 = x(2,n2)
302 x23 = x(3,n2)
303 ELSE
304 ni = l - nrts
305 x11 = xrem(1,ni)
306 x12 = xrem(2,ni)
307 x13 = xrem(3,ni)
308 x21 = xrem(8,ni)
309 x22 = xrem(9,ni)
310 x23 = xrem(10,ni)
311 END IF
312 m1=irectm(1,cand_m(i))
313 m2=irectm(2,cand_m(i))
314
315
316c calcul d'un minorant de la distance
317
318 xmax1 = max(x11,x21)
319 ymax1 = max(x12,x22)
320 zmax1 = max(x13,x23)
321 xmax2 = max(x(1,m1),x(1,m2))
322 ymax2 = max(x(2,m1),x(2,m2))
323 zmax2 = max(x(3,m1),x(3,m2))
324 xmin1 = min(x11,x21)
325 ymin1 = min(x12,x22)
326 zmin1 = min(x13,x23)
327 xmin2 = min(x(1,m1),x(1,m2))
328 ymin2 = min(x(2,m1),x(2,m2))
329 zmin2 = min(x(3,m1),x(3,m2))
330 dd = max(xmin1-xmax2,ymin1-ymax2,zmin1-zmax2,
331 . xmin2-xmax1,ymin2-ymax1,zmin2-zmax1)
332 IF(dd > gap)THEN
333 pene(i) = zero
334 cycle
335 ENDIF
336
337c calcul de la distance^2
338
339 xs12 = x21-x11
340 ys12 = x22-x12
341 zs12 = x23-x13
342 xs2m2 = x(1,m2)-x21
343 ys2m2 = x(2,m2)-x22
344 zs2m2 = x(3,m2)-x23
345 xs2 = xs12*xs12 + ys12*ys12 + zs12*zs12
346 xm12 = x(1,m2)-x(1,m1)
347 ym12 = x(2,m2)-x(2,m1)
348 zm12 = x(3,m2)-x(3,m1)
349 xm2 = xm12*xm12 + ym12*ym12 + zm12*zm12
350 xsm = - (xs12*xm12 + ys12*ym12 + zs12*zm12)
351 xa = xs12*xs2m2 + ys12*ys2m2 + zs12*zs2m2
352 xb = -xm12*xs2m2 - ym12*ys2m2 - zm12*zs2m2
353 det = xm2*xs2 - xsm*xsm
354 det = max(em20,det)
355C
356 alm = (xa*xsm-xb*xs2) / det
357 xs2 = max(xs2,em20)
358 xm2 = max(xm2,em20)
359 alm=min(one,max(zero,alm))
360 als = -(xa + alm*xsm) / xs2
361 als=min(one,max(zero,als))
362 alm = -(xb + als*xsm) / xm2
363 alm=min(one,max(zero,alm))
364
365C PENE = GAP^2 - DIST^2 UTILISE POUR TESTER SI NON NUL
366
367 xx = als*x11 + (one-als)*x21
368 . - alm*x(1,m1) - (one-alm)*x(1,m2)
369 yy = als*x12 + (one-als)*x22
370 . - alm*x(2,m1) - (one-alm)*x(2,m2)
371 zz = als*x13 + (one-als)*x23
372 . - alm*x(3,m1) - (one-alm)*x(3,m2)
373 pene(i) = gap2- xx*xx - yy*yy - zz*zz
374C
375 END DO
376C
377 RETURN
378 END
379C===============================================================================
i11pen3_vox
subroutine i11pen3_vox(jlt, cand_s, cand_m, gapmin, drad, marge, gap_s, gap_m, gap_s_l, gap_m_l, igap, x, irects, irectm, pene, nrts, dgapload)
Definition i11pen3.F:34
i11pen3
subroutine i11pen3(jlt, cand_s, cand_m, gap, x, irects, irectm, pene, nrts)
Definition i11pen3.F:228
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21