C DYNEPA SOURCE CHAT 05/01/12 23:18:27 5004 SUBROUTINE DYNEPA(IPRIG,IPPT,IDDL,NDDL,DT, & ICHT41,ICHT43, & ICHT51,ICHT52,ICHT53,ICHT12,ICHT13,ICHT14, & ICHT61,ICHT62,ICHT63,ICHT64,ICHT65,ICHT66) IMPLICIT INTEGER(I-N) IMPLICIT REAL*8(A-H,O-Z) ************************************************************************ * * D Y N E P A * ----------- * * FONCTION: * --------- * UTILITAIRE DE L'OPERATEUR "DYNE" * CALCUL D'UN PAS DE L'ALGORITHME DE DYNAMIQUE EXPLICITE * * PARAMETRES: (E)=ENTREE (S)=SORTIE (+ = CONTENU DANS UN COMMUN) * ----------- * * IPRIG (E) TABLEAU DE POINTEURS SUR DES RIGIDITES, AVEC : * IPRIG(1) POINTEUR SUR L'OBJET M * IPRIG(2) POINTEUR SUR L'OBJET K * IPRIG(3) POINTEUR SUR L'OBJET A * IPRIG(4) POINTEUR SUR L'OBJET Z * IPPT (E) TABLEAU DE POINTEURS, AVEC : -1 * IPPT(1) POINTEUR SUR LE VECTDOUB CONTENANT M * ( ACTIF ) * IPPT(2) POINTEUR SUR UN SEGMENT MBLO ( ACTIF ) * IPPT(3) POINTEUR SUR UN SEGMENT MBLO ( ACTIF ) * IPPT(4) POINTEUR SUR LE SEGMENT MUCTRA ( ACTIF ) * IPPT(5) POINTEUR SUR LE SEGMENT MMATRI ( ACTIF ) * IPPT(6) POINTEUR SUR LE SEGMENT MZZZ ( ACTIF ) * VCHT1 (E/S) VECTEUR CONTENANT B2 * VCHT3 (E/S) VECTEUR CONTENANT F1 * VEC. (E/S) VECTEURS DE TRAVAIL * VCH1 (E/S) VECTEUR CONTENANT U1 * VCH2 (E/S) VECTEUR CONTENANT V0 * VCH3 (E/S) VECTEUR CONTENANT G0 * VCH4 (E/S) VECTEUR CONTENANT U2 * VCH5 (E/S) VECTEUR CONTENANT V1 * VCH6 (E/S) VECTEUR CONTENANT G1 * * AUTEUR, DATE DE CREATION: * ------------------------- * DENIS ROBERT, LE 21 OCTOBRE 1988. * REVU ET CORRIGE DECEMBRE 89 MP * * LANGAGE: * -------- * ESOPE + FORTRAN77 * ************************************************************************ * * MODULES UTILISES: * -INC CCOPTIO * SEGMENT,MCHT1 REAL*8 VCC1(IDDL) ENDSEGMENT POINTEUR MCHT12.MCHT1,MCHT13.MCHT1,MCHT14.MCHT1 * SEGMENT,MCHT2 REAL*8 VCC2(NDDL) ENDSEGMENT POINTEUR MCHT41.MCHT2,MCHT42.MCHT2,MCHT43.MCHT2 POINTEUR MCHT51.MCHT2,MCHT52.MCHT2,MCHT53.MCHT2 POINTEUR MCHT61.MCHT2,MCHT62.MCHT2,MCHT63.MCHT2 POINTEUR MCHT64.MCHT2,MCHT65.MCHT2,MCHT66.MCHT2 * INTEGER IPRIG(4),IPPT(6) * REAL*8 VEC1(*),VEC2(*),VEC3(*),VEC4(*),VEC5(*),VEC6(*) * REAL*8 VCH1(*),VCH2(*),VCH3(*),VCH4(*),VCH5(*),VCH6(*) * REAL*8 VCHT1(*),VCHT3(*) * PARAMETER ( ZERO=0.D0 ) PARAMETER ( O2=2.D0 ) PARAMETER ( O1=-1.D0 ) * MCHT41=ICHT41 MCHT43=ICHT43 MCHT51=ICHT51 MCHT52=ICHT52 MCHT53=ICHT53 MCHT12=ICHT12 MCHT13=ICHT13 MCHT14=ICHT14 MCHT61=ICHT61 MCHT62=ICHT62 MCHT63=ICHT63 MCHT64=ICHT64 MCHT65=ICHT65 MCHT66=ICHT66 * DT2=DT**2 DT3=DT2/O2 DT4=DT/O2 DT5=O1/DT2 DT6=O1/DT * IF (IPRIG(3).NE.0) THEN * *--- L'OBJET RIGIDITE "A" EXISTE * *--------------------------------------------------------------------* * CALCUL DE : Z L1 = A * M-1 * ( F1 - K * U1 ) - 1/DT2 * B2 * * + A * ( 1/DT2 * U1 + 1/DT * V0 + 1/2 * G0 ) * *--------------------------------------------------------------------* * *--- CALCUL DE : - ( 1/DT2 * U1 + 1/DT * V0 + 1/2 * G0 ) * CALL COLIVE(MCHT61.VCC2,MCHT62.VCC2,NDDL,DT5,DT6, & MCHT51.VCC2) CALL COLIVE(MCHT51.VCC2,MCHT63.VCC2,NDDL,1.D0,-0.5D0, & MCHT52.VCC2) * *--- CALCUL DE : A * - ( 1/DT2 * U1 + 1/DT * V0 + 1/2 * G0 ) * CALL DYNE06(IPPT(3),ICHT52,ICHT12) IF (IERR.NE.0) RETURN * *--- CALCUL DE : - K * U1 * DO 10 I10=1,NDDL MCHT51.VCC2(I10)=ZERO MCHT52.VCC2(I10)=ZERO 10 CONTINUE * END DO * CALL DYNE05(IPRIG(2),ICHT61,ICHT51,NDDL,IPPT(4)) IF (IERR.NE.0) RETURN * *--- CALCUL DE : F1 - K * U1 * CALL COLIVE(MCHT43.VCC2,MCHT51.VCC2,NDDL,1.D0,1.D0, & MCHT52.VCC2) * *--- CALCUL DE : A * M-1 * ( F1 - K * U1 ) * CALL DYNE06(IPPT(2),ICHT52,ICHT13) * *--- CALCUL DE : C1 / NE PAS OUBLIER 1/DT2 * B2 / * CALL COLIVE(MCHT12.VCC1,MCHT13.VCC1,IDDL,-1.D0,1.D0, & MCHT14.VCC1) * DO 15 I15=1,IDDL MCHT13.VCC1(I15)=ZERO 15 CONTINUE * END DO CALL COLIVE(MCHT14.VCC1,MCHT41.VCC2,IDDL,1.D0,DT5, & MCHT13.VCC1) * *--- RESOLUTION : Z L1 = C1 * DO 20 I20=1,IDDL MCHT12.VCC1(I20)=ZERO 20 CONTINUE * END DO * CALL DYNERE(ICHT13,IPPT(6),ICHT12,IDDL) * *--------------------------------------------------------------------* *--- CALCUL DE : G1 = M-1 * ( F1 - K * U1 - AT * L1 ) ---* *--------------------------------------------------------------------* * DO 30 I30=1,NDDL MCHT51.VCC2(I30)=ZERO MCHT53.VCC2(I30)=ZERO 30 CONTINUE * END DO * *--- CALCUL DE : AT * L1 * CALL DYNE08(ICHT12,IDDL,NDDL,IPPT(3),ICHT51) IF (IERR.NE.0) RETURN * *--- CALCUL DE : ( F1 - K * U1 ) - ( AT * L1 ) * CALL COLIVE(MCHT52.VCC2,MCHT51.VCC2,NDDL,1.D0,-1.D0, & MCHT53.VCC2) * ELSE * *--- L'OBJET RIGIDITE "A" N'EXISTE PAS * CALL VTRANS(MCHT53.VCC2,MCHT52.VCC2,NDDL) ENDIF * *--- CALCUL DE : M-1 * 'TOUT LE RESTE' * CALL DYNE10(ICHT53,IPPT(1),ICHT66,NDDL) * *--- CALCUL DE : V1 = V0 + ( DT/2 * ( G0 + G1 ) ) *--- CALCUL DE : U2 = U1 + ( DT * V1 ) + ( DT2/2 * G1 ) * DO 40 I40=1,NDDL MCHT51.VCC2(I40)=ZERO MCHT52.VCC2(I40)=ZERO 40 CONTINUE * END DO * CALL COLIVE(MCHT63.VCC2,MCHT66.VCC2,NDDL,1.D0,1.D0, & MCHT51.VCC2) CALL COLIVE(MCHT62.VCC2,MCHT51.VCC2,NDDL,1.D0,DT4, & MCHT65.VCC2) * CALL COLIVE(MCHT65.VCC2,MCHT66.VCC2,NDDL,DT,DT3, & MCHT52.VCC2) CALL COLIVE(MCHT61.VCC2,MCHT52.VCC2,NDDL,1.D0,1.D0, & MCHT64.VCC2) * END