OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
i14frt.F File Reference
#include "implicit_f.inc"
#include "comlock.inc"
#include "mvsiz_p.inc"
#include "com04_c.inc"
#include "com06_c.inc"
#include "com08_c.inc"
#include "scr07_c.inc"
#include "scr14_c.inc"
#include "scr16_c.inc"
#include "parit_c.inc"
#include "param_c.inc"
#include "scr18_c.inc"
#include "lockon.inc"
#include "lockoff.inc"
#include "vectorize.inc"

Go to the source code of this file.

Functions/Subroutines

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)

Function/Subroutine Documentation

◆ i14frt()

subroutine i14frt ( type(output_), intent(inout) output,
af,
x,
v,
integer ksurf,
type (surf_), dimension(nsurf) igrsurf,
bufsf,
integer nsc,
ksc,
integer nsp,
ksp,
integer, dimension(*) ksi,
integer, dimension(*) impact,
cimp,
nimp,
fric,
integer nfric,
integer, dimension(*) npc,
pld,
gapmin,
stf,
wf,
wst,
de,
ms,
stifn,
fs,
fcont,
fskyi,
integer, dimension(*) isky,
type(h3d_database) h3d_data )

Definition at line 34 of file i14frt.F.

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