1 : C BCOQ2 SOURCE CHAT 06/03/29 21:15:33 5360 2 : SUBROUTINE BCOQ2(B,NSTRS,DJAC,IGAU,IFOU,XEL,NN,T,P,EXCEN,DIM3, 3 : . IARR,XDPGE,YDPGE) 4 : C================================================================== 5 : C CALCUL DE LA MATRICE B DES COQUES @ 2 NOEUDS 6 : C ROUTINE FORTRAN PUR 7 : C GELE JANVIER 87 8 : C================================================================== 9 : C ENTREES 10 : C IGAU=NUMERO DU POINT DE GAUSS 11 : C IFOU=IFOUR DE CCOPTIO 12 : C XEL=COORDONNEES LOCALE DE L'ELEMENT 13 : C NN=NUMERO DU MODE DE FOURIER 14 : C T(IGAU)=POSITION DU POINT DE GAUSS 15 : C P(IGAU)=POIDS DU POINT DE GAUSS 16 : C EXCEN = EXCENTREMENT 17 : C DIM3 = EPAISSEUR DANS L'AUTRE DIMENSION 18 : C XDPGE,YDPGE : COORDONNEE DU POINT AUTOUR DUQUEL 19 : C FAIT LE MOUVEMENT EN DEFO PLAN GENE 20 : C SORTIE 21 : C B(6,8)=MATRICE B AU POINT DE GAUSSPOUR IFOU GT 0 22 : C B(4,6)=MATRICE B AU POINT DE GAUSSPOUR IFOU LE 0 23 : C B(4,9)=MATRICE B AU POINT DE GAUSSPOUR IFOU EQ -3 24 : C DJAC=JACOBIEN AU POINT DE GAUSS=POIGAU*LONG/2 (*R(IGAU), SI 25 : C IFOU EST SUPERIEUR OU EGAL A ZERO) 26 : C IARR=0 SI OK 1 SI LONGUEUR ELEMENT NULLE 27 : C 2 SI R / D INFERIEUR A 10-3 28 : C================================================================== 29 : IMPLICIT INTEGER(I-N) 30 : IMPLICIT REAL*8(A-H,O-Z) 31 : C Include contenant quelques constantes dont XPI, XZERO : 32 : -INC CCREEL 33 : PARAMETER(UNDE=.5D0,UN=1.D0,DEUX=2.D0,TRS=3.D0) 34 : PARAMETER(QUTR=4.D0,SIX=6.D0,DOUZ=12.D0) 35 : DIMENSION B(NSTRS,*),T(*),XEL(3,*),P(*) 36 : C 37 : C ---------------------------------INITIALISATION 38 : IARR=0 39 : DJAC=XZERO 40 : IF(IFOU.GT.0) THEN 41 : LRE=8 41 : CALL ZERO(B,6,8) 42 : ELSE IF(IFOU.LE.0) THEN 43 : LRE=6 43 : CALL ZERO(B,4,6) 44 : IF(IFOU.EQ.-3) THEN 45 : LRE=9 45 : CALL ZERO(B,4,9) 46 : ENDIF 47 : ENDIF 48 : C 49 : D=SQRT((XEL(1,2)-XEL(1,1))**2+(XEL(2,2)-XEL(2,1))**2) 50 : IF(D.EQ.0) THEN 51 : IARR=1 52 : GOTO 4 53 : ENDIF 54 : DD=UN/D 55 : RO=(XEL(1,1)+XEL(1,2))*UNDE 56 : R1=(XEL(2,1)+XEL(2,2))*UNDE 57 : SP=(XEL(1,2)-XEL(1,1))/D 58 : CP=(XEL(2,2)-XEL(2,1))/D 59 : SP2=SP*SP 60 : CP2=CP*CP 61 : SPCP=SP*CP 62 : C X FONCTION FORME NOEUD 2 RRRR RAYON 63 : C Y FONCTION FORME NOEUD 1 D LONGUEUR DD INVERSE LONGUEUR 64 : X=UNDE+UNDE*T(IGAU) 65 : Y=UNDE-UNDE*T(IGAU) 66 : RRRR=RO+UNDE*D*SP*T(IGAU) 67 : C --------------------------------- 68 : C 69 : C TEST D'ERREUR 70 : C 71 : IF(IFOU.GE.0) THEN 72 : IF(ABS(RRRR/D).LE.1.D-03) THEN 73 : IARR=2 74 : GOTO 4 75 : ENDIF 76 : ENDIF 77 : C 78 : C ---------------------------------CALCULS 79 : C 80 : IF(IFOU.LT.0) RRRR =UN 81 : U=X/RRRR 82 : V=Y/RRRR 83 : C 84 : C AXISYMETRIQUE DEFORMATIONS PLANES CONTRAINTES PLANES 85 : C 86 : IF(IFOU.LE.0) THEN 87 : C 88 : C EPSILON S S 89 : C 90 : B(1,1)=-SP*DD 91 : B(1,2)=-CP*DD 92 : B(1,4)=-B(1,1) 93 : B(1,5)=-B(1,2) 94 : C 95 : C KSI S S 96 : C 97 : AUX = SIX*T(IGAU)*DD*DD 98 : B(3,1)= CP*AUX 99 : B(3,2)=-SP*AUX 100 : B(3,3)=(QUTR-SIX*X)*DD 101 : B(3,4)=-B(3,1) 102 : B(3,5)=-B(3,2) 103 : B(3,6)=(DEUX-SIX*X)*DD 104 : C 105 : C MODIFICATION DEFORMATION PLANE GENERALISEE 106 : C 107 : IF (IFOU.EQ.-3) THEN 108 : RRRX=RO+UNDE*D*SP*T(IGAU) 109 : RRRY=R1+UNDE*D*CP*T(IGAU) 110 : B(2,7 )=UN 111 : B(2,8)=XDPGE-RRRX 112 : B(2,9)=RRRY-YDPGE 113 : B(4,8)=CP 114 : B(4,9)=SP 115 : ENDIF 116 : C 117 : C AXISYMETRIQUE 118 : C 119 : IF(IFOU.EQ.0) THEN 120 : C 121 : C EPSILON THETA THETA 122 : C 123 : B(2,1)= V*(SP2+CP2*Y*(UN+DEUX*X)) 124 : B(2,2)= SPCP*U*Y*T(IGAU) 125 : B(2,3)=-D*CP*X*Y*V 126 : B(2,4)= U*(SP2+CP2*X*(TRS-DEUX*X)) 127 : B(2,5)=-B(2,2) 128 : B(2,6)= D*CP*X*Y*U 129 : C 130 : C KSI THETA THETA 131 : C 132 : AUX =-SIX*U*Y*SP*DD 133 : B(4,1)= CP*AUX 134 : B(4,2)=-SP*AUX 135 : B(4,3)= SP*V*(TRS*X-UN) 136 : B(4,4)=-CP*AUX 137 : B(4,5)= SP*AUX 138 : B(4,6)=-SP*U*(TRS*X-DEUX) 139 : ENDIF 140 : C 141 : C FOURIER 142 : C 143 : ELSE IF(IFOU.GT.0) THEN 144 : AN=DBLE(NN) 145 : C 146 : C EPSILON S S 147 : C 148 : B(1,1)=-SP*DD 149 : B(1,2)=-CP*DD 150 : B(1,5)=-B(1,1) 151 : B(1,6)=-B(1,2) 152 : C 153 : C EPSILON THETA THETA 154 : C 155 : B(2,1)= V*(SP2+CP2*Y*(UN+DEUX*X)) 156 : B(2,2)= SPCP*U*Y*T(IGAU) 157 : B(2,3)=+AN*V 158 : B(2,4)=-D*CP*X*Y*V 159 : B(2,5)= U*(SP2+CP2*X*(TRS-DEUX*X)) 160 : B(2,6)=-B(2,2) 161 : B(2,7)=+AN*U 162 : B(2,8)= D*CP*X*Y*U 163 : C 164 : C EPSILON S THETA 165 : C 166 : B(3,1)=-SP*B(2,3) 167 : B(3,2)=-CP*B(2,3) 168 : B(3,3)=-DD-V*SP 169 : B(3,5)=-SP*B(2,7) 170 : B(3,6)=-CP*B(2,7) 171 : B(3,7)=+DD-U*SP 172 : C 173 : C KSI S S 174 : C 175 : AUX= SIX*T(IGAU)*DD*DD 176 : B(4,1)= CP*AUX 177 : B(4,2)=-SP*AUX 178 : B(4,4)=(QUTR-SIX*X)*DD 179 : B(4,5)=-B(4,1) 180 : B(4,6)=-B(4,2) 181 : B(4,8)=(DEUX-SIX*X)*DD 182 : C 183 : C KSI THETA THETA 184 : C 185 : AUX1=SIX*U*Y*SP*DD 186 : AUX=(UN+DEUX*X)*(AN*V)**2+AUX1 187 : B(5,1)=-CP*AUX 188 : B(5,2)= SP*AUX 189 : B(5,3)=-B(2,3)*CP/RRRR 190 : B(5,4)= D*X*(AN*V)**2-SP*V*(UN-TRS*X) 191 : AUX=(TRS-DEUX*X)*(AN*U)**2-AUX1 192 : B(5,5)=-CP*AUX 193 : B(5,6)= SP*AUX 194 : B(5,7)=-B(2,7)*CP/RRRR 195 : B(5,8)=-D*Y*(AN*U)**2-SP*U*(TRS*X-DEUX) 196 : C 197 : C KSI S THETA 198 : C 199 : AUX1=DOUZ*U*Y*AN*DD 200 : AUX=-AN*SP*(DEUX+QUTR*X)*V*V-AUX1 201 : B(6,1)=-CP*AUX 202 : B(6,2)= SP*AUX 203 : B(6,3)= DEUX*(CP/D+SPCP*V)/RRRR 204 : B(6,4)= AN*V*((DEUX-SIX*X)-DEUX*D*X*V*SP) 205 : AUX=AN*SP*(QUTR*X-SIX)*U*U+AUX1 206 : B(6,5)=-CP*AUX 207 : B(6,6)= SP*AUX 208 : B(6,7)= DEUX*(-CP/D+SPCP*U)/RRRR 209 : B(6,8)= AN*U*((SIX*X-QUTR)+DEUX*D*X*V*SP) 210 : ENDIF 211 : IF(IFOU.EQ.0.OR.(IFOU.EQ.1.AND.NN.EQ.0)) THEN 212 : DJAC=D*UNDE*P(IGAU)*RRRR*2*XPI 213 : ELSEIF(IFOU.EQ.1.AND.NN.NE.0) THEN 214 : DJAC=D*UNDE*P(IGAU)*RRRR*XPI 215 : ELSE 216 : DJAC=D*UNDE*P(IGAU)*RRRR*DIM3 217 : ENDIF 218 : * 219 : * ON MODIFIE LA MATRICE B EN CAS D'EXCENTREMENT 220 : * 221 : NSTRS2 = NSTRS / 2 222 : DO 88 IJL=1,NSTRS2 223 : DO 881 IJC=1,LRE 224 : B(IJL,IJC)=B(IJL,IJC)+EXCEN*B(IJL+NSTRS2,IJC) 225 : 881 CONTINUE 226 : 88 CONTINUE 227 : * 228 : 4 CONTINUE 229 : RETURN 230 : END 231 : 232 : 233 : 234 : 235 : 236 : 237 : 238 :
© Cast3M 2003 - All rights reserved.
Disclaimer