C CQ2KSG SOURCE CHAT 05/01/12 22:26:43 5004 SUBROUTINE CQ2KSG(COO,EP,DIM3,IFOU,AN,NBPGAU,XSTRS,RAYON, . RAY,TGAUS,CGAUS,ALF11,ALF12,ALF21,ALF22,A4,S3,AM,A1,A,LRE,REL $ ) *----------------------------------------------------------------------- C CALCUL DE LA MATRICE KSIGMA POUR UNE COQUE A 2 NOEUDS C FEVRIER 86 ; X. DE MAZANCOURT C PROGRAMME REPRIS A 90% SUR LE COQFLA DE INCA *----------------------------------------------------------------------- C ENTREES: C COO:COORDONNEES DE L'ELEMENT C EP:EPAISSEUR DE LA COQUE C IFOU = VALEUR DE IFOUR DE CCOPTIO C AN:NUMERO DE L'HARMONIQUE DE FOURIER C XSTRS:CONTRAINTES C NBPGAU:NOMBRE DE POINTS DE GAUSS C TGAUS:POSITION DES POINTS DE GAUSS C CGAUS:POIDS DES POINTS DE GAUSS C RAYON,RAY,ALF:TABLEAUX DE TRAVAIL C LRE = TAILLE DE LA MATRICE KSIGMA C SORTIES: C REL:MATRICE KSIGMA *----------------------------------------------------------------------- IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) C Include contenant quelques constantes dont XPI : -INC CCREEL PARAMETER(UN=1.D0,UNDEMI=.5D0,UNQUAR=.25D0) DIMENSION A4(8,*),S3(8,*),AM(8,*),A1(8,*) DIMENSION COO(3,*),A(8,*),TGAUS(*),CGAUS(*) DIMENSION ALF11(*),ALF12(*),ALF21(*),ALF22(*),RAYON(*),REL(LRE,*) DIMENSION RAY(*),XSTRS(*),II(6) DATA II/1,2,4,5,6,8/ C --------------------------------PRELIMINAIRES CALL ZERO(REL,LRE,LRE) RAYON(1)=COO(1,1) RAYON(3)=COO(1,2) RAYON(2)=(RAYON(1)+RAYON(3))*UNDEMI DR=RAYON(3)-RAYON(1) DZ=COO(2,2)-COO(2,1) DS=SQRT(DR*DR+DZ*DZ) COSP=DR/DS SINP=DZ/DS XS2=UNDEMI*DS XS4=UNQUAR*DS C -------------------------------- CALL RESO1K(A4,DS,COSP,SINP,S3,A) CALL ZERO(ALF12,5,1) CALL ZERO(ALF22,7,1) C --------------------------------BOUCLE SUR LES POINTS DE GAUSS DO 9 IN=1,NBPGAU C ID1=1+(IN-1)*6 C ID2=2+(IN-1)*6 RAY(IN)=RAYON(2)+XS2*COSP*TGAUS(IN) X=TGAUS(IN)*XS2 X2=X*X X3=X2*X X4=X3*X C ANS=EP*XSTRS(ID1) C ANT=EP*XSTRS(ID2) C C ON PEUT FAIRE VARIER LES CONTRAINTES AVEC LE POINT DE GAUSS C MAIS ON A CHOISI DE PRENDRE CELLES DU 2} POINT DE GAUSS... C IF (IFOU.EQ.1) THEN ANS=EP*XSTRS(7) ANT=EP*XSTRS(8) ELSE ANS=EP*XSTRS(5) ANT=EP*XSTRS(6) ENDIF ALF11(1)=UN ALF11(2)=X ALF11(3)=X2 ALF11(4)=X3 ALF11(5)=X4 DO 12 NC=1,5 ALF21(NC)=ALF11(NC) 12 CONTINUE ALF21(6)=X4*X ALF21(7)=X4*X2 RRRR=RAY(IN) IF(IFOU.EQ.-3) RRRR=UN * * LE CAS DES DEF. PLANES GENERALISEES SERA A TESTER |||| * IF(IFOU.EQ.-1.OR.IFOU.EQ.-2) THEN ANS=ANS*XS4*CGAUS(IN) ANT=0.D0 ELSE ANS=ANS*RRRR*XS4*CGAUS(IN) ANT=(ANT*XS4*CGAUS(IN))/RRRR ENDIF IF (IFOU.EQ.0.OR.(IFOU.EQ.1.AND.AN.EQ.0)) THEN ANS=ANS*2*XPI ANT=ANT*2*XPI ELSEIF (IFOU.EQ.1.AND.AN.NE.0) THEN ANS=ANS*XPI ANT=ANT*XPI ENDIF ANS=ANS*DIM3 ANT=ANT*DIM3 CALL AMUL1K(ALF11,5,ANS,ALF12) CALL AMUL1K(ALF21,7,ANT,ALF22) C --------------------------------------- 9 CONTINUE C C TRANSFERT DE ALF12 (RESP.ALF22) DANS ALF11 (RESP.ALF21) C CALL SHIF1K(ALF12,ALF11,5) CALL SHIF1K(ALF22,ALF21,7) C C CALCUL DE A1,TRANSPOSEE DE A4 C DO 15 NC=1,8 DO 13 NCC=1,8 A1(NCC,NC)=A4(NC,NCC) 13 CONTINUE 15 CONTINUE CALL ZERO(A,8,8) PROD1=COSP*COSP PROD2=SINP*SINP PROD3=SINP*COSP A(1,1)=PROD1*ALF21(1) A(1,2)=PROD1*ALF21(2) A(1,5)=PROD3*ALF21(1) A(1,6)=PROD3*ALF21(2) A(1,7)=PROD3*ALF21(3) A(1,8)=PROD3*ALF21(4) A(2,1)=A(1,2) A(2,2)=ALF11(1)+PROD1*ALF21(3) A(2,5)=PROD3*ALF21(2) A(2,6)=PROD3*ALF21(3) A(2,7)=PROD3*ALF21(4) A(2,8)=PROD3*ALF21(5) A(3,3)=ALF21(1) A(3,4)=ALF21(2) A(4,3)=A(3,4) A(4,4)=ALF11(1)+ALF21(3) A(5,1)=A(1,5) A(5,2)=A(2,5) A(5,5)=PROD2*ALF21(1) A(5,6)=PROD2*ALF21(2) A(5,7)=PROD2*ALF21(3) A(5,8)=PROD2*ALF21(4) A(6,1)=A(1,6) A(6,2)=A(2,6) A(6,5)=A(5,6) A(6,6)=ALF11(1)+A(5,7) A(6,7)=ALF11(2)*2.D0+A(5,8) A(6,8)=ALF11(3)*3.D0+PROD2*ALF21(5) A(7,1)=A(1,7) A(7,2)=A(2,7) A(7,5)=A(5,7) A(7,6)=A(6,7) A(7,7)=ALF11(3)*4.D0+PROD2*ALF21(5) A(7,8)=ALF11(4)*6.D0+PROD2*ALF21(6) A(8,1)=A(1,8) A(8,2)=A(2,8) A(8,5)=A(5,8) A(8,6)=A(6,8) A(8,7)=A(7,8) A(8,8)=ALF11(5)*9.D0+PROD2*ALF21(7) 20 CONTINUE CALL MULMAT (S3,A,A4,8,8,8) CALL MULMAT (A,A1,S3,8,8,8) CALL ZERO(AM,8,8) AM(1,1)=ALF21(1) AM(3,3)=ALF21(1) AM(5,5)=ALF21(1) AM(1,2)=ALF21(2) AM(2,1)=ALF21(2) AM(3,4)=ALF21(2) AM(4,3)=ALF21(2) AM(5,6)=ALF21(2) AM(6,5)=ALF21(2) AM(2,2)=ALF21(3) AM(4,4)=ALF21(3) AM(6,6)=ALF21(3) AM(5,7)=ALF21(3) AM(7,5)=ALF21(3) AM(5,8)=ALF21(4) AM(8,5)=ALF21(4) AM(6,7)=ALF21(4) AM(7,6)=ALF21(4) AM(7,7)=ALF21(5) AM(6,8)=ALF21(5) AM(8,6)=ALF21(5) AM(7,8)=ALF21(6) AM(8,7)=ALF21(6) AM(8,8)=ALF21(7) CALL MULMAT (S3,AM,A4,8,8,8) CALL MULMAT (AM,A1,S3,8,8,8) AN2=AN*AN DO 30 NC=1,8 DO 31 NCC=1,8 A(NCC,NC)=A(NCC,NC)+AN2*AM(NCC,NC) 31 CONTINUE 30 CONTINUE PROD1=2.D0*COSP PROD2=2.D0*SINP CALL ZERO(AM,8,8) AM(1,3)=PROD1*ALF21(1) AM(1,4)=PROD1*ALF21(2) AM(2,3)=AM(1,4) AM(2,4)=PROD1*ALF21(3) AM(3,1)=AM(1,3) AM(3,2)=AM(2,3) AM(3,5)=PROD2*ALF21(1) AM(3,6)=PROD2*ALF21(2) AM(3,7)=PROD2*ALF21(3) AM(3,8)=PROD2*ALF21(4) AM(4,1)=AM(1,4) AM(4,2)=AM(2,4) AM(4,5)=AM(3,6) AM(4,6)=AM(3,7) AM(4,7)=AM(3,8) AM(7,4)=AM(4,7) AM(4,8)=PROD2*ALF21(5) AM(5,3)=AM(3,5) AM(5,4)=AM(4,5) AM(6,3)=AM(3,6) AM(6,4)=AM(4,6) AM(7,3)=AM(3,7) AM(8,3)=AM(3,8) AM(8,4)=AM(4,8) 19 CONTINUE CALL MULMAT (S3,AM,A4,8,8,8) CALL MULMAT (AM,A1,S3,8,8,8) C C CALCUL FINAL DE LA MATRICE KSIGMA C IF(IFOU.EQ.1) THEN DO 32 NC=1,8 DO 33 NCC=1,8 REL(NCC,NC) =(A(NCC,NC)+AN*AM(NCC,NC))*2.D0 33 CONTINUE 32 CONTINUE ELSE IF(IFOU.EQ.0.OR.IFOU.EQ.-1.OR.IFOU.EQ.-2) THEN DO 82 NC=1,6 IC=II(NC) DO 83 NCC=1,6 ICC=II(NCC) REL(NCC,NC) =(A(ICC,IC)+AN*AM(ICC,IC))*2.D0 83 CONTINUE 82 CONTINUE ENDIF RETURN END