OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
i3fri3.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!|| i3fri3 ../engine/source/interfaces/inter3d/i3fri3.F
25!||--- called by ------------------------------------------------------
26!|| intvo3 ../engine/source/interfaces/inter3d/intvo3.F
27!||--- uses -----------------------------------------------------
28!|| h3d_mod ../engine/share/modules/h3d_mod.F
29!|| output_mod ../common_source/modules/output/output_mod.F90
30!||====================================================================
31 SUBROUTINE i3fri3(OUTPUT,LFT ,LLT ,NFT ,
32 2 X ,E ,IRECT ,MSR ,NSV ,
33 3 IRTL ,NTY ,CST ,IRTLO ,FRIC0 ,
34 4 FRIC ,IMAST ,FSAV ,FSKYI ,ISKY ,
35 5 FCONT ,H3D_DATA,N1 ,N2 ,N3 ,
36 6 IX1 ,IX2 ,IX3 ,IX4 ,H1 ,
37 7 H2 ,H3 ,H4 ,SSC ,TTC ,
38 8 XFACE ,STIF ,XP ,YP ,ZP ,
39 9 FNI )!,FXI ,FYI ,FZI ,FX1 ,
40! 1 FX2 ,FX3 ,FX4 ,FY1 ,FY2 ,
41! 2 FY3 ,FY4 ,FZ1 ,FZ2 ,FZ3 ,
42! 3 FZ4 )
43C-----------------------------------------------
44C M o d u l e s
45C-----------------------------------------------
46 USE h3d_mod
47 USE output_mod, ONLY: output_
48C-----------------------------------------------
49C I m p l i c i t T y p e s
50C-----------------------------------------------
51#include "implicit_f.inc"
52#include "comlock.inc"
53C-----------------------------------------------
54C G l o b a l P a r a m e t e r s
55C-----------------------------------------------
56#include "mvsiz_p.inc"
57C-----------------------------------------------
58C C o m m o n B l o c k s
59C-----------------------------------------------
60#include "scr07_c.inc"
61#include "scr14_c.inc"
62#include "scr16_c.inc"
63#include "com06_c.inc"
64#include "com08_c.inc"
65#include "parit_c.inc"
66#include "scr18_c.inc"
67C-----------------------------------------------
68C D u m m y A r g u m e n t s
69C-----------------------------------------------
70 TYPE(output_), INTENT(inout) :: OUTPUT
71 INTEGER NTY, IMAST,LFT, LLT, NFT
72 my_real
73 . FRIC
74 INTEGER IRECT(4,*), MSR(*), NSV(*), IRTL(*), IRTLO(*), ISKY(*)
75 my_real
76 . X(3,*), E(*), CST(2,*), FRIC0(3,*), FSAV(*),
77 . FSKYI(LSKYI,NFSKYI),FCONT(3,*)
78 TYPE(h3d_database) :: H3D_DATA
79 my_real, DIMENSION(MVSIZ), INTENT(IN) :: n1,n2,n3
80 INTEGER, DIMENSION(MVSIZ), INTENT(IN):: IX1,IX2,IX3,IX4
81 my_real, DIMENSION(MVSIZ), INTENT(IN) :: H1,H2,H3,H4
82 my_real, DIMENSION(MVSIZ), INTENT(IN) :: ssc,ttc,xface,stif
83 my_real, DIMENSION(MVSIZ), INTENT(IN) :: xp,yp,zp,fni
84C-----------------------------------------------
85C L o c a l V a r i a b l e s
86C-----------------------------------------------
87 INTEGER I, IL, LOLD, JJ, NN, J3,
88 . j2, j1, ig, i3, i2, i1
89 INTEGER NISKYL
90 my_real
91 . h(4), xx1(4), xx2(4), xx3(4),ss0, tt0, xc,
92 . yc, zc, xc0, yc0, zc0, sp, sm, tp, tm, ansx, ansy, ansz, fmax,
93 . stf, fti, fn, tn1, tn2, tn3, tn, dtm, econvt
94 my_real, DIMENSION(MVSIZ) :: fxi,fyi,fzi
95 my_real, DIMENSION(MVSIZ) :: fx1,fx2,fx3,fx4
96 my_real, DIMENSION(MVSIZ) :: fy1,fy2,fy3,fy4
97 my_real, DIMENSION(MVSIZ) :: fz1,fz2,fz3,fz4
98C-----------------------------------------------
99C Source Lines
100C-----------------------------------------------
101 DO 300 i=lft,llt
102 fxi(i)=zero
103 fyi(i)=zero
104 fzi(i)=zero
105 il=i+nft
106 econvt = zero
107 IF(xface(i)==zero) THEN
108C-------------------------------
109C POINT NON IMPACTE
110C-------------------------------
111 irtlo(il)=0
112 fric0(1,il)=zero
113 fric0(2,il)=zero
114 fric0(3,il)=zero
115 ELSE
116C
117 lold=iabs(irtlo(il))
118 IF(lold==0)THEN
119C-------------------------------
120C POINT NON IMPACTE PRECEDEMENT
121C-------------------------------
122 irtlo(il)=irtl(il)*xface(i)
123 cst(1,il)=ssc(i)
124 cst(2,il)=ttc(i)
125 ELSE
126C-------------------------------
127C POINT IMPACTE PRECEDEMENT
128C-------------------------------
129 ss0=cst(1,il)
130 tt0=cst(2,il)
131 fxi(i)=fric0(1,il)
132 fyi(i)=fric0(2,il)
133 fzi(i)=fric0(3,il)
134C
135 xc=xp(i)
136 yc=yp(i)
137 zc=zp(i)
138 DO jj=1,4
139 nn=msr(irect(jj,lold))
140 xx1(jj)=x(1,nn)
141 xx2(jj)=x(2,nn)
142 xx3(jj)=x(3,nn)
143 ENDDO
144 xc0=zero
145 yc0=zero
146 zc0=zero
147 sp=one+ss0
148 sm=one-ss0
149 tp= fourth*(one+tt0)
150 tm= fourth*(one-tt0)
151 h(1)=tm*sm
152 h(2)=tm*sp
153 h(3)=tp*sp
154 h(4)=tp*sm
155 DO jj=1,4
156 xc0=xc0+h(jj)*xx1(jj)
157 yc0=yc0+h(jj)*xx2(jj)
158 zc0=zc0+h(jj)*xx3(jj)
159 ENDDO
160 ansx= (xc-xc0)
161 ansy= (yc-yc0)
162 ansz= (zc-zc0)
163 econvt = econvt + half*(fxi(i)*ansx+fyi(i)*ansy+fzi(i)*ansz)
164C
165 fmax= -min(fric*fni(i),zero)
166C
167 stf=em01*stif(i)
168 fxi(i)=fxi(i) + ansx*stf
169 fyi(i)=fyi(i) + ansy*stf
170 fzi(i)=fzi(i) + ansz*stf
171 fti=sqrt(fxi(i)*fxi(i)+fyi(i)*fyi(i)+fzi(i)*fzi(i))
172C
173 fn=fxi(i)*n1(i)+fyi(i)*n2(i)+fzi(i)*n3(i)
174 tn1=fxi(i)-n1(i)*fn
175 tn2=fyi(i)-n2(i)*fn
176 tn3=fzi(i)-n3(i)*fn
177 tn=sqrt(tn1*tn1+tn2*tn2+tn3*tn3)
178 IF(tn/=zero)THEN
179 tn1=tn1/tn
180 tn2=tn2/tn
181 tn3=tn3/tn
182 ELSE
183 tn3=zero
184 tn=sqrt(n1(i)*n1(i)+n2(i)*n2(i))
185 IF(tn/=zero)THEN
186 tn2=-n1(i)/tn
187 tn1=n2(i)/tn
188 ELSE
189 tn2=zero
190 tn1=one
191 ENDIF
192 ENDIF
193C
194 IF(fti>fmax)THEN
195C-------------------------------
196C POINT GLISSANT
197C-------------------------------
198 fxi(i)=tn1*fmax
199 fyi(i)=tn2*fmax
200 fzi(i)=tn3*fmax
201 fric0(1,il)=fxi(i)
202 fric0(2,il)=fyi(i)
203 fric0(3,il)=fzi(i)
204 irtlo(il)=irtl(il)*xface(i)
205 cst(1,il)=ssc(i)
206 cst(2,il)=ttc(i)
207 ELSE
208C-------------------------------
209C POINT NON GLISSANT
210C-------------------------------
211 fxi(i)=tn1*fti
212 fyi(i)=tn2*fti
213 fzi(i)=tn3*fti
214 ENDIF
215 econvt = econvt + half*(fxi(i)*ansx+fyi(i)*ansy+fzi(i)*ansz)
216 ENDIF
217 ENDIF
218C
219 300 CONTINUE
220C
221C
222!$OMP ATOMIC
223 econtv = econtv + econvt ! Frictional contact energy
224!$OMP ATOMIC
225 fsav(27) = fsav(27) + econvt
226C---------------------------------
227C SAUVEGARDE DE L'IMPULSION TOTALE
228C---------------------------------
229 dtm=imast*dt12
230 DO 350 i=lft,llt
231 fsav(4)=fsav(4)+fxi(i)*dtm
232 fsav(5)=fsav(5)+fyi(i)*dtm
233 fsav(6)=fsav(6)+fzi(i)*dtm
234 350 CONTINUE
235C
236C
237 DO 400 i=lft,llt
238 fx1(i)=fxi(i)*h1(i)
239 fy1(i)=fyi(i)*h1(i)
240 fz1(i)=fzi(i)*h1(i)
241C
242 fx2(i)=fxi(i)*h2(i)
243 fy2(i)=fyi(i)*h2(i)
244 fz2(i)=fzi(i)*h2(i)
245C
246 fx3(i)=fxi(i)*h3(i)
247 fy3(i)=fyi(i)*h3(i)
248 fz3(i)=fzi(i)*h3(i)
249C
250 fx4(i)=fxi(i)*h4(i)
251 fy4(i)=fyi(i)*h4(i)
252 fz4(i)=fzi(i)*h4(i)
253 400 CONTINUE
254C
255 IF(iparit==0)THEN
256C---------------------------------
257C FRICTION MAIN
258C---------------------------------
259 DO i=lft,llt
260 il=i+nft
261 j3=3*ix1(i)
262 j2=j3-1
263 j1=j2-1
264 e(j1)=e(j1)+fx1(i)
265 e(j2)=e(j2)+fy1(i)
266 e(j3)=e(j3)+fz1(i)
267C
268 j3=3*ix2(i)
269 j2=j3-1
270 j1=j2-1
271 e(j1)=e(j1)+fx2(i)
272 e(j2)=e(j2)+fy2(i)
273 e(j3)=e(j3)+fz2(i)
274C
275 j3=3*ix3(i)
276 j2=j3-1
277 j1=j2-1
278 e(j1)=e(j1)+fx3(i)
279 e(j2)=e(j2)+fy3(i)
280 e(j3)=e(j3)+fz3(i)
281C
282 j3=3*ix4(i)
283 j2=j3-1
284 j1=j2-1
285 e(j1)=e(j1)+fx4(i)
286 e(j2)=e(j2)+fy4(i)
287 e(j3)=e(j3)+fz4(i)
288C---------------------------------
289C FRICTION SECOND C---------------------------------
290 ig=nsv(il)
291 i3=3*ig
292 i2=i3-1
293 i1=i2-1
294 e(i1)=e(i1)-fxi(i)
295 e(i2)=e(i2)-fyi(i)
296 e(i3)=e(i3)-fzi(i)
297C
298 ENDDO
299C
300 ELSE
301C
302#include "lockon.inc"
303 niskyl = nisky
304 nisky = nisky + 5 * llt
305#include "lockoff.inc"
306C
307 IF(kdtint==0)THEN
308 DO i=lft,llt
309 niskyl = niskyl + 1
310 fskyi(niskyl,1)=fx1(i)
311 fskyi(niskyl,2)=fy1(i)
312 fskyi(niskyl,3)=fz1(i)
313 fskyi(niskyl,4)=zero
314 isky(niskyl) = ix1(i)
315 niskyl = niskyl + 1
316 fskyi(niskyl,1)=fx2(i)
317 fskyi(niskyl,2)=fy2(i)
318 fskyi(niskyl,3)=fz2(i)
319 fskyi(niskyl,4)=zero
320 isky(niskyl) = ix2(i)
321 niskyl = niskyl + 1
322 fskyi(niskyl,1)=fx3(i)
323 fskyi(niskyl,2)=fy3(i)
324 fskyi(niskyl,3)=fz3(i)
325 fskyi(niskyl,4)=zero
326 isky(niskyl) = ix3(i)
327 niskyl = niskyl + 1
328 fskyi(niskyl,1)=fx4(i)
329 fskyi(niskyl,2)=fy4(i)
330 fskyi(niskyl,3)=fz4(i)
331 fskyi(niskyl,4)=zero
332 isky(niskyl) = ix4(i)
333 niskyl = niskyl + 1
334 fskyi(niskyl,1)=-fxi(i)
335 fskyi(niskyl,2)=-fyi(i)
336 fskyi(niskyl,3)=-fzi(i)
337 fskyi(niskyl,4)=zero
338 il=i+nft
339 isky(niskyl) = nsv(il)
340 ENDDO
341 ELSE
342 DO i=lft,llt
343 niskyl = niskyl + 1
344 fskyi(niskyl,1)=fx1(i)
345 fskyi(niskyl,2)=fy1(i)
346 fskyi(niskyl,3)=fz1(i)
347 fskyi(niskyl,4)=zero
348 fskyi(niskyl,5)=zero
349 isky(niskyl) = ix1(i)
350 niskyl = niskyl + 1
351 fskyi(niskyl,1)=fx2(i)
352 fskyi(niskyl,2)=fy2(i)
353 fskyi(niskyl,3)=fz2(i)
354 fskyi(niskyl,4)=zero
355 fskyi(niskyl,5)=zero
356 isky(niskyl) = ix2(i)
357 niskyl = niskyl + 1
358 fskyi(niskyl,1)=fx3(i)
359 fskyi(niskyl,2)=fy3(i)
360 fskyi(niskyl,3)=fz3(i)
361 fskyi(niskyl,4)=zero
362 fskyi(niskyl,5)=zero
363 isky(niskyl) = ix3(i)
364 niskyl = niskyl + 1
365 fskyi(niskyl,1)=fx4(i)
366 fskyi(niskyl,2)=fy4(i)
367 fskyi(niskyl,3)=fz4(i)
368 fskyi(niskyl,4)=zero
369 fskyi(niskyl,5)=zero
370 isky(niskyl) = ix4(i)
371 niskyl = niskyl + 1
372 fskyi(niskyl,1)=-fxi(i)
373 fskyi(niskyl,2)=-fyi(i)
374 fskyi(niskyl,3)=-fzi(i)
375 fskyi(niskyl,4)=zero
376 fskyi(niskyl,5)=zero
377 il=i+nft
378 isky(niskyl) = nsv(il)
379 ENDDO
380 ENDIF
381 ENDIF
382C
383 IF(anim_v(4)+outp_v(4)+h3d_data%N_VECT_CONT>0.AND.
384 . ((tt>=output%TANIM .AND. tt<=output%TANIM_STOP).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP).OR.
385 . (manim>=4.AND.manim<=15).OR. h3d_data%MH3D /= 0))THEN
386#include "lockon.inc"
387 DO i=1,llt
388 fcont(1,ix1(i)) =fcont(1,ix1(i)) + fx1(i)
389 fcont(2,ix1(i)) =fcont(2,ix1(i)) + fy1(i)
390 fcont(3,ix1(i)) =fcont(3,ix1(i)) + fz1(i)
391 fcont(1,ix2(i)) =fcont(1,ix2(i)) + fx2(i)
392 fcont(2,ix2(i)) =fcont(2,ix2(i)) + fy2(i)
393 fcont(3,ix2(i)) =fcont(3,ix2(i)) + fz2(i)
394 fcont(1,ix3(i)) =fcont(1,ix3(i)) + fx3(i)
395 fcont(2,ix3(i)) =fcont(2,ix3(i)) + fy3(i)
396 fcont(3,ix3(i)) =fcont(3,ix3(i)) + fz3(i)
397 fcont(1,ix4(i)) =fcont(1,ix4(i)) + fx4(i)
398 fcont(2,ix4(i)) =fcont(2,ix4(i)) + fy4(i)
399 fcont(3,ix4(i)) =fcont(3,ix4(i)) + fz4(i)
400 fcont(1,nsv(i+nft))=fcont(1,nsv(i+nft))- fxi(i)
401 fcont(2,nsv(i+nft))=fcont(2,nsv(i+nft))- fyi(i)
402 fcont(3,nsv(i+nft))=fcont(3,nsv(i+nft))- fzi(i)
403 ENDDO
404#include "lockoff.inc"
405 ENDIF
406C
407 RETURN
408 END
409C 054 +++
410!||====================================================================
411!|| i5fri3 ../engine/source/interfaces/inter3d/i3fri3.F
412!||--- called by ------------------------------------------------------
413!|| intvo3 ../engine/source/interfaces/inter3d/intvo3.F
414!||--- uses -----------------------------------------------------
415!|| h3d_mod ../engine/share/modules/h3d_mod.F
416!|| output_mod ../common_source/modules/output/output_mod.F90
417!||====================================================================
418 SUBROUTINE i5fri3(OUTPUT,LFT ,LLT ,NFT ,IPARI ,
419 2 X ,E ,IRECT ,MSR ,NSV ,
420 3 IRTL ,NTY ,CST ,IRTLO ,FRIC0 ,
421 4 FRIC ,IMAST ,FSAV ,FSKYI ,ISKY ,
422 5 FCONT ,V ,CF ,FROT_P,FREQ ,
423 6 FTSAV ,FTCONT,H3D_DATA,N1 ,N2 ,
424 7 N3 ,IX1 ,IX2 ,IX3 ,IX4 ,
425 8 XP ,YP ,ZP ,SSC ,TTC ,
426 9 XFACE ,STIF ,H1 ,H2 ,H3 ,
427 1 H4 ,AREA ,FNI)
428C-----------------------------------------------
429C M o d u l e s
430C-----------------------------------------------
431 USE h3d_mod
432 USE output_mod, ONLY: output_
433C-----------------------------------------------
434C I m p l i c i t T y p e s
435C-----------------------------------------------
436#include "implicit_f.inc"
437#include "comlock.inc"
438C-----------------------------------------------
439C G l o b a l P a r a m e t e r s
440C-----------------------------------------------
441#include "mvsiz_p.inc"
442C-----------------------------------------------
443C C o m m o n B l o c k s
444C-----------------------------------------------
445#include "scr03_c.inc"
446#include "scr07_c.inc"
447#include "scr14_c.inc"
448#include "scr16_c.inc"
449#include "com06_c.inc"
450#include "com08_c.inc"
451#include "parit_c.inc"
452#include "scr18_c.inc"
453#include "impl1_c.inc"
454C-----------------------------------------------
455C D u m m y A r g u m e n t s
456C-----------------------------------------------
457 TYPE(output_), INTENT(inout) :: OUTPUT
458 INTEGER IPARI(*), NTY, IMAST, LFT, LLT, NFT
459 my_real FRIC
460 INTEGER IRECT(4,*), MSR(*), NSV(*), IRTL(*), IRTLO(*), ISKY(*)
461 my_real X(3,*), E(*), CST(2,*), FRIC0(3,*), FSAV(*),
462 . FSKYI(LSKYI,NFSKYI),
463 . V(3,*), FCONT(3,*),CF(*), FREQ, FROT_P(*),FTSAV(*),
464 . FTCONT(3,*)
465 TYPE(H3D_DATABASE) :: H3D_DATA
466 INTEGER, DIMENSION(MVSIZ), INTENT(IN) :: IX1,IX2,IX3,IX4
467 my_real, DIMENSION(MVSIZ), INTENT(IN) :: N1,N2,N3
468 my_real, DIMENSION(MVSIZ), INTENT(IN) :: XP,YP,ZP
469 my_real, DIMENSION(MVSIZ), INTENT(IN) :: SSC,TTC,XFACE,STIF
470 my_real, DIMENSION(MVSIZ), INTENT(IN) :: H1,H2,H3,H4,AREA,FNI
471C-----------------------------------------------
472C L o c a l V a r i a b l e s
473C-----------------------------------------------
474 INTEGER I, IL, LOLD, JJ, NN, J3,J2, J1, IG, I3, I2, I1, K, IFQ, MFROT, NISKYL
475 my_real
476 . h(4), xx1(4), xx2(4), xx3(4),
477 . fxi(mvsiz), fyi(mvsiz),
478 . fzi(mvsiz), fx1(mvsiz), fx2(mvsiz), fx3(mvsiz), fx4(mvsiz), fy1(mvsiz), fy2(mvsiz),
479 . fy3(mvsiz), fy4(mvsiz), fz1(mvsiz), fz2(mvsiz), fz3(mvsiz), fz4(mvsiz),
480 . ss0, tt0, xc, econvt, alpha, alphi,
481 . yc, zc, xc0, yc0, zc0, sp, sm, tp, tm, ansx, ansy, ansz, fmax,
482 . stf, fti, fn, tn1, tn2, tn3, tn, dtm, xmu, vx,vy,vz,vv,v2,p,
483 . vv1,vv2,dmu,aa
484 INTEGER IRTLO_I(MVSIZ)
485 my_real fric0_i(3,mvsiz)
486C-----------------------------------------------
487C S o u r c e L i n e s
488C-----------------------------------------------
489C------for implicit run---save previous values-
490 IF (inconv/=1) THEN
491 DO i=lft,llt
492 il=i+nft
493 irtlo_i(i)=irtlo(il)
494 fric0_i(1,i)=fric0(1,il)
495 fric0_i(2,i)=fric0(2,il)
496 fric0_i(3,i)=fric0(3,il)
497 ENDDO
498 ENDIF
499 DO 300 i=lft,llt
500 fxi(i)=zero
501 fyi(i)=zero
502 fzi(i)=zero
503 il=i+nft
504 econvt = zero
505 IF(xface(i)==zero) THEN
506C-------------------------------
507C Point without Impact
508C-------------------------------
509 irtlo(il)=0
510 fric0(1,il)=zero
511 fric0(2,il)=zero
512 fric0(3,il)=zero
513 ELSE
514C
515 lold=iabs(irtlo(il))
516 IF(lold==0)THEN
517C-------------------------------
518C Point not impacted previously
519C-------------------------------
520 irtlo(il)=irtl(il)*xface(i)
521 cst(1,il)=ssc(i)
522 cst(2,il)=ttc(i)
523 ELSE
524C-------------------------------
525C Point previously impacted
526C-------------------------------
527 ss0=cst(1,il)
528 tt0=cst(2,il)
529 fxi(i)=fric0(1,il)
530 fyi(i)=fric0(2,il)
531 fzi(i)=fric0(3,il)
532C
533 xc=xp(i)
534 yc=yp(i)
535 zc=zp(i)
536 DO jj=1,4
537 nn=msr(irect(jj,lold))
538 xx1(jj)=x(1,nn)
539 xx2(jj)=x(2,nn)
540 xx3(jj)=x(3,nn)
541 ENDDO
542 xc0=zero
543 yc0=zero
544 zc0=zero
545 sp=one+ss0
546 sm=one-ss0
547 tp=fourth*(one + tt0)
548 tm=fourth*(one - tt0)
549 h(1)=tm*sm
550 h(2)=tm*sp
551 h(3)=tp*sp
552 h(4)=tp*sm
553 DO jj=1,4
554 xc0=xc0+h(jj)*xx1(jj)
555 yc0=yc0+h(jj)*xx2(jj)
556 zc0=zc0+h(jj)*xx3(jj)
557 ENDDO
558 ansx= (xc-xc0)
559 ansy= (yc-yc0)
560 ansz= (zc-zc0)
561 econvt = econvt + half*(fxi(i)*ansx+fyi(i)*ansy+fzi(i)*ansz)
562 stf=em01*stif(i)
563 fxi(i)=fxi(i) + ansx*stf
564 fyi(i)=fyi(i) + ansy*stf
565 fzi(i)=fzi(i) + ansz*stf
566Cfriction filtering
567 IF (codvers>=44) THEN
568 ifq = ipari(31)
569 IF (ifq>0) THEN
570 IF (ifq==3) freq = max(one,freq*dt12)
571 alpha = freq
572 alphi = one - alpha
573 k = 3*(il-1)
574 IF (fni(i)==0) THEN
575 ftsav(k+1) = zero
576 ftsav(k+2) = zero
577 ftsav(k+3) = zero
578 ELSE
579 fxi(i)= alpha*fxi(i) + alphi*ftsav(k+1)
580 fyi(i)= alpha*fyi(i) + alphi*ftsav(k+2)
581 fzi(i)= alpha*fzi(i) + alphi*ftsav(k+3)
582 IF (inconv==1) THEN
583 ftsav(k+1) = fxi(i)
584 ftsav(k+2) = fyi(i)
585 ftsav(k+3) = fzi(i)
586 ENDIF
587 ENDIF
588 ENDIF
589 ENDIF
590 fti=sqrt(fxi(i)*fxi(i)+fyi(i)*fyi(i)+fzi(i)*fzi(i))
591C
592 fn =fxi(i)*n1(i)+fyi(i)*n2(i)+fzi(i)*n3(i)
593 tn1=fxi(i)-n1(i)*fn
594 tn2=fyi(i)-n2(i)*fn
595 tn3=fzi(i)-n3(i)*fn
596 tn =sqrt(tn1*tn1+tn2*tn2+tn3*tn3)
597 IF(tn/=zero)THEN
598 tn1=tn1/tn
599 tn2=tn2/tn
600 tn3=tn3/tn
601 ELSE
602 tn3=zero
603 tn=sqrt(n1(i)*n1(i)+n2(i)*n2(i))
604 IF(tn/=zero)THEN
605 tn2=-n1(i)/tn
606 tn1= n2(i)/tn
607 ELSE
608 tn2=zero
609 tn1=one
610 ENDIF
611 ENDIF
612C
613 p = -fni(i)/area(i)
614 xmu = fric
615 IF(codvers>=44) THEN
616C---------------------------------
617C NEW FRICTION MODELS
618C---------------------------------
619 mfrot = ipari(30)
620 IF(mfrot>=1) THEN
621 ig=nsv(il)
622 vx = v(1,ig) - h(1)*v(1,ix1(i)) - h(2)*v(1,ix2(i))
623 . - h(3)*v(1,ix3(i)) - h(4)*v(1,ix4(i))
624 vy = v(2,ig) - h(1)*v(2,ix1(i)) - h(2)*v(2,ix2(i))
625 . - h(3)*v(2,ix3(i)) - h(4)*v(2,ix4(i))
626 vz = v(3,ig) - h(1)*v(3,ix1(i)) - h(2)*v(3,ix2(i))
627 . - h(3)*v(3,ix3(i)) - h(4)*v(3,ix4(i))
628C--- tangantial velocities
629 aa = n1(i)*vx + n2(i)*vy + n3(i)*vz
630 vx = vx - n1(i)*aa
631 vy = vy - n2(i)*aa
632 vz = vz - n3(i)*aa
633 v2 = vx*vx + vy*vy + vz*vz
634 vv = sqrt(max(em30, v2))
635 ENDIF
636 IF(mfrot==1)THEN
637 xmu = fric + (frot_p(1) + frot_p(4)*p ) * p
638 . +(frot_p(2) + frot_p(3)*p) * vv + frot_p(5)*v2
639 ELSEIF(mfrot==2)THEN
640C--- Darmstad's law
641 xmu = fric
642 . + frot_p(1)*exp(frot_p(2)*vv)*p**2
643 . + frot_p(3)*exp(frot_p(4)*vv)*p
644 . + frot_p(5)*exp(frot_p(6)*vv)
645 ELSEIF(mfrot==3)THEN
646C--- Renard's law
647 IF(vv>=0.AND.vv<=frot_p(5)) THEN
648 dmu = frot_p(3)-frot_p(1)
649 vv1 = vv / frot_p(5)
650 xmu = frot_p(1)+ dmu*vv1*(two-vv1)
651 ELSEIF(vv>frot_p(5).AND.vv<frot_p(6)) THEN
652 dmu = frot_p(4)-frot_p(3)
653 vv1 = (vv - frot_p(5))/(frot_p(6)-frot_p(5))
654 xmu = frot_p(3)+ dmu * (three-two*vv1)*vv1**2
655 ELSE
656 dmu = frot_p(2)-frot_p(4)
657 vv2 = (vv - frot_p(6))**2
658 xmu = frot_p(2) - dmu / (one + dmu*vv2)
659 ENDIF
660 ENDIF
661 fmax = -min(xmu*fni(i),zero)
662 IF(fti>fmax)THEN
663C--- sliding
664 fxi(i)=tn1*fmax
665 fyi(i)=tn2*fmax
666 fzi(i)=tn3*fmax
667 fric0(1,il)=fxi(i)
668 fric0(2,il)=fyi(i)
669 fric0(3,il)=fzi(i)
670 irtlo(il)=irtl(il)*xface(i)
671 cst(1,il)=ssc(i)
672 cst(2,il)=ttc(i)
673 ELSE
674C--- tied
675 fxi(i)=tn1*fti
676 fyi(i)=tn2*fti
677 fzi(i)=tn3*fti
678 ENDIF
679 ELSE
680C-----------------------------------------------------
681C CODVERS < 44
682C-----------------------------------------------------
683 ig=nsv(il)
684 vx = v(1,ig) - h(1)*v(1,ix1(i)) - h(2)*v(1,ix2(i))
685 . - h(3)*v(1,ix3(i)) - h(4)*v(1,ix4(i))
686 vy = v(2,ig) - h(1)*v(2,ix1(i)) - h(2)*v(2,ix2(i))
687 . - h(3)*v(2,ix3(i)) - h(4)*v(2,ix4(i))
688 vz = v(3,ig) - h(1)*v(3,ix1(i)) - h(2)*v(3,ix2(i))
689 . - h(3)*v(3,ix3(i)) - h(4)*v(3,ix4(i))
690 vv = vx*vx + vy*vy + vz*vz
691 xmu = fric + ( cf(1) + cf(4)*p ) * p
692 . + ( cf(2) + cf(3)*p ) * sqrt(vv) + cf(5)*vv
693 fmax= -min(xmu*fni(i),zero)
694 IF(fti>fmax)THEN
695C------------------------------------------------
696C SLIDING POINT
697C-------------------------------
698 fxi(i)=tn1*fmax
699 fyi(i)=tn2*fmax
700 fzi(i)=tn3*fmax
701 fric0(1,il)=fxi(i)
702 fric0(2,il)=fyi(i)
703 fric0(3,il)=fzi(i)
704 irtlo(il)=irtl(il)*xface(i)
705 cst(1,il)=ssc(i)
706 cst(2,il)=ttc(i)
707 ELSE
708C-------------------------------
709C NON SLIDING POINTS
710C-------------------------------
711 fxi(i)=tn1*fti
712 fyi(i)=tn2*fti
713 fzi(i)=tn3*fti
714 ENDIF
715 ENDIF
716 econvt = econvt + half*(fxi(i)*ansx+fyi(i)*ansy+fzi(i)*ansz)
717 ENDIF
718 ENDIF
719C
720 300 CONTINUE
721C
722 IF (inconv==1) THEN
723!$OMP ATOMIC
724 econtv = econtv + econvt ! Frictional contact energy
725!$OMP ATOMIC
726 fsav(27) = fsav(27) + econvt
727C---------------------------------
728C TOTAL IMPULSE BACKUP
729C---------------------------------
730 dtm=imast*dt12
731 DO i=lft,llt
732 fsav(4)=fsav(4)+fxi(i)*dtm
733 fsav(5)=fsav(5)+fyi(i)*dtm
734 fsav(6)=fsav(6)+fzi(i)*dtm
735 ENDDO
736 END IF !(INCONV==1) THEN
737C
738 DO i=lft,llt
739 fx1(i)=fxi(i)*h1(i)
740 fy1(i)=fyi(i)*h1(i)
741 fz1(i)=fzi(i)*h1(i)
742C
743 fx2(i)=fxi(i)*h2(i)
744 fy2(i)=fyi(i)*h2(i)
745 fz2(i)=fzi(i)*h2(i)
746C
747 fx3(i)=fxi(i)*h3(i)
748 fy3(i)=fyi(i)*h3(i)
749 fz3(i)=fzi(i)*h3(i)
750C
751 fx4(i)=fxi(i)*h4(i)
752 fy4(i)=fyi(i)*h4(i)
753 fz4(i)=fzi(i)*h4(i)
754 ENDDO
755C
756 IF(iparit==0)THEN
757C---------------------------------
758C FRICTION MAIN
759C---------------------------------
760 DO i=lft,llt
761 il=i+nft
762 j3=3*ix1(i)
763 j2=j3-1
764 j1=j2-1
765 e(j1)=e(j1)+fx1(i)
766 e(j2)=e(j2)+fy1(i)
767 e(j3)=e(j3)+fz1(i)
768C
769 j3=3*ix2(i)
770 j2=j3-1
771 j1=j2-1
772 e(j1)=e(j1)+fx2(i)
773 e(j2)=e(j2)+fy2(i)
774 e(j3)=e(j3)+fz2(i)
775C
776 j3=3*ix3(i)
777 j2=j3-1
778 j1=j2-1
779 e(j1)=e(j1)+fx3(i)
780 e(j2)=e(j2)+fy3(i)
781 e(j3)=e(j3)+fz3(i)
782C
783 j3=3*ix4(i)
784 j2=j3-1
785 j1=j2-1
786 e(j1)=e(j1)+fx4(i)
787 e(j2)=e(j2)+fy4(i)
788 e(j3)=e(j3)+fz4(i)
789C---------------------------------
790C FRICTION SECOND C---------------------------------
791 ig=nsv(il)
792 i3=3*ig
793 i2=i3-1
794 i1=i2-1
795 e(i1)=e(i1)-fxi(i)
796 e(i2)=e(i2)-fyi(i)
797 e(i3)=e(i3)-fzi(i)
798C
799 ENDDO
800C
801 ELSE
802C
803#include "lockon.inc"
804 niskyl = nisky
805 nisky = nisky + 5 * llt
806#include "lockoff.inc"
807C
808 IF(kdtint==0)THEN
809 DO 440 i=lft,llt
810 niskyl = niskyl + 1
811 fskyi(niskyl,1)=fx1(i)
812 fskyi(niskyl,2)=fy1(i)
813 fskyi(niskyl,3)=fz1(i)
814 fskyi(niskyl,4)=zero
815 isky(niskyl) = ix1(i)
816 niskyl = niskyl + 1
817 fskyi(niskyl,1)=fx2(i)
818 fskyi(niskyl,2)=fy2(i)
819 fskyi(niskyl,3)=fz2(i)
820 fskyi(niskyl,4)=zero
821 isky(niskyl) = ix2(i)
822 niskyl = niskyl + 1
823 fskyi(niskyl,1)=fx3(i)
824 fskyi(niskyl,2)=fy3(i)
825 fskyi(niskyl,3)=fz3(i)
826 fskyi(niskyl,4)=zero
827 isky(niskyl) = ix3(i)
828 niskyl = niskyl + 1
829 fskyi(niskyl,1)=fx4(i)
830 fskyi(niskyl,2)=fy4(i)
831 fskyi(niskyl,3)=fz4(i)
832 fskyi(niskyl,4)=zero
833 isky(niskyl) = ix4(i)
834 niskyl = niskyl + 1
835 fskyi(niskyl,1)=-fxi(i)
836 fskyi(niskyl,2)=-fyi(i)
837 fskyi(niskyl,3)=-fzi(i)
838 fskyi(niskyl,4)=zero
839 il=i+nft
840 isky(niskyl) = nsv(il)
841 440 CONTINUE
842 ELSE
843 DO i=lft,llt
844 niskyl = niskyl + 1
845 fskyi(niskyl,1)=fx1(i)
846 fskyi(niskyl,2)=fy1(i)
847 fskyi(niskyl,3)=fz1(i)
848 fskyi(niskyl,4)=zero
849 fskyi(niskyl,5)=zero
850 isky(niskyl) = ix1(i)
851 niskyl = niskyl + 1
852 fskyi(niskyl,1)=fx2(i)
853 fskyi(niskyl,2)=fy2(i)
854 fskyi(niskyl,3)=fz2(i)
855 fskyi(niskyl,4)=zero
856 fskyi(niskyl,5)=zero
857 isky(niskyl) = ix2(i)
858 niskyl = niskyl + 1
859 fskyi(niskyl,1)=fx3(i)
860 fskyi(niskyl,2)=fy3(i)
861 fskyi(niskyl,3)=fz3(i)
862 fskyi(niskyl,4)=zero
863 fskyi(niskyl,5)=zero
864 isky(niskyl) = ix3(i)
865 niskyl = niskyl + 1
866 fskyi(niskyl,1)=fx4(i)
867 fskyi(niskyl,2)=fy4(i)
868 fskyi(niskyl,3)=fz4(i)
869 fskyi(niskyl,4)=zero
870 fskyi(niskyl,5)=zero
871 isky(niskyl) = ix4(i)
872 niskyl = niskyl + 1
873 fskyi(niskyl,1)=-fxi(i)
874 fskyi(niskyl,2)=-fyi(i)
875 fskyi(niskyl,3)=-fzi(i)
876 fskyi(niskyl,4)=zero
877 fskyi(niskyl,5)=zero
878 il=i+nft
879 isky(niskyl) = nsv(il)
880 ENDDO
881 ENDIF
882 ENDIF
883C------for implicit run---restore previous values-
884 IF (inconv/=1) THEN
885 DO i=lft,llt
886 il=i+nft
887 irtlo(il)=irtlo_i(i)
888 fric0(1,il)=fric0_i(1,i)
889 fric0(2,il)=fric0_i(2,i)
890 fric0(3,il)=fric0_i(3,i)
891 ENDDO
892 RETURN
893 ENDIF
894C
895 IF(anim_v(4)+outp_v(4)+h3d_data%N_VECT_CONT>0 .AND.
896 . ((tt>=output%TANIM .AND. tt<=output%TANIM_STOP).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP).OR.
897 . (manim>=4.AND.manim<=15).OR. h3d_data%MH3D /= 0))THEN
898#include "lockon.inc"
899 DO i=1,llt
900 fcont(1,ix1(i)) =fcont(1,ix1(i)) + fx1(i)
901 fcont(2,ix1(i)) =fcont(2,ix1(i)) + fy1(i)
902 fcont(3,ix1(i)) =fcont(3,ix1(i)) + fz1(i)
903 fcont(1,ix2(i)) =fcont(1,ix2(i)) + fx2(i)
904 fcont(2,ix2(i)) =fcont(2,ix2(i)) + fy2(i)
905 fcont(3,ix2(i)) =fcont(3,ix2(i)) + fz2(i)
906 fcont(1,ix3(i)) =fcont(1,ix3(i)) + fx3(i)
907 fcont(2,ix3(i)) =fcont(2,ix3(i)) + fy3(i)
908 fcont(3,ix3(i)) =fcont(3,ix3(i)) + fz3(i)
909 fcont(1,ix4(i)) =fcont(1,ix4(i)) + fx4(i)
910 fcont(2,ix4(i)) =fcont(2,ix4(i)) + fy4(i)
911 fcont(3,ix4(i)) =fcont(3,ix4(i)) + fz4(i)
912 fcont(1,nsv(i+nft))=fcont(1,nsv(i+nft))- fxi(i)
913 fcont(2,nsv(i+nft))=fcont(2,nsv(i+nft))- fyi(i)
914 fcont(3,nsv(i+nft))=fcont(3,nsv(i+nft))- fzi(i)
915 ENDDO
916#include "lockoff.inc"
917 ENDIF
918C
919 IF(anim_v(12)+outp_v(12)+h3d_data%N_VECT_PCONT>0.AND.
920 . ((tt>=output%TANIM .AND. tt<=output%TANIM_STOP).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP) .OR.
921 . (manim>=4.AND.manim<=15).OR.h3d_data%MH3D/=0))THEN
922#include "lockon.inc"
923 DO i=1,llt
924 ftcont(1,ix1(i)) =ftcont(1,ix1(i)) + fx1(i)
925 ftcont(2,ix1(i)) =ftcont(2,ix1(i)) + fy1(i)
926 ftcont(3,ix1(i)) =ftcont(3,ix1(i)) + fz1(i)
927 ftcont(1,ix2(i)) =ftcont(1,ix2(i)) + fx2(i)
928 ftcont(2,ix2(i)) =ftcont(2,ix2(i)) + fy2(i)
929 ftcont(3,ix2(i)) =ftcont(3,ix2(i)) + fz2(i)
930 ftcont(1,ix3(i)) =ftcont(1,ix3(i)) + fx3(i)
931 ftcont(2,ix3(i)) =ftcont(2,ix3(i)) + fy3(i)
932 ftcont(3,ix3(i)) =ftcont(3,ix3(i)) + fz3(i)
933 ftcont(1,ix4(i)) =ftcont(1,ix4(i)) + fx4(i)
934 ftcont(2,ix4(i)) =ftcont(2,ix4(i)) + fy4(i)
935 ftcont(3,ix4(i)) =ftcont(3,ix4(i)) + fz4(i)
936 ftcont(1,nsv(i+nft))=ftcont(1,nsv(i+nft))- fxi(i)
937 ftcont(2,nsv(i+nft))=ftcont(2,nsv(i+nft))- fyi(i)
938 ftcont(3,nsv(i+nft))=ftcont(3,nsv(i+nft))- fzi(i)
939 ENDDO
940#include "lockoff.inc"
941 ENDIF
942C
943 RETURN
944 END
alpha
#define alpha
Definition eval.h:35
i3fri3
subroutine i3fri3(output, lft, llt, nft, x, e, irect, msr, nsv, irtl, nty, cst, irtlo, fric0, fric, imast, fsav, fskyi, isky, fcont, h3d_data, n1, n2, n3, ix1, ix2, ix3, ix4, h1, h2, h3, h4, ssc, ttc, xface, stif, xp, yp, zp, fni)
Definition i3fri3.F:40
i5fri3
subroutine i5fri3(output, lft, llt, nft, ipari, x, e, irect, msr, nsv, irtl, nty, cst, irtlo, fric0, fric, imast, fsav, fskyi, isky, fcont, v, cf, frot_p, freq, ftsav, ftcont, h3d_data, n1, n2, n3, ix1, ix2, ix3, ix4, xp, yp, zp, ssc, ttc, xface, stif, h1, h2, h3, h4, area, fni)
Definition i3fri3.F:428
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
h3d_mod
Definition h3d_mod.F:311
h3d_mod::h3d_database
Definition h3d_mod.F:438