OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
func_inters.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!|| func_inters ../starter/source/tools/curve/func_inters.F
25!||--- called by ------------------------------------------------------
26!|| law36_upd ../starter/source/materials/mat/mat036/law36_upd.F
27!|| law58_upd ../starter/source/materials/mat/mat058/law58_upd.F
28!|| law88_upd ../starter/source/materials/mat/mat088/law88_upd.F
29!||--- uses -----------------------------------------------------
30!||====================================================================
31 SUBROUTINE func_inters(TITR,MAT_ID,FUNC1,FUNC2,FAC1,FAC2,NPC,PLD,XINT,YINT)
32C-----------------------------------------------
33C M o d u l e s
34C-----------------------------------------------
35 USE names_and_titles_mod , ONLY : nchartitle
36C-----------------------------------------------
37C I m p l i c i t T y p e s
38C-----------------------------------------------
39#include "implicit_f.inc"
40C-----------------------------------------------
41C C o m m o n B l o c k s
42C-----------------------------------------------
43#include "tabsiz_c.inc"
44C-----------------------------------------------
45C D u m m y A r g u m e n t s
46C-----------------------------------------------
47 CHARACTER(LEN=NCHARTITLE) :: TITR
48 INTEGER :: MAT_ID,FUNC1,FUNC2
49 my_real :: xint,yint,fac1,fac2
50 INTEGER ,DIMENSION(SNPC) :: NPC
51 my_real ,DIMENSION(STF) :: pld
52C-----------------------------------------------
53 INTENT(IN) :: titr,func1,func2,mat_id,npc,pld,fac1,fac2
54 INTENT(INOUT) :: xint,yint
55C-----------------------------------------------
56C L o c a l V a r i a b l e s
57C-----------------------------------------------
58 INTEGER :: I,J,ID,NP1,NP2,J1,K,K1,FOUND
59 my_real :: s1,s2,t1,t2,x1,x2,y1,y2,ax,bx,ay,by,cx,cy,dm,alpha,beta
60C=======================================================================
61c Check common points between 2 curves => intersection
62 xint = zero
63 yint = zero
64 found = 0
65 np1 = (npc(func1+1)-npc(func1)) / 2 ! nb of points of the 1st curve
66 np2 = (npc(func2+1)-npc(func2)) / 2 ! nb of points of the 2nd curve
67 DO j=2,np1+1
68 j1=2*(j-2)
69 s1=pld(npc(func1)+j1)
70 t1=pld(npc(func1)+j1+1)*fac1
71 DO k=2,np2+1
72 k1=2*(k-2)
73 x1=pld(npc(func2)+k1)
74 y1=pld(npc(func2)+k1+1)*fac2
75 IF (x1 == s1 .AND. y1 == t1 .AND. x1> zero) THEN
76 xint = x1
77 yint = y1
78 found = 1
79 EXIT
80 ENDIF
81 ENDDO
82 IF (found == 1) EXIT
83 ENDDO
84c Check intersection of curve segments
85 IF (found == 0) THEN
86 DO j=2,np1
87 j1=2*(j-2)
88 s1=pld(npc(func1)+j1)
89 s2=pld(npc(func1)+j1+2)
90 t1=pld(npc(func1)+j1+1)*fac1
91 t2=pld(npc(func1)+j1+3)*fac1
92 DO k=2,np2
93 k1=2*(k-2)
94 x1=pld(npc(func2)+k1)
95 x2=pld(npc(func2)+k1+2)
96 y1=pld(npc(func2)+k1+1)*fac2
97 y2=pld(npc(func2)+k1+3)*fac2
98 IF (x2 < s1 .or. s2 < x1) cycle
99 ax = x2 - x1
100 ay = y2 - y1
101 bx = s1 - s2
102 by = t1 - t2
103 dm = ay*bx - ax*by
104 IF (dm /= zero) THEN ! check if segments are not parallel
105 cx = s1 - x1
106 cy = t1 - y1
107 alpha = (bx * cy - by * cx) / dm
108 beta = (ax * cy - ay * cx) / dm
109 IF (alpha > zero .and. alpha < one .and.
110 . beta < zero .and. beta >-one) THEN
111 xint = x1 + alpha * ax
112 yint = y1 + alpha * ay
113 found = 1
114 EXIT
115 ENDIF
116 ENDIF
117 ENDDO
118 IF (found == 1) EXIT
119 ENDDO
120 END IF
121c-----------
122 RETURN
123 END
124!||====================================================================
125!|| func_inters_shear ../starter/source/tools/curve/func_inters.F
126!||--- called by ------------------------------------------------------
127!|| law58_upd ../starter/source/materials/mat/mat058/law58_upd.F
128!||--- uses -----------------------------------------------------
129!|| message_mod ../starter/share/message_module/message_mod.F
130!|| table_mod ../starter/share/modules1/table_mod.f
131!||====================================================================
132 SUBROUTINE func_inters_shear(TITR,MAT_ID ,FUNC ,FUND, FAC1 ,FAC2,
133 . NPC,PLD,XINT1 ,YINT1 ,XINT2 ,YINT2 )
134C-----------------------------------------------
135C M o d u l e s
136C-----------------------------------------------
137 USE message_mod
138 USE table_mod
139 USE names_and_titles_mod , ONLY : nchartitle
140C-----------------------------------------------
141C I m p l i c i t T y p e s
142C-----------------------------------------------
143#include "implicit_f.inc"
144C-----------------------------------------------
145C D u m m y A r g u m e n t s
146C-----------------------------------------------
147 CHARACTER(LEN=NCHARTITLE) :: TITR
148 !INTEGER , DIMENSION(NFUNC) :: IFUNC,FUNC_ID
149 INTEGER FUNC,FUND,NPC(*)
150 INTEGER :: MAT_ID
151 my_real
152 . xint1 ,yint1 ,xint2 ,yint2,fac1,fac2,pld(*)
153! TYPE(TTABLE) TABLE(*)
154C-----------------------------------------------
155 INTENT(IN) :: titr,func,fund,mat_id,npc,pld,fac1,fac2
156 INTENT(INOUT) :: xint1 ,yint1 ,xint2 ,yint2
157C-----------------------------------------------
158C L o c a l V a r i a b l e s
159C-----------------------------------------------
160 INTEGER I,J,ID,NP1,NP2,J1,K,K1
161 my_real
162 . s1,s2,t1,t2,x1,x2,y1,y2,sx,ty,dydx,dtds
163C=======================================================================
164 !direction SHEAR, kfunc(3) kfunc(6)
165 xint1 = zero
166 yint1 = zero
167 xint2 = zero
168 yint2 = zero
169 sx = zero
170 ty = zero
171 np1 = (npc(func+1)-npc(func)) / 2
172 np2 = (npc(fund+1)-npc(fund)) / 2
173 DO j=2,np1+1
174 j1=2*(j-2)
175 s1=pld(npc(func)+j1)
176 t1=pld(npc(func)+j1+1)*fac1
177 DO k=2,np2+1
178 k1=2*(k-2)
179 x1=pld(npc(fund)+k1)
180 y1=pld(npc(fund)+k1+1)*fac2
181 IF(x1 == s1 .AND. y1 == t1 .AND.x1> zero)THEN
182 xint1 = x1
183 yint1 = y1
184 GOTO 361
185 ENDIF
186 ENDDO
187 ENDDO
188 DO j=2,np1
189 j1=2*(j-2)
190 s1=pld(npc(func)+j1)
191 s2=pld(npc(func)+j1+2)
192 t1=pld(npc(func)+j1+1)*fac1
193 t2=pld(npc(func)+j1+3)*fac1
194 DO k=2,np2
195 k1=2*(k-2)
196 x1=pld(npc(fund)+k1)
197 x2=pld(npc(fund)+k1+2)
198 y1=pld(npc(fund)+k1+1)*fac2
199 y2=pld(npc(fund)+k1+3)*fac2
200 IF(x1>zero.AND.x2>zero.AND.s1>zero.AND.s2>zero)THEN
201 IF (y2>=t1 .AND. y1<=t2 .AND. x2>=s1 .AND. x1<=s2) THEN
202 dydx = (y2-y1) / (x2-x1)
203 dtds = (t2-t1) / (s2-s1)
204 IF (dydx > dtds) THEN
205 !SX = (DTDS*S2-DYDX*X2-T2+Y2) / (DTDS-DYDX)
206 !TY = T2 + DTDS*(SX - S2)
207 sx = (t1-y1-dtds*s1+dydx*x1) / (dydx-dtds)
208 ty = t1 + dtds*(sx - s1)
209 IF (ty>=y1 .AND. ty<=y2 .AND. sx>=x1 .AND. sx<=x2.AND.sx/=zero)THEN
210 xint1 = sx
211 yint1 = ty
212 GOTO 361
213 ENDIF
214 ENDIF
215 ENDIF
216 ENDIF
217 ENDDO
218 ENDDO
219 361 CONTINUE
220 DO j=2,np1+1
221 j1=2*(j-2)
222 s1=pld(npc(func)+j1)
223 t1=pld(npc(func)+j1+1)*fac1
224 DO k=2,np2+1
225 k1=2*(k-2)
226 x1=pld(npc(fund)+k1)
227 y1=pld(npc(fund)+k1+1)*fac2
228 IF(x1 == s1 .AND. y1 == t1 .AND.x1 < zero)THEN
229 xint2 = x1
230 yint2 = y1
231 GOTO 362
232 ENDIF
233 ENDDO
234 ENDDO
235 DO j=2,np1
236 j1=2*(j-2)
237 s1=pld(npc(func)+j1)
238 s2=pld(npc(func)+j1+2)
239 t1=pld(npc(func)+j1+1)*fac1
240 t2=pld(npc(func)+j1+3)*fac1
241 DO k=2,np2
242 k1=2*(k-2)
243 x1=pld(npc(fund)+k1)
244 x2=pld(npc(fund)+k1+2)
245 y1=pld(npc(fund)+k1+1)*fac2
246 y2=pld(npc(fund)+k1+3)*fac2
247 IF(x1<zero.AND.x2<zero.AND.s1<zero.AND.s2<zero)THEN
248 IF (y2>=t1 .AND. y1<=t2 .AND. x2>=s1 .AND. x1<=s2) THEN
249 dydx = (y2-y1) / (x2-x1)
250 dtds = (t2-t1) / (s2-s1)
251 IF (dydx > dtds) THEN
252 !SX = (DTDS*S2-DYDX*X2-T2+Y2) / (DTDS-DYDX)
253 !TY = T2 + DTDS*(SX - S2)
254 sx = (t1-y1-dtds*s1+dydx*x1) / (dydx-dtds)
255 ty = t1 + dtds*(sx - s1)
256 IF (ty>=y1 .AND. ty<=y2 .AND. sx>=x1 .AND. sx<=x2.AND.sx/=zero)THEN
257 xint2 = sx
258 yint2 = ty
259 GOTO 362
260 ENDIF
261 ENDIF
262 ENDIF
263 ENDIF
264 ENDDO
265 ENDDO
266 362 CONTINUE
267c-----------
268 RETURN
269 END
270!||====================================================================
271!|| func_inters_c ../starter/source/tools/curve/func_inters.F
272!||--- called by ------------------------------------------------------
273!|| law88_upd ../starter/source/materials/mat/mat088/law88_upd.f
274!||--- uses -----------------------------------------------------
275!|| message_mod ../starter/share/message_module/message_mod.F
276!|| table_mod ../starter/share/modules1/table_mod.F
277!||====================================================================
278 SUBROUTINE func_inters_c(TITR,MAT_ID ,FUNC,FUND,FAC1,FAC2,NPC,PLD,XINC,YINC )
279C-----------------------------------------------
280C M o d u l e s
281C-----------------------------------------------
282 USE message_mod
283 USE table_mod
284 USE names_and_titles_mod , ONLY : nchartitle
285C-----------------------------------------------
286C I m p l i c i t T y p e s
287C-----------------------------------------------
288#include "implicit_f.inc"
289C-----------------------------------------------
290C D u m m y A r g u m e n t s
291C-----------------------------------------------
292 CHARACTER(LEN=NCHARTITLE) :: TITR
293 INTEGER FUNC,FUND,NPC(*)
294 INTEGER :: MAT_ID
295 my_real xinc,yinc,fac1,fac2,pld(*)
296C-----------------------------------------------
297 INTENT(IN) :: titr,func,fund,mat_id,npc,pld,fac1,fac2
298 INTENT(INOUT) :: xinc,yinc
299C-----------------------------------------------
300C L o c a l V a r i a b l e s
301C-----------------------------------------------
302 INTEGER I,J,ID,NP1,NP2,J1,K,K1
303 my_real
304 . s1,s2,t1,t2,x1,x2,y1,y2,sx,ty,dydx,dtds,det,b1,b2,x,y
305C=======================================================================
306 !NFUNC = IPM(10,IMAT)+IPM(6,IMAT)in updmat.f
307 !IFUNC => IPM(10+1:10+NFUNC,IMAT) in updmat.f
308! Calculation intersection
309 xinc = zero
310 yinc = zero
311 np1 = (npc(func+1)-npc(func)) / 2
312 np2 = (npc(fund+1)-npc(fund)) / 2
313 DO j=2,np1+1
314 j1=2*(j-2)
315 s1=pld(npc(func)+j1)
316 t1=pld(npc(func)+j1+1)*fac1
317 DO k=2,np2+1
318 k1=2*(k-2)
319 x1=pld(npc(fund)+k1)
320 y1=pld(npc(fund)+k1+1)*fac2
321 IF(x1 == s1 .AND. y1 == t1 .AND.x1 < zero)THEN
322 xinc = x1
323 yinc = y1
324 GOTO 350
325 ENDIF
326 ENDDO
327 ENDDO
328 DO j=2,np1
329 j1=2*(j-2)
330 s2=pld(npc(func)+j1)
331 s1=pld(npc(func)+j1+2)
332 t2=pld(npc(func)+j1+1)*fac1
333 t1=pld(npc(func)+j1+3)*fac1
334 IF(s1 < zero .OR. s2 < zero) THEN
335 DO k=2,np2
336 k1=2*(k-2)
337 x2=pld(npc(fund)+k1)
338 x1=pld(npc(fund)+k1+2)
339 y2=pld(npc(fund)+k1+1)*fac2
340 y1=pld(npc(fund)+k1+3)*fac2
341 IF(x1 < zero .OR. x2 < zero) THEN
342 dydx = (y2-y1) / (x2-x1)
343 dtds = (t2-t1) / (s2-s1)
344 det = dtds - dydx
345 IF(det /= zero ) THEN
346 b1 = y1 - dydx*x1
347 b2 = t1 - dtds*s1
348 x = (b1 - b2) / det
349 y = (-dydx*b2 + b1*dtds)/det
350 IF(x <= x1 .AND. x >= x2 .AND. x <= s1 .AND. x >= s2 .AND.
351 . y <= y1 .AND. y >= y2 .AND. y <= t1 .AND. y >= t2 ) THEN
352 xinc = x
353 yinc = y
354 GOTO 350
355 ENDIF
356 ENDIF
357 ENDIF
358 ENDDO ! K
359 ENDIF ! S1, S2
360 ENDDO
361 350 CONTINUE
362c-----------
363 RETURN
364 END
365!||====================================================================
366!|| table_inters ../starter/source/tools/curve/func_inters.F
367!||--- called by ------------------------------------------------------
368!|| law119_upd ../starter/source/materials/mat/mat119/law119_upd.F
369!||--- uses -----------------------------------------------------
370!|| table_mod ../starter/share/modules1/table_mod.F
371!||====================================================================
372 SUBROUTINE table_inters(TABLE,FUNC1,FUNC2,FAC1,FAC2,XINT,YINT)
373C-----------------------------------------------
374C M o d u l e s
375C-----------------------------------------------
376 USE table_mod
377C-----------------------------------------------
378C I m p l i c i t T y p e s
379C-----------------------------------------------
380#include "implicit_f.inc"
381C-----------------------------------------------
382C C o m m o n B l o c k s
383C-----------------------------------------------
384#include "com04_c.inc"
385C-----------------------------------------------
386C D u m m y A r g u m e n t s
387C-----------------------------------------------
388 INTEGER :: FUNC1,FUNC2
389 my_real :: xint,yint,fac1,fac2
390 TYPE(ttable), DIMENSION(NTABLE) :: TABLE
391C-----------------------------------------------
392 INTENT(IN) :: func1,func2,fac1,fac2
393 INTENT(INOUT) :: xint,yint
394C-----------------------------------------------
395C L o c a l V a r i a b l e s
396C-----------------------------------------------
397 INTEGER :: J,K,NP1,NP2,NDIM,FOUND
398 my_real :: s1,s2,t1,t2,x1,x2,y1,y2,ax,bx,ay,by,cx,cy,dm,alpha,beta
399C=======================================================================
400c Check common points between 2 curves => intersection
401 ndim = table(func1)%NDIM
402 np1 = SIZE(table(func1)%X(1)%VALUES)
403 np2 = SIZE(table(func2)%X(1)%VALUES)
404 xint = zero
405 yint = zero
406 found = 0
407 DO j=2,np1
408 s1 = table(func1)%X(1)%VALUES(j)
409 t1 = table(func1)%Y%VALUES(j)*fac1
410 DO k=2,np2
411 x1 = table(func2)%X(1)%VALUES(k)
412 y1 = table(func2)%Y%VALUES(k)*fac2
413 IF (s1 > zero .and. x1 == s1 .and. y1 == t1) THEN
414 xint = x1
415 yint = y1
416 found = 1
417 EXIT
418 ENDIF
419 ENDDO
420 IF (found == 1) EXIT
421 ENDDO
422c Check intersection of curve segments
423 IF (found == 0) THEN
424 DO j=2,np1
425 s1 = table(func1)%X(1)%VALUES(j-1)
426 s2 = table(func1)%X(1)%VALUES(j)
427 t1 = table(func1)%Y%VALUES(j-1)*fac1
428 t2 = table(func1)%Y%VALUES(j)*fac1
429 DO k=2,np2
430 x1 = table(func2)%X(1)%VALUES(k-1)
431 x2 = table(func2)%X(1)%VALUES(k)
432 y1 = table(func2)%Y%VALUES(k-1)*fac2
433 y2 = table(func2)%Y%VALUES(k)*fac2
434 IF (x2 < s1 .or. s2 < x1) cycle
435 ax = x2 - x1
436 ay = y2 - y1
437 bx = s1 - s2
438 by = t1 - t2
439 dm = ay*bx - ax*by
440 IF (dm /= zero) THEN ! check if segments are not parallel
441 cx = s1 - x1
442 cy = t1 - y1
443 alpha = (bx * cy - by * cx) / dm
444 beta = (ax * cy - ay * cx) / dm
445 IF (alpha >= zero .and. alpha < one .and.
446 . beta <= zero .and. beta >-one .and.
447 . s1 > zero) THEN
448 xint = x1 + alpha * ax
449 yint = y1 + alpha * ay
450 found = 1
451 EXIT
452 ENDIF
453 ENDIF
454 ENDDO
455 IF (found == 1) EXIT
456 ENDDO
457 END IF
458c-----------
459 RETURN
460 END
my_real
#define my_real
Definition cppsort.cpp:32
alpha
#define alpha
Definition eval.h:35
func_inters
subroutine func_inters(titr, mat_id, func1, func2, fac1, fac2, npc, pld, xint, yint)
Definition func_inters.F:32
table_inters
subroutine table_inters(table, func1, func2, fac1, fac2, xint, yint)
Definition func_inters.F:373
func_inters_shear
subroutine func_inters_shear(titr, mat_id, func, fund, fac1, fac2, npc, pld, xint1, yint1, xint2, yint2)
Definition func_inters.F:134
func_inters_c
subroutine func_inters_c(titr, mat_id, func, fund, fac1, fac2, npc, pld, xinc, yinc)
Definition func_inters.F:279
law88_upd
subroutine law88_upd(iout, titr, uparam, npc, pld, nfunc, ifunc, mat_id, func_id, pm, nfunct)
Definition law88_upd.F:40
message_mod
Definition message_mod.F:1249
names_and_titles_mod
Definition names_and_titles_mod.F:997
names_and_titles_mod::nchartitle
integer, parameter nchartitle
Definition names_and_titles_mod.F:1003
table_mod
Definition table_mod.F:112
starter
program starter
Definition starter.F:39
table_mod::ttable
Definition table_mod.F:125