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