1 : C TUFIR1 SOURCE CHAT 05/01/13 03:56:11 5004 2 : SUBROUTINE TUFIR1(REG,RAYON,EPAI,V1,V2,ANDEG,D,REL,BPSS,KERRE) 3 : C======================================================================= 4 : C CALCUL DE LA MATRICE DE RIGIDITE DE L'ELEMENT DE TUYAU FISSURE 5 : C DANS LE REPERE GLOBAL 6 : C ENTREES: 7 : C WORK = TABLEAU CONTENANT LES CARACTERISTIQUES GEOMETRIQUES DE 8 : C L'ELEMENT 9 : C D(2,2)= MATRICE DE HOOKE 10 : C TRAVAIL 11 : C REL(12,12) = MATRICE DE RIGIDITE LOCALE 12 : C V1(3),V2(3) = VECTEURS ORIENTANT LE TUYAU FISSURE 13 : C BPSS(3,3) = MATRICE DE PASSAGE 14 : C SORTIE: 15 : C REG(12,12)= MATRICE DE RIGIDITE GLOBALE 16 : C KERRE = CODE D'ERREUR CORRESPONDANT A LA POSITION RELATIVE DES 17 : C DEUX VECTEURS V1 ET V2 DEFINISSANT L'ELEMENT TUFI: 18 : C 0 SI PAS DE PROBLEME, 1 SINON. 19 : C 1 SI V2 EST NUL OU SI V1 ET V2 SONT COLINEAIRES. 20 : C 2 SI V1 EST NUL. 21 : C 3 SI LE TUYAU FISSURE EST TRP EPAI R/T INF A 2 22 : C======================================================================= 23 : IMPLICIT INTEGER(I-N) 24 : IMPLICIT REAL*8(A-H,O-Z) 25 : C Include contenant quelques constantes dont XPI : 26 : -INC CCREEL 27 : DIMENSION REL(12,*),V1(*),V2(*),D(2,*),BPSS(3,*) 28 : DIMENSION REG(12,*),BV1(3),BV2(3) 29 : C 30 : DATA EPS/1.D-10/ 31 : DATA UN,UNS2/1.D0,.5D0/ 32 : C 33 : C D6 = COEFFICIENT INTERVENANT POUR LES TERMES DE PENALISATION. 34 : C 35 : DATA DEUX,D6/2.D0,1.D06/ 36 : C 37 : KERRE = 0 38 : TETA1 = ANDEG*UNS2 39 : YOU = D(1,1) 40 : XNU = D(2,2) 41 : CISA = YOU /(DEUX *(UN + XNU)) 42 : C======================================================================= 43 : C LES CARACTERISTIQUES REPRESENTENT: 44 : C EPAI = EPAISSEUR DU TUYAU 45 : C RAYON = RAYON EXTERIEUR DU TUYAU 46 : C TETA1 = DEMI-ANGLE D'OUVERTURE DE LA FISSURE EN DEGRE 47 : C V1 = VECTEUR ORIENTANT LES NOEUDS 1 ET 2 (AXE DU TUYAU) 48 : C V2 = VECTEUR ORIENTANT L'AXE DE LA FISSURE 49 : C YOU = MODULE D'YOUNG 50 : C CISA = MODULE DE CISAILLEMENT 51 : C TABLEAU DE TRAVAIL: 52 : C BPSS = MATRICE DE PASSAGE 53 : C ON OBTIENT: 54 : C REG(12,12) = MATRICE DE RIGIDITE AXES GLOBAUX 55 : C 56 : C======================================================================= 57 : C 58 : C CONVERSION DE TETA1 EN RADIAN 59 : C 60 : TETA = (TETA1 * XPI)/180.D0 61 : C 62 : C ON MET DANS 'RAYON' LE RAYON MOYEN DU TUYAU. 63 : C 64 : RAYMO = RAYON - EPAI*UNS2 65 : C 66 : C ON CALCULE A 67 : C 68 : RSURT= RAYMO / EPAI 69 : IF(RSURT.LE.10.0D0.AND.RSURT.GE.4.9D0) THEN 70 : AXX= (.125D0*RSURT - .25D0 ) **.25D0 71 : ELSE IF(RSURT.GT.10.D0.AND.RSURT.LE.35.D0) THEN 72 : AXX= (.4D0*RSURT - 3.D0 ) **.25D0 73 : ELSE 74 : KERRE=3 75 : ENDIF 76 : IF(KERRE.EQ.3) GOTO 666 77 : C 78 : C NORMALISATION DES VECTEURS V1 ET V2. 79 : C 80 : XNV1 = XZERO 81 : XNV2 = XZERO 82 : DO 5 I=1,3 83 : XNV1 = XNV1 + (V1(I)*V1(I)) 84 : XNV2 = XNV2 + (V2(I)*V2(I)) 85 : 5 CONTINUE 86 : C 87 : C MISE A ZERO DE REG 88 : C 88 : CALL ZERO(REG,12,12) 89 : IF (XNV1.LT.EPS) THEN 90 : KERRE = 2 91 : GOTO 666 92 : ENDIF 93 : IF (XNV2.LT.EPS) THEN 94 : KERRE = 1 95 : GOTO 666 96 : ENDIF 97 : XNV1 = UN/SQRT(XNV1) 98 : XNV2 = UN/SQRT(XNV2) 99 : DO 10 I=1,3 100 : BV1(I) = V1(I)*XNV1 101 : BV2(I) = V2(I)*XNV2 102 : 10 CONTINUE 103 : C 104 : C ORTHOGONALISATION ET RENORMALISATION DE V2. 105 : C 106 : PS = (BV1(1)*BV2(1))+(BV1(2)*BV2(2))+(BV1(3)*BV2(3)) 107 : C 108 : C TEST SUR LA COLINEARITE DE V1 ET V2. 109 : C 110 : IF(ABS(PS).GE.(.99D0)) THEN 111 : KERRE = 1 112 : GOTO 666 113 : ELSE 114 : XNV2 = UN/SQRT(UN-(PS**2)) 115 : DO 15 I=1,3 116 : BV2(I) = XNV2*(BV2(I)-(PS*BV1(I))) 117 : 15 CONTINUE 118 : ENDIF 119 : C 120 : C CALCUL DE LA MATRICE DE PASSAGE 121 : C LOCAL = BPSS * GLOBAL 122 : C BV1 REPRESENTE L'AXE LOCAL DES Z 123 : C BV2 REPRESENTE L'AXE LOCAL DES Y 124 : C 125 : DO 20 I=1,3 126 : BPSS(3,I) =BV1(I) 127 : BPSS(2,I) =BV2(I) 128 : 20 CONTINUE 129 : BPSS(1,1) = (BPSS(2,2)*BPSS(3,3))-(BPSS(2,3)*BPSS(3,2)) 130 : BPSS(1,2) = (BPSS(2,3)*BPSS(3,1))-(BPSS(2,1)*BPSS(3,3)) 131 : BPSS(1,3) = (BPSS(2,1)*BPSS(3,2))-(BPSS(2,2)*BPSS(3,1)) 132 : C 133 : C CALCUL DES TERMES DE LA MATRICE DE RIGIDITE REPERE LOCAL 134 : C 135 : TESPI = TETA/XPI 135 : CALL ZERO(REL,12,12) 136 : C 137 : C ON TESTE SUR L'IMPORTANCE DE L'ANGLE D'OUVERTURE DE LA FISSURE. SI 138 : C CE DERNIER EST INFERIEUR A UN DEGRE, ON ANNULE TOUS LES COEFFI- 139 : C CIENTS DE COUPLAGE ET ON PENALISE TOUS LES DEGRES DE LIBERTE. DANS 140 : C LE CAS CONTRAIRE, SEULS SONT PERMIS LA ROTATION AUTOUR DE L'AXE 141 : C LOCAL DES X ET L'ALLONGEMENT SUIVANT L'AXE LOCAL DES Z ,CES DEUX 142 : C DEGRES DE LIBERTE ETANT COUPLES. 143 : C 144 : IF (TETA1.LE.(0.5D0)) THEN 145 : REL(3,3) = DEUX * YOU * XPI * EPAI * RAYMO * D6 146 : REL(4,4) = YOU * XPI * EPAI * (RAYMO**2) * D6 147 : ELSE 148 : C 148 : CALL TUFIFP(TESPI,AXX,FP,FM,FMP,FOP,FOM) 149 : DENOM = FP*FM - FMP*FMP 150 : COEF = (YOU * XPI * EPAI)/(DENOM * (TETA**2)) 151 : REL(3,3) = COEF * DEUX * FM 152 : REL(3,4) = -COEF * RAYMO * FMP 153 : REL(4,3) = REL(3,4) 154 : REL(4,4) = COEF * (RAYMO**2) * FP * UNS2 155 : ENDIF 156 : C 157 : C INITIALISATION DES TERMES DE PENALISATION 158 : C 159 : REL(1,1) = DEUX * CISA * XPI * EPAI * RAYMO * D6 160 : REL(2,2) = REL(1,1) 161 : REL(5,5) = YOU * XPI * EPAI * (RAYMO**2) * D6 162 : REL(6,6) = DEUX * CISA * XPI * EPAI * (RAYMO**2) * D6 163 : C 164 : C CALCUL DE LA MATRICE (12,12) A PARTIR DE LA MATRICE (6,6) 165 : C 166 : DO 30 I=1,6 167 : DO 301 J=1,6 168 : REL(I+6,J+6) = REL(I,J) 169 : REL(I+6,J ) = -REL(I,J) 170 : REL(I ,J+6) = -REL(I,J) 171 : 301 CONTINUE 172 : 30 CONTINUE 172 : CALL POUROT(REG,12,BPSS,REL) 173 : 666 CONTINUE 174 : RETURN 175 : END 176 : 177 : 178 : 179 : 180 :
© Cast3M 2003 - All rights reserved.
Disclaimer