OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
i14frt.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!|| i14frt ../engine/source/interfaces/int14/i14frt.F
25!||--- called by ------------------------------------------------------
26!|| i14cmp ../engine/source/interfaces/int14/i14cmp.F
27!||--- calls -----------------------------------------------------
28!|| ninterp ../engine/source/interfaces/int14/ninterp.F
29!||--- uses -----------------------------------------------------
30!|| anim_mod ../common_source/modules/output/anim_mod.F
31!|| groupdef_mod ../common_source/modules/groupdef_mod.F
32!|| h3d_mod ../engine/share/modules/h3d_mod.F
33!|| output_mod ../common_source/modules/output/output_mod.F90
34!||====================================================================
35 SUBROUTINE i14frt(OUTPUT,AF ,X ,V ,KSURF ,IGRSURF,
36 2 BUFSF ,NSC ,KSC ,NSP ,KSP ,
37 3 KSI ,IMPACT ,CIMP ,NIMP ,FRIC ,
38 4 NFRIC ,NPC ,PLD ,GAPMIN ,STF ,
39 5 WF ,WST ,DE ,MS ,STIFN ,
40 6 FS ,FCONT ,FSKYI ,ISKY ,H3D_DATA)
41C-----------------------------------------------
42C M o d u l e s
43C-----------------------------------------------
44 USE h3d_mod
45 USE groupdef_mod
46 USE anim_mod
47 USE output_mod
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 "com04_c.inc"
61#include "com06_c.inc"
62#include "com08_c.inc"
63#include "scr07_c.inc"
64#include "scr14_c.inc"
65#include "scr16_c.inc"
66#include "parit_c.inc"
67#include "param_c.inc"
68#include "scr18_c.inc"
69C-----------------------------------------------
70C D u m m y A r g u m e n t s
71C-----------------------------------------------
72 TYPE(output_) ,INTENT(INOUT) :: OUTPUT
73 INTEGER NSC,NSP, KSURF, KSI(*), IMPACT(*),
74 . NFRIC, NPC(*),ISKY(*)
75C REAL
76 my_real
77 . af(*), x(3,*), v(3,*) , bufsf(*), ksp(*),
78 . ksc(*), fric, cimp(3,*),nimp(3,*), ms(*),
79 . stifn(*), fs(nthvki),fcont(3,*),fskyi(lskyi,nfskyi), pld(*),
80 . wf(*), wst(*),gapmin , stf , de
81 TYPE(h3d_database) :: H3D_DATA
82 TYPE (SURF_) , DIMENSION(NSURF) :: IGRSURF
83C-----------------------------------------------
84C L o c a l V a r i a b l e s
85C-----------------------------------------------
86 INTEGER ADRBUF, I, IL, J3, J2, J1, IN, I3, I2, I1
87 INTEGER NISKYL
88 INTEGER NDEB, NREST
89 INTEGER DGR
90 INTEGER IPT,NPT,II,JJ
91 my_real
92 . a, b, c, an, bn, cn, rot(9),
93 . x1, x2, x3, n1, n2, n3, n,
94 . xpvn1, ypvn1, zpvn1, sgnxn, sgnyn, sgnzn,
95 . xi, yi, zi,
96 . fxn, fyn, fzn, fxt, fyt, fzt, ftp, ftpa, fmax, fn, tn,
97 . cost, dist, psca, e, fx, fy, fz,
98 . fn1, fn2, fn3, ft1, ft2, ft3, am1, am2, am3,
99 . dd, dt1inv
100 my_real
101 . xpv(3,mvsiz),nv(3,mvsiz) ,tv(3,mvsiz),
102 . fnpv(mvsiz) ,kfric(mvsiz),xq(3,mvsiz),
103 . fv(3,mvsiz) ,st(mvsiz)
104C-----------------------------------------------
105 adrbuf=igrsurf(ksurf)%IAD_BUFR
106C-----------------------------------------------
107 a =bufsf(adrbuf+1)
108 b =bufsf(adrbuf+2)
109 c =bufsf(adrbuf+3)
110 dgr=bufsf(adrbuf+36)
111 an=a**dgr
112 bn=b**dgr
113 cn=c**dgr
114 an=1./an
115 bn=1./bn
116 cn=1./cn
117 DO i=1,9
118 rot(i)=bufsf(adrbuf+7+i-1)
119 END DO
120C-----------------------------------------------
121 IF (dt1==zero) THEN
122 dt1inv=zero
123 ELSE
124 dt1inv=one/dt1
125 END IF
126C-----------------------------------------------
127 IF (iparit/=0) THEN
128#include "lockon.inc"
129 niskyl = nisky
130 nisky = nisky + nsc + nsp
131#include "lockoff.inc"
132 END IF
133C-----------------------------------------------
134 de=zero
135C-----------------------------------------------
136C POINTS JUSTE IMPACTES.
137C-----------------------------------------------
138 ndeb =0
139 nrest=nsc
140 50 CONTINUE
141C-----------------------------------------------
142 DO i=1,min(mvsiz,nrest)
143 il=ksc(i+ndeb)
144 in=ksi(il)
145 fxn=wf(in)*nimp(1,il)
146 fyn=wf(in)*nimp(2,il)
147 fzn=wf(in)*nimp(3,il)
148C------------------------------------------------------------
149C RETOUR DES FORCES EN GLOBAL
150C------------------------------------------------------------
151 fn1 =rot(1)*fxn+rot(4)*fyn+rot(7)*fzn
152 fn2 =rot(2)*fxn+rot(5)*fyn+rot(8)*fzn
153 fn3 =rot(3)*fxn+rot(6)*fyn+rot(9)*fzn
154 fv(1,i) =fn1
155 fv(2,i) =fn2
156 fv(3,i) =fn3
157C---------------------------------
158 fs(1)=fs(1)-fn1*dt1
159 fs(2)=fs(2)-fn2*dt1
160 fs(3)=fs(3)-fn3*dt1
161C---------------------------------
162 ENDDO
163C----------------------------------------------------------
164C MOMENTS ET STABILITE.
165C----------------------------------------------------------
166 IF(kdtint==0)THEN
167 DO i=1,min(mvsiz,nrest)
168 il=ksc(i+ndeb)
169 in=ksi(il)
170 IF(impact(il)>0)THEN
171 st(i)=stf+two*wst(in)*dt1inv
172 ELSE
173 st(i)=zero
174 ENDIF
175 ENDDO
176 ELSE
177 DO i=1,min(mvsiz,nrest)
178 il=ksc(i+ndeb)
179 in=ksi(il)
180 IF(impact(il)>0)THEN
181 st(i)=stf+wst(in)*dt1inv
182 ELSE
183 st(i)=zero
184 ENDIF
185 ENDDO
186 ENDIF
187C----------------------------------------------------------
188 DO i=1,min(mvsiz,nrest)
189 il=ksc(i+ndeb)
190 in=ksi(il)
191C---------------------------------
192C CALCUL DES MOMENTS MS ^ F (REP. GLOBAL).
193C les forces sont transmises au (meme) point
194C que l'on a recu !!!
195C---------------------------------
196 x1=x(1,in)-bufsf(adrbuf+16)
197 x2=x(2,in)-bufsf(adrbuf+17)
198 x3=x(3,in)-bufsf(adrbuf+18)
199 am1=x2*fv(3,i)-x3*fv(2,i)
200 am2=x3*fv(1,i)-x1*fv(3,i)
201 am3=x1*fv(2,i)-x2*fv(1,i)
202C---------------------------------
203C Assemblage des FORCES, moments et rigidites d'interface
204C au nd main.
205C---------------------------------
206 bufsf(adrbuf+25)=bufsf(adrbuf+25)-fv(1,i)
207 bufsf(adrbuf+26)=bufsf(adrbuf+26)-fv(2,i)
208 bufsf(adrbuf+27)=bufsf(adrbuf+27)-fv(3,i)
209 bufsf(adrbuf+28)=bufsf(adrbuf+28)-am1
210 bufsf(adrbuf+29)=bufsf(adrbuf+29)-am2
211 bufsf(adrbuf+30)=bufsf(adrbuf+30)-am3
212 bufsf(adrbuf+31)=bufsf(adrbuf+31)+st(i)
213C STIFRM += [ ST2 * |(CIMP-M)^N|**2 ] + [ ST1 * |(CIMP-M)^T|**2 ]
214C est majore par :
215 dd = x1**2+x2**2+x3**2
216 bufsf(adrbuf+32)=bufsf(adrbuf+32)+dd*st(i)
217C---------------------------------
218 ENDDO
219C----------------------------------------------------------
220 IF(kdtint/=0)THEN
221 DO i=1,min(mvsiz,nrest)
222 il=ksc(i+ndeb)
223 in=ksi(il)
224 IF(impact(il)>0)THEN
225 st(i)=stf
226 ELSE
227 st(i)=zero
228 ENDIF
229 ENDDO
230 ENDIF
231C---------------------------------
232C Assemblage des FORCES aux noeuds seconds.
233C---------------------------------
234 IF (iparit==0) THEN
235C-----------------------------------------------
236#include "vectorize.inc"
237 DO i=1,min(mvsiz,nrest)
238 il=ksc(i+ndeb)
239 in=ksi(il)
240 i3=3*in
241 i2=i3-1
242 i1=i2-1
243 af(i1)=af(i1)+fv(1,i)
244 af(i2)=af(i2)+fv(2,i)
245 af(i3)=af(i3)+fv(3,i)
246C Stabilite
247 stifn(in)=stifn(in)+st(i)
248 ENDDO
249 ELSE
250C-----------------------------------------------
251 IF(kdtint==0)THEN
252 DO i=1,min(mvsiz,nrest)
253 il=ksc(i+ndeb)
254 in=ksi(il)
255 niskyl = niskyl + 1
256 fskyi(niskyl,1)=fv(1,i)
257 fskyi(niskyl,2)=fv(2,i)
258 fskyi(niskyl,3)=fv(3,i)
259C Stabilite
260 fskyi(niskyl,4)=st(i)
261 isky(niskyl) =in
262 ENDDO
263 ELSE
264 DO i=1,min(mvsiz,nrest)
265 il=ksc(i+ndeb)
266 in=ksi(il)
267 niskyl = niskyl + 1
268 fskyi(niskyl,1)=fv(1,i)
269 fskyi(niskyl,2)=fv(2,i)
270 fskyi(niskyl,3)=fv(3,i)
271C Stabilite
272 fskyi(niskyl,4)=st(i)
273 fskyi(niskyl,5)=wst(in)
274C 2C dans modsti FSKYI(NISKYL,5)=2.*WST(IN)
275 isky(niskyl) =in
276 ENDDO
277 ENDIF
278
279 ENDIF
280C------------------------------------------------------------
281C ANIM (FORCES DE CONTACT).
282C------------------------------------------------------------
283 IF(anim_v(4)+outp_v(4)+h3d_data%N_VECT_CONT >0.AND.
284 . ((tt>=tanim .AND. tt<=tanim_stop).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP) .OR.
285 . (manim>=4.AND.manim<=15).OR. h3d_data%MH3D /= 0))THEN
286#include "lockon.inc"
287#include "vectorize.inc"
288 DO i=1,min(mvsiz,nrest)
289 il=ksc(i+ndeb)
290 in=ksi(il)
291 fcont(1,in) =fcont(1,in) + fv(1,i)
292 fcont(2,in) =fcont(2,in) + fv(2,i)
293 fcont(3,in) =fcont(3,in) + fv(3,i)
294 ENDDO
295#include "lockoff.inc"
296 ENDIF
297C---------------------------------
298C Pour Travail des forces sur noeuds seconds
299C 1ere partie : ici
300C 2eme partie : apres calcul de DT2.
301C---------------------------------
302 DO i=1,min(mvsiz,nrest)
303 il=ksc(i+ndeb)
304 in=ksi(il)
305 de=de+fv(1,i)*v(1,in)+fv(2,i)*v(2,in)+fv(3,i)*v(3,in)
306 ENDDO
307C---------------------------------
308C Groupe suivant.
309C---------------------------------
310 IF (nrest-mvsiz>0) THEN
311 nrest=nrest-mvsiz
312 ndeb =ndeb +mvsiz
313 GOTO 50
314 END IF
315C-------------------------------
316C POINTS PRECEDEMMENT IMPACTES.
317C-------------------------------
318 ndeb =0
319 nrest=nsp
320 100 CONTINUE
321C-------------------------------
322C Passage au repere local :
323C-------------------------------
324 DO i=1,min(mvsiz,nrest)
325 il=ksp(i+ndeb)
326 in=ksi(il)
327 x1=x(1,in)-bufsf(adrbuf+4)
328 x2=x(2,in)-bufsf(adrbuf+5)
329 x3=x(3,in)-bufsf(adrbuf+6)
330 xpv(1,i)=rot(1)*x1+rot(2)*x2+rot(3)*x3
331 xpv(2,i)=rot(4)*x1+rot(5)*x2+rot(6)*x3
332 xpv(3,i)=rot(7)*x1+rot(8)*x2+rot(9)*x3
333 ENDDO
334C-------------------------------
335C Extraction des infos deja calculees.
336C-------------------------------
337 DO i=1,min(mvsiz,nrest)
338 il=ksp(i+ndeb)
339 in=ksi(il)
340 fnpv(i) =wf(in)
341 nv(1,i) =nimp(1,il)
342 nv(2,i) =nimp(2,il)
343 nv(3,i) =nimp(3,il)
344 ENDDO
345C----------------------------------------------------------
346C Coefficient de Frottement :
347C MU=F(FORCE NORMALE LOCALE).
348C----------------------------------------------------------
349 IF (nfric==0) THEN
350 DO i=1,min(mvsiz,nrest)
351 kfric(i)=fric
352 ENDDO
353 ELSE
354 CALL ninterp(nfric,npc,pld,min(mvsiz,nrest),fnpv,kfric)
355 DO i=1,min(mvsiz,nrest)
356 kfric(i)=fric*kfric(i)
357 ENDDO
358 ENDIF
359C----------------------------------------------------------
360C Frottement :
361C----------------------------------------------------------
362 DO i=1,min(mvsiz,nrest)
363 il=ksp(i+ndeb)
364 in=ksi(il)
365C------------------------------
366 fxt=zero
367 fyt=zero
368 fzt=zero
369C-------------------------------
370 xi=cimp(1,il)
371 yi=cimp(2,il)
372 zi=cimp(3,il)
373C------------------------------
374 fmax= max(kfric(i)*fnpv(i),zero)
375C Force tangente K(XI-XP) et Norme.
376 fxt=(xi-xpv(1,i))*stf
377 fyt=(yi-xpv(2,i))*stf
378 fzt=(zi-xpv(3,i))*stf
379 ftp=sqrt(fxt*fxt+fyt*fyt+fzt*fzt)
380C-----------------------------------------------
381 fxn=fnpv(i)*nv(1,i)
382 fyn=fnpv(i)*nv(2,i)
383 fzn=fnpv(i)*nv(3,i)
384C-----------------------------------------------
385 fn=fxt*nv(1,i)+fyt*nv(2,i)+fzt*nv(3,i)
386 tv(1,i)=fxt-nv(1,i)*fn
387 tv(2,i)=fyt-nv(2,i)*fn
388 tv(3,i)=fzt-nv(3,i)*fn
389 tn=sqrt(tv(1,i)*tv(1,i)+tv(2,i)*tv(2,i)+tv(3,i)*tv(3,i))
390 tn=max(em20,tn)
391 tv(1,i)=tv(1,i)/tn
392 tv(2,i)=tv(2,i)/tn
393 tv(3,i)=tv(3,i)/tn
394C-----------------------------------------------
395C Force de Frottement.
396 ftpa=min(tn,fmax)
397 fxt=tv(1,i)*ftpa
398 fyt=tv(2,i)*ftpa
399 fzt=tv(3,i)*ftpa
400C-------------------------------
401C proj. sur tangente.
402 cost =(xpv(1,i)-xi)*tv(1,i)
403 . +(xpv(2,i)-yi)*tv(2,i)
404 . +(xpv(3,i)-zi)*tv(3,i)
405 xq(1,i) =xi+cost*tv(1,i)
406 xq(2,i) =yi+cost*tv(2,i)
407 xq(3,i) =zi+cost*tv(3,i)
408C-------------------------------
409C glissement du point sur la tgte.
410 xq(1,i) =xq(1,i)+ftpa*tv(1,i)/max(em20,stf)
411 xq(2,i) =xq(2,i)+ftpa*tv(2,i)/max(em20,stf)
412 xq(3,i) =xq(3,i)+ftpa*tv(3,i)/max(em20,stf)
413C---------------------------------
414 fv(1,i)=fxt+fxn
415 fv(2,i)=fyt+fyn
416 fv(3,i)=fzt+fzn
417C---------------------------------
418 ENDDO
419C----------------------------------------------------------
420C Frottement Suite :
421C----------------------------------------------------------
422#include "vectorize.inc"
423 DO i=1,min(mvsiz,nrest)
424 il=ksp(i+ndeb)
425 in=ksi(il)
426C-------------------------------
427C nouvelle normale.
428C-------------------------------
429 xpvn1 =xq(1,i)**(dgr-1)
430 sgnxn=-one
431 IF (xpvn1*xq(1,i)>=zero) sgnxn=one
432 ypvn1 =xq(2,i)**(dgr-1)
433 sgnyn=-one
434 IF (ypvn1*xq(2,i)>=zero) sgnyn=one
435 zpvn1 =xq(3,i)**(dgr-1)
436 sgnzn=-one
437 IF (zpvn1*xq(3,i)>=zero) sgnzn=one
438 n1 =sgnxn*xpvn1*an
439 n2 =sgnyn*ypvn1*bn
440 n3 =sgnzn*zpvn1*cn
441 n =n1*n1+n2*n2+n3*n3
442 n =sqrt(n)
443C-------------------------------
444C nouveau point d'impact : projection hors ellips.
445C la position de l'impact et le calcul de la normale s'affinent
446C en l'absence de glissement.
447C-------------------------------
448 e =n1*xq(1,i)+n2*xq(2,i)+n3*xq(3,i)
449 dist =(e-exp((dgr-1)*log(max(e,em20))/dgr))
450 . / max(exp((dgr-1)*log(em20)/dgr),n)
451 n1 =n1/max(em20,n)
452 n2 =n2/max(em20,n)
453 n3 =n3/max(em20,n)
454 cimp(1,il)=xq(1,i)-dist*n1
455 cimp(2,il)=xq(2,i)-dist*n2
456 cimp(3,il)=xq(3,i)-dist*n3
457C memorisation des composantes de la normale.
458 nimp(1,il)=n1
459 nimp(2,il)=n2
460 nimp(3,il)=n3
461C-------------------------------
462C rabattre F // PI apres glissement
463C ROMPT LA RELATION FTOT=F(PENETRATION MAX.)
464C-------------------------------
465 fx=fv(1,i)
466 fy=fv(2,i)
467 fz=fv(3,i)
468C-------------------------------
469 xi=cimp(1,il)+gapmin*nimp(1,il)
470 yi=cimp(2,il)+gapmin*nimp(2,il)
471 zi=cimp(3,il)+gapmin*nimp(3,il)
472 psca = (xi-xpv(1,i))*fx
473 . +(yi-xpv(2,i))*fy
474 . +(zi-xpv(3,i))*fz
475 dist = (xi-xpv(1,i))**2
476 . +(yi-xpv(2,i))**2
477 . +(zi-xpv(3,i))**2
478 fx=psca/dist*(xi-xpv(1,i))
479 fy=psca/dist*(yi-xpv(2,i))
480 fz=psca/dist*(zi-xpv(3,i))
481C decomposition pour sorties :
482 psca=fx*nv(1,i)+fy*nv(2,i)+fz*nv(3,i)
483 fxn =psca*nv(1,i)
484 fyn =psca*nv(2,i)
485 fzn =psca*nv(3,i)
486 psca=fx*tv(1,i)+fy*tv(2,i)+fz*tv(3,i)
487 fxt =psca*tv(1,i)
488 fyt =psca*tv(2,i)
489 fzt =psca*tv(3,i)
490C------------------------------------------------------------
491C RETOUR DES FORCES EN GLOBAL
492C------------------------------------------------------------
493 fn1 =rot(1)*fxn+rot(4)*fyn+rot(7)*fzn
494 fn2 =rot(2)*fxn+rot(5)*fyn+rot(8)*fzn
495 fn3 =rot(3)*fxn+rot(6)*fyn+rot(9)*fzn
496 ft1 =rot(1)*fxt+rot(4)*fyt+rot(7)*fzt
497 ft2 =rot(2)*fxt+rot(5)*fyt+rot(8)*fzt
498 ft3 =rot(3)*fxt+rot(6)*fyt+rot(9)*fzt
499 fv(1,i) =fn1+ft1
500 fv(2,i) =fn2+ft2
501 fv(3,i) =fn3+ft3
502C---------------------------------
503 fs(1)=fs(1)-fn1*dt1
504 fs(2)=fs(2)-fn2*dt1
505 fs(3)=fs(3)-fn3*dt1
506 fs(4)=fs(4)-ft1*dt1
507 fs(5)=fs(5)-ft2*dt1
508 fs(6)=fs(6)-ft3*dt1
509C---------------------------------
510 ENDDO
511C----------------------------------------------------------
512C MOMENTS ET STABILITE.
513C----------------------------------------------------------
514 IF(kdtint==0)THEN
515 DO i=1,min(mvsiz,nrest)
516 il=ksp(i+ndeb)
517 in=ksi(il)
518 st(i)=stf+two*wst(in)*dt1inv
519 ENDDO
520 ELSE
521 DO i=1,min(mvsiz,nrest)
522 il=ksp(i+ndeb)
523 in=ksi(il)
524 st(i)=stf+wst(in)*dt1inv
525 ENDDO
526 ENDIF
527C----------------------------------------------------------
528 DO i=1,min(mvsiz,nrest)
529 il=ksp(i+ndeb)
530 in=ksi(il)
531C---------------------------------
532C CALCUL DES MOMENTS MS ^ F (REP. GLOBAL).
533C les forces sont transmises au (meme) point
534C que l'on a recu !!!
535C---------------------------------
536 x1=x(1,in)-bufsf(adrbuf+16)
537 x2=x(2,in)-bufsf(adrbuf+17)
538 x3=x(3,in)-bufsf(adrbuf+18)
539 am1=x2*fv(3,i)-x3*fv(2,i)
540 am2=x3*fv(1,i)-x1*fv(3,i)
541 am3=x1*fv(2,i)-x2*fv(1,i)
542C---------------------------------
543C Assemblage des FORCES, moments et rigidites d'interface
544C au nd main.
545C---------------------------------
546 bufsf(adrbuf+25)=bufsf(adrbuf+25)-fv(1,i)
547 bufsf(adrbuf+26)=bufsf(adrbuf+26)-fv(2,i)
548 bufsf(adrbuf+27)=bufsf(adrbuf+27)-fv(3,i)
549 bufsf(adrbuf+28)=bufsf(adrbuf+28)-am1
550 bufsf(adrbuf+29)=bufsf(adrbuf+29)-am2
551 bufsf(adrbuf+30)=bufsf(adrbuf+30)-am3
552 bufsf(adrbuf+31)=bufsf(adrbuf+31)+st(i)
553C STIFRM += [ ST2 * |(CIMP-M)^N|**2 ] + [ ST1 * |(CIMP-M)^T|**2 ]
554C est majore par :
555 dd = x1**2+x2**2+x3**2
556 bufsf(adrbuf+32)=bufsf(adrbuf+32)+dd*st(i)
557C---------------------------------
558 ENDDO
559C---------------------------------
560C Assemblage des FORCES aux noeuds seconds.
561C---------------------------------
562 IF (iparit==0) THEN
563C-----------------------------------------------
564#include "vectorize.inc"
565 DO i=1,min(mvsiz,nrest)
566 il=ksp(i+ndeb)
567 in=ksi(il)
568 i3=3*in
569 i2=i3-1
570 i1=i2-1
571 af(i1)=af(i1)+fv(1,i)
572 af(i2)=af(i2)+fv(2,i)
573 af(i3)=af(i3)+fv(3,i)
574C Stabilite
575 stifn(in)=stifn(in)+st(i)
576 ENDDO
577 ELSE
578C-----------------------------------------------
579 IF(kdtint==0)THEN
580 DO i=1,min(mvsiz,nrest)
581 il=ksp(i+ndeb)
582 in=ksi(il)
583 niskyl = niskyl + 1
584 fskyi(niskyl,1)=fv(1,i)
585 fskyi(niskyl,2)=fv(2,i)
586 fskyi(niskyl,3)=fv(3,i)
587C Stabilite
588 fskyi(niskyl,4)=st(i)
589 isky(niskyl) =in
590 ENDDO
591 ELSE
592 DO i=1,min(mvsiz,nrest)
593 il=ksp(i+ndeb)
594 in=ksi(il)
595 niskyl = niskyl + 1
596 fskyi(niskyl,1)=fv(1,i)
597 fskyi(niskyl,2)=fv(2,i)
598 fskyi(niskyl,3)=fv(3,i)
599C Stabilite
600 fskyi(niskyl,4)=stf
601 fskyi(niskyl,5)=wst(in)
602C 2C dans modsti FSKYI(NISKYL,5)=2.*WST(IN)
603 isky(niskyl) =in
604 ENDDO
605 ENDIF
606
607 ENDIF
608C------------------------------------------------------------
609C ANIM (FORCES DE CONTACT).
610C------------------------------------------------------------
611 IF(anim_v(4)+outp_v(4)+h3d_data%N_VECT_CONT>0.AND.
612 . (tt>=tanim .AND. tt<=tanim_stop.OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP).OR.
613 . (manim>=4.AND.manim<=15).OR. h3d_data%MH3D /= 0))THEN
614#include "lockon.inc"
615#include "vectorize.inc"
616 DO i=1,min(mvsiz,nrest)
617 il=ksp(i+ndeb)
618 in=ksi(il)
619 fcont(1,in) =fcont(1,in) + fv(1,i)
620 fcont(2,in) =fcont(2,in) + fv(2,i)
621 fcont(3,in) =fcont(3,in) + fv(3,i)
622 ENDDO
623#include "lockoff.inc"
624 ENDIF
625C---------------------------------
626C Pour Travail des forces sur noeuds seconds
627C 1ere partie : ici
628C 2eme partie : apres calcul de DT2.
629C---------------------------------
630 DO i=1,min(mvsiz,nrest)
631 il=ksp(i+ndeb)
632 in=ksi(il)
633 de=de+fv(1,i)*v(1,in)+fv(2,i)*v(2,in)+fv(3,i)*v(3,in)
634 ENDDO
635C---------------------------------
636C Groupe suivant.
637C---------------------------------
638 IF (nrest-mvsiz>0) THEN
639 nrest=nrest-mvsiz
640 ndeb =ndeb +mvsiz
641 GOTO 100
642 END IF
643C---------------------------------
644C Working force at interface (Madymo)
645C---------------------------------
646 fs(7)=fs(7)+de*dt1*half
647 IF (igrsurf(ksurf)%TYPE==100) THEN
648C Madymo Ellipsoids.
649!$OMP ATOMIC
650 output%TH%WFEXT=output%TH%WFEXT+de*dt1*half
651 ENDIF
652C------------------------------------------------------------
653 RETURN
654 END
i14frt
subroutine i14frt(output, af, x, v, ksurf, igrsurf, bufsf, nsc, ksc, nsp, ksp, ksi, impact, cimp, nimp, fric, nfric, npc, pld, gapmin, stf, wf, wst, de, ms, stifn, fs, fcont, fskyi, isky, h3d_data)
Definition i14frt.F:41
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
anim_mod
Definition anim_mod.F:96
groupdef_mod
Definition groupdef_mod.F:662
h3d_mod
Definition h3d_mod.F:308
ninterp
subroutine ninterp(ifunc, npc, pld, npoint, xv, yv)
Definition ninterp.F:32
h3d_mod::h3d_database
Definition h3d_mod.F:435