Sources Esope | Cast3M

dyforc
C DYFORC    SOURCE    BP208322  20/09/18    21:16:20     10718          
        SUBROUTINE DYFORC(XPTB,XVALB,IPALB,IPLIB,XPALB,NLIAB,NPLB,IND
     &,IND2,PDTS2,I,iannul,IESC,IROLE,XCHPFB)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
*--------------------------------------------------------------------*
*                                                                    *
*     Opérateur DYNE : algorithme de Fu - de Vogelaere               *
*     ________________________________________________               *
*                                                                    *
*     Calcul des forces de choc sur base B pour les liaisons de      *
*     type POINT_PLAN.                                               *
*                                                                    *
*     Paramètres:                                                    *
*                                                                    *
* e   XPTB    Tableau des déplacements des points                    *
* e   IPALB   Renseigne sur la liaison.                              *
* e   IPLIB   Tableau contenant les numéros "DYNE" de la liaison.    *
* e   XPALB   Tableau contenant les paramètres de la liaison.        *
* es  XVALB   Tableau contenant les variables internes de liaisons.  *
* es  XCHPFB  Tableau contenant les valeurs des futurs chpoints      *
* e   NLIAB   Nombre de liaisons sur la base B.                      *
* e   NPLB    Nombre total de points intervenant dans les liaisons.  *
* e   IND     Indice du pas.                                         *
* e   I       numéro de la liaison.                                  *
* e   IESC    numéro DYNE du noeud esclave                           *
* e   IROLE   indique qui est le maitre et l'esclave                 *
*             = 0 premiere passe , =1 deuxième passe                 *
*                                                                    *
*                                                                    *
*     Auteur, date de création:                                      *
*                                                                    *
*     Samuel DURAND      : le 19 Octobre 1996  : Création            *
*                                                                    *
*--------------------------------------------------------------------*
*
      INTEGER IPALB(NLIAB,*),IPLIB(NLIAB,*)
      INTEGER ICAND(2)
      REAL*8  XPALB(NLIAB,*),XPTB(NPLB,2,*),XVALB(NLIAB,4,*)
      REAL*8  XXXN(3),XCHPFB(2,NLIAB,4,*)
*
        XEPSI = 0.0001
        IGP    = IPALB(I,2)
        IDIM   = IPALB(I,3)
*       ILOCAL = IPALB(I,23)
        ITYP = IPALB(I,1)
        XPOID = 0.D0
*
        IF (ITYP.EQ.35) THEN
           ID1 = 6
        ELSE
           ID1 = 7
        ENDIF
*
        IF (IROLE.EQ.0) THEN
           KMAI = 0
           IMAI = ID1 +4*IDIM
           IBUT = IMAI + IPALB(I,21)*IDIM
           KBUT = IPALB(I,21)
           LMAI = ID1 + (IPALB(I,21)+IPALB(I,22)+4)*IDIM
           LBUT = LMAI +IPALB(I,21)*IDIM
           JMAI = 26 + IPALB(I,21)+IPALB(I,22)
           JBUT = JMAI + IPALB(I,21)
           MBUT = LBUT +IPALB(I,22)*IDIM +IPALB(I,21)
           NBUT = MBUT +IPALB(I,21)+IPALB(I,22)
        ELSE
           KMAI =  IPALB(I,21)
           IBUT = ID1 + 4*IDIM
           IMAI = IBUT + IPALB(I,21)*IDIM
           KBUT = 0
           LBUT = ID1 + (IPALB(I,21)+IPALB(I,22)+4)*IDIM
           LMAI = LBUT +IPALB(I,21)*IDIM
           JBUT = 26 + 2*(IPALB(I,21)+IPALB(I,22))
           JMAI = JBUT +IPALB(I,21)
           MBUT = LMAI +IPALB(I,22)*IDIM
           NBUT = MBUT +IPALB(I,21)+IPALB(I,22)
        ENDIF
        ID2 = ID1 +  IDIM
        ID4 = ID1 +3*IDIM
        CALL DYCAND(IPALB,IPLIB,XPALB,XPTB,IND,I,NLIAB,NPLB,
     &                  XXXN,XDEP,XPOID,ICAND,IESC,IROLE)
            XVITN =0.D0
            IDCAN1 = IMAI + (ICAND(1)-1)*IDIM
            IDCAN2 = IMAI + (ICAND(2)-1)*IDIM
            DO 80 ID=1,IDIM
*
*   Déplacement du point de contact au pas courant
                XDE2 =(XPTB(IPLIB(I,KMAI+ICAND(2)),1,ID))*(1-XPOID)
     &             +  (XPTB(IPLIB(I,KMAI+ICAND(1)),1,ID))*XPOID
                XVALB(I,IND,3+ID) = XVALB(I,IND,3+ID)+
     &          XPALB(I,IDCAN1+ID)*XPOID +XPALB(I,IDCAN2+ID)*(1-XPOID)
cbp,2020-09*   Déplacement du point de contact au pas précédent
cbp,2020-09                XDM2 =(XPTB(IPLIB(I,KMAI+ICAND(2)),IND2,ID))*(1-XPOID)
cbp,2020-09     &           +  (XPTB(IPLIB(I,KMAI+ICAND(1)),IND2,ID))*XPOID
*   Vitesse du point de contact au pas courant
                XVI2 =(XPTB(IPLIB(I,KMAI+ICAND(2)),2,ID))*(1-XPOID)
     &             +  (XPTB(IPLIB(I,KMAI+ICAND(1)),2,ID))*XPOID
cbp,2020-09*   Déplacement de la butee
cbp,2020-09                XDB = (XPTB(IPLIB(I,KBUT+IESC),IND,ID)
cbp,2020-09     &                -XPTB(IPLIB(I,KBUT+IESC),IND2,ID))
*   Vitesse de la butee
                XVB = XPTB(IPLIB(I,KBUT+IESC),2,ID)
*   Vitesse de glissement
cbp,2020-09     XPALB(I,ID2+ID)=(XDB+XDM2-XDE2)/PDTS2
                XPALB(I,ID2+ID)=XVB-XVI2
*   Vitesse normale
cbp,2020-09     XVITN = XVITN + ((XDB+XDM2-XDE2)*XXXN(ID)/PDTS2)
                XVITN = XVITN + (XVB-XVI2)*XXXN(ID)
 80             CONTINUE
           DO 82 ID=1,IDIM
*   Vitesse tangentielle
                XPALB(I,ID2+ID) = XPALB(I,ID2+ID) - (XVITN*XXXN(ID))
 82        CONTINUE
*   Calcul de la force de choc
                XPALB(I,1) = XPALB(I,MBUT+IESC)
            IF (ITYP.EQ.35) THEN
                CALL DYCHE4(XDEP,IDIM,IGP,XPALB,NLIAB,I,XFN,XFT,XPUS,
     &iannul)
            ELSE
                XPALB(I,7) = XPALB(I,NBUT+IESC)
                CALL DYCHA4(XDEP,XVITN,IDIM,IGP,XPALB,NLIAB,I,XFN,XFT
     &,XPUS,iannul)
            ENDIF
            XVALB(I,IND,3) = XVALB(I,IND,3) + XVITN
            IPALB(I,2) = IGP
            IF (IGP.EQ.1) THEN
                PS =0.D0
                DO 84 ID=1,IDIM
                    PS=PS+(XPALB(I,ID2+ID)*XPALB(I,ID2+ID))
 84             CONTINUE
                XVITT = SQRT(PS)
            ELSE
                XVITT = 0.D0
            ENDIF
            XVALB(I,IND,11) = XVALB(I,IND,11)+XVITT
            IF (XDEP.GT.0) THEN
              DO 90 ID=1,IDIM
                XPALB(I,ID+LMAI+(ICAND(1)-1)*IDIM) =XPALB(I,ID+LMAI+
     &(ICAND(1)-1)*IDIM)-XPOID * (XFN*XXXN(ID) + XPALB(I,ID4+ID))
                XPALB(I,ID+LMAI+(ICAND(2)-1)*IDIM) =XPALB(I,ID+LMAI+
     &(ICAND(2)-1)*IDIM)- ( 1.-XPOID )*(XFN*XXXN(ID) + XPALB(I,ID4+ID))
                XPALB(I,ID+LBUT+(IESC-1)*IDIM) =XPALB(I,ID+LBUT+
     &(IESC-1)*IDIM) + XFN * XXXN(ID) +  XPALB(I,ID4+ID)
  90          CONTINUE
                XCHPFB(1,I,IND,IPLIB(I,KBUT+IESC))=
     &XCHPFB(1,I,IND,IPLIB(I,KBUT+IESC))+XFN
                XCHPFB(2,I,IND,IPLIB(I,KBUT+IESC))=
     &XCHPFB(2,I,IND,IPLIB(I,KBUT+IESC))+ABS(XFT)
                XCHPFB(1,I,IND,IPLIB(I,KMAI+ICAND(1)))=
     &XCHPFB(1,I,IND,IPLIB(I,KMAI+ICAND(1)))+XPOID *XFN
                XCHPFB(2,I,IND,IPLIB(I,KMAI+ICAND(1)))=
     &XCHPFB(2,I,IND,IPLIB(I,KMAI+ICAND(1)))+XPOID *ABS(XFT)
                XCHPFB(1,I,IND,IPLIB(I,KMAI+ICAND(2)))=
     &XCHPFB(1,I,IND,IPLIB(I,KMAI+ICAND(2)))+( 1.-XPOID )*XFN
                XCHPFB(2,I,IND,IPLIB(I,KMAI+ICAND(2)))=
     &XCHPFB(2,I,IND,IPLIB(I,KMAI+ICAND(2)))+( 1.-XPOID )*ABS(XFT)
              IF (XPOID.GT.XEPSI) THEN
                IPALB(I,JMAI+ICAND(1))=1
              ENDIF
              IF ((1.0-XPOID).GT.XEPSI) THEN
                IPALB(I,JMAI+ICAND(2))=1
              ENDIF
                IPALB(I,JBUT+IESC) = 1
                XVALB(I,IND,1) = XVALB(I,IND,1)+XFN
                XVALB(I,IND,10) = XVALB(I,IND,10)+ ABS(XFT)
                XVALB(I,IND,12) = XVALB(I,IND,12)+ XPUS
            ENDIF
*
      END