OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
i25cor3e.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!|| i25cor3e ../starter/source/interfaces/inter3d1/i25cor3e.F
25!||--- called by ------------------------------------------------------
26!|| inint3 ../starter/source/interfaces/inter3d1/inint3.F
27!||--- uses -----------------------------------------------------
28!|| tri7box ../starter/share/modules1/tri7box.F
29!||====================================================================
30 SUBROUTINE i25cor3e(
31 1 JLT ,LEDGE ,IRECT ,X ,
32 2 CAND_S ,CAND_M ,
33 4 XXS1 ,XXS2 ,XYS1 ,XYS2 ,XZS1 ,
34 5 XZS2 ,XXM1 ,XXM2 ,XYM1 ,XYM2 ,
35 6 XZM1 ,XZM2 ,EX ,EY ,EZ ,
36 7 FX ,FY ,FZ ,
37 8 N1 ,N2 ,M1 ,M2 ,NEDGE ,
38 9 GAPE ,GAPVE ,
39 A IEDGE ,ADMSR ,LBOUND ,EDG_BISECTOR ,
40 B VTX_BISECTOR ,ITAB ,IGAP0 ,IGAP ,
41 C GAP_E_L )
42C-----------------------------------------------
43C M o d u l e s
44C-----------------------------------------------
45 USE tri7box
46C-----------------------------------------------
47C I m p l i c i t T y p e s
48C-----------------------------------------------
49#include "implicit_f.inc"
50C-----------------------------------------------
51C G l o b a l P a r a m e t e r s
52C-----------------------------------------------
53#include "mvsiz_p.inc"
54C-----------------------------------------------
55C C o m m o n B l o c k s
56C-----------------------------------------------
57#include "param_c.inc"
58C-----------------------------------------------
59C D u m m y A r g u m e n t s
60C-----------------------------------------------
61 INTEGER LEDGE(NLEDGE,*), IRECT(4,*), CAND_M(*), CAND_S(*), ADMSR(4,*),
62 . LBOUND(*), JLT, NEDGE, IEDGE, ITAB(*), IGAP0, IGAP,
63 . N1(MVSIZ), N2(MVSIZ),
64 . M1(MVSIZ), M2(MVSIZ)
65C REAL
66 my_real
67 . X(3,*),
68 . XXS1(MVSIZ), XXS2(MVSIZ), XYS1(MVSIZ), XYS2(MVSIZ),
69 . XZS1(MVSIZ), XZS2(MVSIZ), XXM1(MVSIZ), XXM2(MVSIZ),
70 . XYM1(MVSIZ), XYM2(MVSIZ), XZM1(MVSIZ), XZM2(MVSIZ),
71 . GAPE(*),GAPVE(MVSIZ),GAP_E_L(NEDGE),
72 . ex(mvsiz), ey(mvsiz), ez(mvsiz), fx(mvsiz), fy(mvsiz), fz(mvsiz)
73 real*4 edg_bisector(3,4,*), vtx_bisector(3,2,*)
74C-----------------------------------------------
75C L o c a l V a r i a b l e s
76C-----------------------------------------------
77 INTEGER I ,NN, J, JRM, K, KRM, I1, J1, I2, J2,
78 . IE, JE, SOL_EDGE, SH_EDGE, IM(MVSIZ), IS(MVSIZ)
79 my_real
80 . AAA, DX, DY, DZ, DD, NNI, NI2, INVCOS, GAPE_M(MVSIZ), GAPE_S(MVSIZ)
81C-----------------------------------------------
82 DO I=1,jlt
83 IF(cand_s(i)<=nedge) THEN
84
85 i1=ledge(1,cand_s(i))
86 j1=ledge(2,cand_s(i))
87 n1(i)=ledge(5,cand_s(i))
88 n2(i)=ledge(6,cand_s(i))
89
90 IF(irect(j1,i1)==n1(i).AND.irect(mod(j1,4)+1,i1)==n2(i))THEN
91 is(i)= 1
92 ELSEIF(irect(j1,i1)==n2(i).AND.irect(mod(j1,4)+1,i1)==n1(i))THEN
93 is(i)=-1
94 ELSE
95 print *,'i25cor3e - internal problem',cand_s(i),n1(i),n2(i),
96 . irect(j1,i1),irect(mod(j1,4)+1,i1)
97 END IF
98
99 i2=ledge(1,cand_m(i))
100 j2=ledge(2,cand_m(i))
101 m1(i)=ledge(5,cand_m(i))
102 m2(i)=ledge(6,cand_m(i))
103
104 IF(irect(j2,i2)==m1(i).AND.irect(mod(j2,4)+1,i2)==m2(i))THEN
105 im(i)= 1
106 ELSEIF(irect(j2,i2)==m2(i).AND.irect(mod(j2,4)+1,i2)==m1(i))THEN
107 im(i)=-1
108 ELSE
109 print *,'i25cor3e - internal problem',cand_m(i),m1(i),m2(i),
110 . irect(j2,i2),irect(mod(j2,4)+1,i2)
111 END IF
112
113 xxs1(i) = x(1,n1(i))
114 xys1(i) = x(2,n1(i))
115 xzs1(i) = x(3,n1(i))
116 xxs2(i) = x(1,n2(i))
117 xys2(i) = x(2,n2(i))
118 xzs2(i) = x(3,n2(i))
119 xxm1(i) = x(1,m1(i))
120 xym1(i) = x(2,m1(i))
121 xzm1(i) = x(3,m1(i))
122 xxm2(i) = x(1,m2(i))
123 xym2(i) = x(2,m2(i))
124 xzm2(i) = x(3,m2(i))
125 END IF
126 END DO
127
128 DO i=1,jlt
129 gape_m(i)=gape(cand_m(i))
130 IF(cand_s(i)<=nedge) THEN
131 gape_s(i)=gape(cand_s(i))
132 ELSE ! TBD
133 END IF
134 gapve(i)=gape_m(i)+gape_s(i)
135 END DO
136 IF(igap==3) THEN
137 DO i=1,jlt
138 gape_m(i)=min(gape_m(i),gap_e_l(cand_m(i)))
139 IF(cand_s(i)<=nedge) THEN
140 gape_s(i)=min(gape_s(i),gap_e_l(cand_s(i)))
141 gapve(i)=min(gap_e_l(cand_m(i))+gap_e_l(cand_s(i)),gapve(i))
142 ENDIF
143 ENDDO
144 ENDIF
145C
146 sol_edge=iedge/10 ! solids
147 sh_edge =iedge-10*sol_edge ! shells
148
149 ex(1:jlt)=zero
150 ey(1:jlt)=zero
151 ez(1:jlt)=zero
152 fx(1:jlt)=zero
153 fy(1:jlt)=zero
154 fz(1:jlt)=zero
155 IF(sh_edge/=0)THEN
156 DO i=1,jlt
157C
158 ie=ledge(1,cand_m(i))
159 je=ledge(2,cand_m(i))
160 ex(i) = edg_bisector(1,je,ie)
161 ey(i) = edg_bisector(2,je,ie)
162 ez(i) = edg_bisector(3,je,ie)
163C
164 IF(iabs(ledge(7,cand_s(i)))/=1)THEN
165 ie=ledge(1,cand_s(i))
166 je=ledge(2,cand_s(i))
167 fx(i) = edg_bisector(1,je,ie)
168 fy(i) = edg_bisector(2,je,ie)
169 fz(i) = edg_bisector(3,je,ie)
170 END IF
171
172 END DO
173 END IF
174
175 DO i=1,jlt
176C
177 IF(ledge(3,cand_m(i))/=0)cycle
178C
179 ie=ledge(1,cand_m(i))
180 je=ledge(2,cand_m(i))
181 IF(im(i)==1)THEN
182 i1=admsr(je,ie)
183 i2=admsr(mod(je,4)+1,ie)
184 ELSE ! IM(I)==-1
185 i2=admsr(je,ie)
186 i1=admsr(mod(je,4)+1,ie)
187 END IF
188C
189 ex(i) = edg_bisector(1,je,ie)
190 ey(i) = edg_bisector(2,je,ie)
191 ez(i) = edg_bisector(3,je,ie)
192C
193C DX, DY, DZ Direction for shifting
194 dx = vtx_bisector(1,1,i1)+vtx_bisector(1,2,i1)
195 dy = vtx_bisector(2,1,i1)+vtx_bisector(2,2,i1)
196 dz = vtx_bisector(3,1,i1)+vtx_bisector(3,2,i1)
197C
198 nni = ex(i)*dx + ey(i)*dy + ez(i)*dz
199 ni2 = dx*dx + dy*dy + dz*dz
200
201 IF(nni < zero)THEN
202 dx=dx-two*nni*ex(i)
203 dy=dy-two*nni*ey(i)
204 dz=dz-two*nni*ez(i)
205 nni=-nni
206 END IF
207
208 IF(two*nni*nni < ni2)THEN
209c scharp angle bound nodal normal to 45 from edge normal
210 aaa = sqrt(max(zero,ni2-nni*nni)) - nni
211 dx = dx + aaa*ex(i)
212 dy = dy + aaa*ey(i)
213 dz = dz + aaa*ez(i)
214 ENDIF
215 dd=one/max(em20,sqrt(dx*dx+dy*dy+dz*dz))
216 dx = dx*dd
217 dy = dy*dd
218 dz = dz*dd
219 invcos = one / max(em20,ex(i)*dx + ey(i)*dy + ez(i)*dz)
220 dx = dx*invcos
221 dy = dy*invcos
222 dz = dz*invcos
223C
224 xxm1(i) = xxm1(i)-gape_m(i)*dx
225 xym1(i) = xym1(i)-gape_m(i)*dy
226 xzm1(i) = xzm1(i)-gape_m(i)*dz
227C
228C DX, DY, DZ Direction for shifting
229 dx = vtx_bisector(1,1,i2)+vtx_bisector(1,2,i2)
230 dy = vtx_bisector(2,1,i2)+vtx_bisector(2,2,i2)
231 dz = vtx_bisector(3,1,i2)+vtx_bisector(3,2,i2)
232C
233 nni = ex(i)*dx + ey(i)*dy + ez(i)*dz
234 ni2 = dx*dx + dy*dy + dz*dz
235
236 IF(nni < zero)THEN
237 dx=dx-two*nni*ex(i)
238 dy=dy-two*nni*ey(i)
239 dz=dz-two*nni*ez(i)
240 nni=-nni
241 END IF
242
243 IF(two*nni*nni < ni2)THEN
244c scharp angle bound nodal normal to 45 from edge normal
245 aaa = sqrt(max(zero,ni2-nni*nni)) - nni
246 dx = dx + aaa*ex(i)
247 dy = dy + aaa*ey(i)
248 dz = dz + aaa*ez(i)
249 ENDIF
250 dd=one/max(em20,sqrt(dx*dx+dy*dy+dz*dz))
251 dx = dx*dd
252 dy = dy*dd
253 dz = dz*dd
254 invcos = one / max(em20,ex(i)*dx + ey(i)*dy + ez(i)*dz)
255 dx = dx*invcos
256 dy = dy*invcos
257 dz = dz*invcos
258C
259 xxm2(i) = xxm2(i)-gape_m(i)*dx
260 xym2(i) = xym2(i)-gape_m(i)*dy
261 xzm2(i) = xzm2(i)-gape_m(i)*dz
262C
263 END DO
264
265 IF(igap0/=0) THEN
266 DO i=1,jlt
267C
268 IF(ledge(3,cand_s(i))/=0)cycle
269C
270 ie=ledge(1,cand_s(i))
271 je=ledge(2,cand_s(i))
272 IF(is(i)==1)THEN
273 i1=admsr(je,ie)
274 i2=admsr(mod(je,4)+1,ie)
275 ELSE ! IS(I)==-1
276 i2=admsr(je,ie)
277 i1=admsr(mod(je,4)+1,ie)
278 END IF
279C
280 fx(i) = edg_bisector(1,je,ie)
281 fy(i) = edg_bisector(2,je,ie)
282 fz(i) = edg_bisector(3,je,ie)
283C
284C DX, DY, DZ Direction for shifting
285 dx = vtx_bisector(1,1,i1)+vtx_bisector(1,2,i1)
286 dy = vtx_bisector(2,1,i1)+vtx_bisector(2,2,i1)
287 dz = vtx_bisector(3,1,i1)+vtx_bisector(3,2,i1)
288C
289 nni = fx(i)*dx + fy(i)*dy + fz(i)*dz
290 ni2 = dx*dx + dy*dy + dz*dz
291
292 IF(nni < zero)THEN
293 dx=dx-two*nni*fx(i)
294 dy=dy-two*nni*fy(i)
295 dz=dz-two*nni*fz(i)
296 nni=-nni
297 END IF
298
299 IF(two*nni*nni < ni2)THEN
300c scharp angle bound nodal normal to 45 from edge normal
301 aaa = sqrt(max(zero,ni2-nni*nni)) - nni
302 dx = dx + aaa*fx(i)
303 dy = dy + aaa*fy(i)
304 dz = dz + aaa*fz(i)
305 ENDIF
306 dd=one/max(em20,sqrt(dx*dx+dy*dy+dz*dz))
307 dx = dx*dd
308 dy = dy*dd
309 dz = dz*dd
310 invcos = one / max(em20,fx(i)*dx + fy(i)*dy + fz(i)*dz)
311 dx = dx*invcos
312 dy = dy*invcos
313 dz = dz*invcos
314C
315 xxs1(i) = xxs1(i)-gape_s(i)*dx
316 xys1(i) = xys1(i)-gape_s(i)*dy
317 xzs1(i) = xzs1(i)-gape_s(i)*dz
318C
319C DX, DY, DZ Direction for shifting
320 dx = vtx_bisector(1,1,i2)+vtx_bisector(1,2,i2)
321 dy = vtx_bisector(2,1,i2)+vtx_bisector(2,2,i2)
322 dz = vtx_bisector(3,1,i2)+vtx_bisector(3,2,i2)
323C
324 nni = fx(i)*dx + fy(i)*dy + fz(i)*dz
325 ni2 = dx*dx + dy*dy + dz*dz
326
327 IF(nni < zero)THEN
328 dx=dx-two*nni*fx(i)
329 dy=dy-two*nni*fy(i)
330 dz=dz-two*nni*fz(i)
331 nni=-nni
332 END IF
333
334 IF(two*nni*nni < ni2)THEN
335c scharp angle bound nodal normal to 45 from edge normal
336 aaa = sqrt(max(zero,ni2-nni*nni)) - nni
337 dx = dx + aaa*fx(i)
338 dy = dy + aaa*fy(i)
339 dz = dz + aaa*fz(i)
340 ENDIF
341 dd=one/max(em20,sqrt(dx*dx+dy*dy+dz*dz))
342 dx = dx*dd
343 dy = dy*dd
344 dz = dz*dd
345 invcos = one / max(em20,fx(i)*dx + fy(i)*dy + fz(i)*dz)
346 dx = dx*invcos
347 dy = dy*invcos
348 dz = dz*invcos
349C
350 xxs2(i) = xxs2(i)-gape_s(i)*dx
351 xys2(i) = xys2(i)-gape_s(i)*dy
352 xzs2(i) = xzs2(i)-gape_s(i)*dz
353C
354 END DO
355 END IF
356C
357 RETURN
358 END
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
i25cor3e
subroutine i25cor3e(jlt, ledge, irect, x, cand_s, cand_m, xxs1, xxs2, xys1, xys2, xzs1, xzs2, xxm1, xxm2, xym1, xym2, xzm1, xzm2, ex, ey, ez, fx, fy, fz, n1, n2, m1, m2, nedge, gape, gapve, iedge, admsr, lbound, edg_bisector, vtx_bisector, itab, igap0, igap, gap_e_l)
Definition i25cor3e.F:42