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_vox1 ../starter/source/interfaces/inter3d1/i11pen3.F
25!||--- called by ------------------------------------------------------
26!|| i11sto_vox1 ../starter/source/interfaces/inter3d1/i11sto.F
27!||====================================================================
28 SUBROUTINE i11pen3_vox1(JLT ,CAND_S ,CAND_M ,GAPMIN ,DRAD ,
29 . MARGE ,GAP_S ,GAP_M ,GAP_S_L,GAP_M_L,
30 . IGAP ,X ,IRECTS ,IRECTM ,PENE ,
31 . DGAPLOAD )
32C-----------------------------------------------
33C I m p l i c i t T y p e s
34C-----------------------------------------------
35#include "implicit_f.inc"
36C-----------------------------------------------
37C G l o b a l P a r a m e t e r s
38C-----------------------------------------------
39#include "mvsiz_p.inc"
40C-----------------------------------------------
41C D u m m y A r g u m e n t s
42C-----------------------------------------------
43 INTEGER JLT, IGAP
44 INTEGER IRECTS(2,*), IRECTM(2,*),CAND_S(*),CAND_M(*)
45 my_real
46 . gapmin, drad, marge
47 my_real , INTENT(IN) :: dgapload
48 my_real
49 . x(3,*), pene(mvsiz),
50 . gap_s(*), gap_m(*), gap_s_l(*), gap_m_l(*)
51C-----------------------------------------------
52C L o c a l V a r i a b l e s
53C-----------------------------------------------
54 INTEGER I, IG,N1,N2,M1,M2
55 my_real
56 . XS12,YS12,ZS12,XM12,YM12,ZM12,XA,XB,
57 . xs2,xm2,xsm,xs2m2,ys2,ym2,ysm,ys2m2,zs2,zm2,zsm,zs2m2,
58 . xx,yy,zz,als,alm,det,
59 . gap2, gapv(mvsiz)
60C-----------------------------------------------
61C
62 IF(igap==0)THEN
63 DO i=1,jlt
64 gapv(i)=max(drad,gapmin+dgapload)+marge
65 ENDDO
66 ELSE
67 DO i=1,jlt
68 gapv(i)=gap_s(cand_s(i))+gap_m(cand_m(i))
69 IF(igap == 3)
70 . gapv(i)=min(gap_s_l(cand_s(i))+gap_m_l(cand_m(i)),gapv(i))
71 gapv(i)=max(drad,max(gapmin,gapv(i))+dgapload)+marge
72 ENDDO
73 ENDIF
74C--------------------------------------------------------
75C
76C--------------------------------------------------------
77C F = [A*X1+(1-A)*X2-B*X3-(1-B)*X4]^2 + [..Y..]^2 + [..Z..]^2
78C DF/DA = 0 = (X1-X2)(A(X1-X2)+X2-X4 +B(X4-X3))+...
79C DF/DA = 0 = A(X1-X2)^2 +X2-X4 + B(X1-X2)(X4-X3))+...
80C DF/DA = 0 = A[(X1-X2)^2 + (Y1-Y2)^2 + (Z1-Z2)^2]
81C + B[(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
82C + (X1-X2)(X2-X4) + (Y1-Y2)(Y2-Y4) + (Z1-Z2)(Z2-Z4)
83C DF/DB = 0 = (X4-X3)(A(X1-X2)+X2-X4 +B(X4-X3))+...
84C DF/DB = 0 = B[(X4-X3)^2 + (Y4-Y3)^2 + (Z4-Z3)^2]
85C + A[(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
86C + (X4-X3)(X2-X4) + (Y4-Y3)(Y2-Y4) + (Z4-Z3)(Z2-Z4)
87C XS2 = [(X1-X2)^2 + (Y1-Y2)^2 + (Z1-Z2)^2]
88C XM2 = [(X4-X3)^2 + (Y4-Y3)^2 + (Z4-Z3)^2]
89C XSM = [(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
90C XA = (X1-X2)(X2-X4) + (Y1-Y2)(Y2-Y4) + (Z1-Z2)(Z2-Z4)
91C XB = (X4-X3)(X2-X4) + (Y4-Y3)(Y2-Y4) + (Z4-Z3)(Z2-Z4)
92C A XS2 + B XSM + XA = 0
93C A XSM + B XM2 + XB = 0
94C
95C A = -(XA + B XSM)/XS2
96C -(XA + B XSM)*XSM + B XM2*XS2 + XB*XS2 = 0
97C -B XSM*XSM + B XM2*XS2 + XB*XS2-XA*XSM = 0
98C B*(XM2*XS2 - XSM*XSM) = -XB*XS2+XA*XSM
99C B = (XA*XSM-XB*XS2) / (XM2*XS2 - XSM*XSM)
100C A = (XB*XSM-XA*XM2) / (XM2*XS2 - XSM*XSM)
101 DO i=1,jlt
102 n1=irects(1,cand_s(i))
103 n2=irects(2,cand_s(i))
104 m1=irectm(1,cand_m(i))
105 m2=irectm(2,cand_m(i))
106 xs12 = x(1,n2)-x(1,n1)
107 ys12 = x(2,n2)-x(2,n1)
108 zs12 = x(3,n2)-x(3,n1)
109 xs2 = xs12*xs12 + ys12*ys12 + zs12*zs12
110 xm12 = x(1,m2)-x(1,m1)
111 ym12 = x(2,m2)-x(2,m1)
112 zm12 = x(3,m2)-x(3,m1)
113 xm2 = xm12*xm12 + ym12*ym12 + zm12*zm12
114 xsm = - (xs12*xm12 + ys12*ym12 + zs12*zm12)
115 xs2m2 = x(1,m2)-x(1,n2)
116 ys2m2 = x(2,m2)-x(2,n2)
117 zs2m2 = x(3,m2)-x(3,n2)
118 xa = xs12*xs2m2 + ys12*ys2m2 + zs12*zs2m2
119 xb = -xm12*xs2m2 - ym12*ys2m2 - zm12*zs2m2
120 det = xm2*xs2 - xsm*xsm
121 det = max(em20,det)
122C
123 alm = (xa*xsm-xb*xs2) / det
124 xs2 = max(xs2,em20)
125 xm2 = max(xm2,em20)
126 alm=min(one,max(zero,alm))
127 als = -(xa + alm*xsm) / xs2
128 als=min(one,max(zero,als))
129 alm = -(xb + als*xsm) / xm2
130 alm=min(one,max(zero,alm))
131
132C !!!!!!!!!!!!!!!!!!!!!!!
133C PENE = GAP^2 - DIST^2 UTILISE POUR TESTER SI NON NUL
134C!!!!!!!!!!!!!!!!!!!!!!!!
135 xx = als*x(1,n1) + (1.-als)*x(1,n2)
136 . - alm*x(1,m1) - (1.-alm)*x(1,m2)
137 yy = als*x(2,n1) + (1.-als)*x(2,n2)
138 . - alm*x(2,m1) - (1.-alm)*x(2,m2)
139 zz = als*x(3,n1) + (1.-als)*x(3,n2)
140 . - alm*x(3,m1) - (1.-alm)*x(3,m2)
141 gap2=gapv(i)*gapv(i)
142 pene(i) = max(zero,gap2- xx*xx - yy*yy - zz*zz)
143C
144 ENDDO
145C
146 RETURN
147 END
148!||====================================================================
149!|| i11pen3 ../starter/source/interfaces/inter3d1/i11pen3.F
150!||--- called by ------------------------------------------------------
151!|| i11sto ../starter/source/interfaces/inter3d1/i11sto.F
152!||====================================================================
153 SUBROUTINE i11pen3(JLT ,CAND_N,CAND_E,GAP ,X,
154 . IRECTS ,IRECTM ,PENE )
155C-----------------------------------------------
156C I m p l i c i t T y p e s
157C-----------------------------------------------
158#include "implicit_f.inc"
159C-----------------------------------------------
160C G l o b a l P a r a m e t e r s
161C-----------------------------------------------
162#include "mvsiz_p.inc"
163C-----------------------------------------------
164C D u m m y A r g u m e n t s
165C-----------------------------------------------
166 INTEGER JLT
167 INTEGER IRECTS(2,*), IRECTM(2,*),CAND_N(*),CAND_E(*)
168 my_real
169 . GAP
170 my_real
171 . x(3,*), pene(mvsiz)
172C-----------------------------------------------
173C L o c a l V a r i a b l e s
174C-----------------------------------------------
175 INTEGER I, IG,N1,N2,M1,M2
176 my_real
177 . XS12,YS12,ZS12,XM12,YM12,ZM12,XA,XB,
178 . XS2,XM2,XSM,XS2M2,YS2,YM2,YSM,YS2M2,ZS2,ZM2,ZSM,ZS2M2,
179 . xx,yy,zz,als,alm,det,
180 . gap2
181C-----------------------------------------------
182 gap2=gap*gap
183C
184C--------------------------------------------------------
185C
186C--------------------------------------------------------
187C F = [A*X1+(1-A)*X2-B*X3-(1-B)*X4]^2 + [..Y..]^2 + [..Z..]^2
188C DF/DA = 0 = (X1-X2)(A(X1-X2)+X2-X4 +B(X4-X3))+...
189C DF/DA = 0 = A(X1-X2)^2 +X2-X4 + B(X1-X2)(X4-X3))+...
190C DF/DA = 0 = A[(X1-X2)^2 + (Y1-Y2)^2 + (Z1-Z2)^2]
191C + B[(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
192C + (X1-X2)(X2-X4) + (Y1-Y2)(Y2-Y4) + (Z1-Z2)(Z2-Z4)
193C DF/DB = 0 = (X4-X3)(A(X1-X2)+X2-X4 +B(X4-X3))+...
194C DF/DB = 0 = B[(X4-X3)^2 + (Y4-Y3)^2 + (Z4-Z3)^2]
195C + A[(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
196C + (X4-X3)(X2-X4) + (Y4-Y3)(Y2-Y4) + (Z4-Z3)(Z2-Z4)
197C XS2 = [(X1-X2)^2 + (Y1-Y2)^2 + (Z1-Z2)^2]
198C XM2 = [(X4-X3)^2 + (Y4-Y3)^2 + (Z4-Z3)^2]
199C XSM = [(X1-X2)(X4-X3) + (Y1-Y2)(Y4-Y3) + (Z1-Z2)(Z4-Z3)]
200C XA = (X1-X2)(X2-X4) + (Y1-Y2)(Y2-Y4) + (Z1-Z2)(Z2-Z4)
201C XB = (X4-X3)(X2-X4) + (Y4-Y3)(Y2-Y4) + (Z4-Z3)(Z2-Z4)
202C A XS2 + B XSM + XA = 0
203C A XSM + B XM2 + XB = 0
204C
205C A = -(XA + B XSM)/XS2
206C -(XA + B XSM)*XSM + B XM2*XS2 + XB*XS2 = 0
207C -B XSM*XSM + B XM2*XS2 + XB*XS2-XA*XSM = 0
208C B*(XM2*XS2 - XSM*XSM) = -XB*XS2+XA*XSM
209C B = (XA*XSM-XB*XS2) / (XM2*XS2 - XSM*XSM)
210C A = (XB*XSM-XA*XM2) / (XM2*XS2 - XSM*XSM)
211 DO i=1,jlt
212 n1=irects(1,cand_n(i))
213 n2=irects(2,cand_n(i))
214 m1=irectm(1,cand_e(i))
215 m2=irectm(2,cand_e(i))
216 xs12 = x(1,n2)-x(1,n1)
217 ys12 = x(2,n2)-x(2,n1)
218 zs12 = x(3,n2)-x(3,n1)
219 xs2 = xs12*xs12 + ys12*ys12 + zs12*zs12
220 xm12 = x(1,m2)-x(1,m1)
221 ym12 = x(2,m2)-x(2,m1)
222 zm12 = x(3,m2)-x(3,m1)
223 xm2 = xm12*xm12 + ym12*ym12 + zm12*zm12
224 xsm = - (xs12*xm12 + ys12*ym12 + zs12*zm12)
225 xs2m2 = x(1,m2)-x(1,n2)
226 ys2m2 = x(2,m2)-x(2,n2)
227 zs2m2 = x(3,m2)-x(3,n2)
228 xa = xs12*xs2m2 + ys12*ys2m2 + zs12*zs2m2
229 xb = -xm12*xs2m2 - ym12*ys2m2 - zm12*zs2m2
230 det = xm2*xs2 - xsm*xsm
231 det = max(em20,det)
232C
233 alm = (xa*xsm-xb*xs2) / det
234 xs2 = max(xs2,em20)
235 xm2 = max(xm2,em20)
236 alm=min(one,max(zero,alm))
237 als = -(xa + alm*xsm) / xs2
238 als=min(one,max(zero,als))
239 alm = -(xb + als*xsm) / xm2
240 alm=min(one,max(zero,alm))
241
242C !!!!!!!!!!!!!!!!!!!!!!!
243C PENE = GAP^2 - DIST^2 UTILISE POUR TESTER SI NON NUL
244C!!!!!!!!!!!!!!!!!!!!!!!!
245 xx = als*x(1,n1) + (1.-als)*x(1,n2)
246 . - alm*x(1,m1) - (1.-alm)*x(1,m2)
247 yy = als*x(2,n1) + (1.-als)*x(2,n2)
248 . - alm*x(2,m1) - (1.-alm)*x(2,m2)
249 zz = als*x(3,n1) + (1.-als)*x(3,n2)
250 . - alm*x(3,m1) - (1.-alm)*x(3,m2)
251 pene(i) = max(zero,gap2- xx*xx - yy*yy - zz*zz)
252C
253 ENDDO
254C
255 RETURN
256 END
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
i11pen3
subroutine i11pen3(jlt, cand_n, cand_e, gap, x, irects, irectm, pene)
Definition i11pen3.F:155
i11pen3_vox1
subroutine i11pen3_vox1(jlt, cand_s, cand_m, gapmin, drad, marge, gap_s, gap_m, gap_s_l, gap_m_l, igap, x, irects, irectm, pene, dgapload)
Definition i11pen3.F:32