1 : C JACOBP SOURCE CHAT 05/01/13 00:48:29 5004 2 : SUBROUTINE JACOBP(XEL,SHP,IDIM,NBBB,NBNO,DJAC) 3 : *----------------------------------------------------------------------- 4 : * 5 : * SOUS-PROGRAMME JACOBI ADAPTE AUX MILIEUX POREUX 6 : * 7 : * CALCULE LES DERIVES DES FONCTIONS DE FORME 8 : * DANS LA GEOMETRIE DEFORMEE A PARTIR DES DERIVEES 9 : * DANS LA GEOMETRIE REDUITE AINSI QUE LE JACOBIEN 10 : * ENTREES 11 : * XEL(3,NBBB)=COORDONNEES LOCALES DE L ELEMENT 12 : * SHP(6,NBBB)=DERIVEES PAR RAPPORT A LA GEOMETRIE DE REFERENCE 13 : * IDIM=DIMENSION 14 : * NBBB=NOMBRE DE NOEUDS 15 : * NBNO=NOMBRE DE FONCTIONS DE FORME 16 : * SORTIES 17 : * SHP(6,NBNO)=DERIVEES PAR RAPPORT A LA GEOMTRIE DE L ELEMENT 18 : * DJAC =JACOBIEN ( MIS A 0 SI NON INVERSIBLE ) 19 : * 20 : *----------------------------------------------------------------------- 21 : IMPLICIT INTEGER(I-N) 22 : IMPLICIT REAL*8 (A-H,O-Z) 23 : DIMENSION XEL(3,*),SHP(6,*) 24 : DIMENSION D(3,3),DINV(3,3),V(3) 25 : DXDQSI=0.D0 26 : DXDETA=0.D0 27 : DYDQSI=0.D0 28 : DYDETA=0.D0 29 : C 30 : C CAS MONODIMENSIONEL 31 : C 32 : IF(IDIM.NE.1) GOTO 20 33 : DO 400 I=1,NBBB 34 : DXDQSI=DXDQSI+SHP(2,I)*XEL(1,I) 35 : 400 CONTINUE 36 : DJAC=DXDQSI 37 : XXXX = DJAC 38 : IF(DJAC.NE.0.D0) XXXX=1.D0/DJAC 39 : DO 410 I=1,NBNO 40 : SHP(2,I)=SHP(2,I)*XXXX 41 : 410 CONTINUE 42 : GOTO 666 43 : C 44 : C CAS 2 DIMENSIONS 45 : C 46 : 20 CONTINUE 47 : IF (IDIM.NE.2) GOTO 30 48 : DO 100 I=1,NBBB 49 : DXDQSI=DXDQSI+SHP(2,I)*XEL(1,I) 50 : DXDETA=DXDETA+SHP(3,I)*XEL(1,I) 51 : DYDQSI=DYDQSI+SHP(2,I)*XEL(2,I) 52 : DYDETA=DYDETA+SHP(3,I)*XEL(2,I) 53 : 100 CONTINUE 54 : DJAC=DXDQSI*DYDETA-DXDETA*DYDQSI 55 : XXXX = DJAC 56 : IF(DJAC.NE.0.D0) XXXX=1.D0/DJAC 57 : DQSIDX= DYDETA*XXXX 58 : DQSIDY=-DXDETA*XXXX 59 : DETADX=-DYDQSI*XXXX 60 : DETADY= DXDQSI*XXXX 61 : DO 200 I=1,NBNO 62 : TT=SHP(2,I)*DQSIDX+SHP(3,I)*DETADX 63 : SHP(3,I)=SHP(2,I)*DQSIDY+SHP(3,I)*DETADY 64 : SHP(2,I)=TT 65 : 200 CONTINUE 66 : GOTO 666 67 : C 68 : C CAS TRIDIMENSIONEL 69 : C 70 : 30 IF (IDIM.NE.3) GOTO 666 71 : DO 310 I=1,3 72 : DO 310 J=1,3 73 : D(I,J)=0.D0 74 : 310 CONTINUE 75 : DO 300 I=1,NBBB 76 : C 77 : D(1,1)=D(1,1)+SHP(2,I)*XEL(1,I) 78 : D(2,1)=D(2,1)+SHP(3,I)*XEL(1,I) 79 : D(3,1)=D(3,1)+SHP(4,I)*XEL(1,I) 80 : C 81 : D(1,2)=D(1,2)+SHP(2,I)*XEL(2,I) 82 : D(2,2)=D(2,2)+SHP(3,I)*XEL(2,I) 83 : D(3,2)=D(3,2)+SHP(4,I)*XEL(2,I) 84 : C 85 : D(1,3)=D(1,3)+SHP(2,I)*XEL(3,I) 86 : D(2,3)=D(2,3)+SHP(3,I)*XEL(3,I) 87 : D(3,3)=D(3,3)+SHP(4,I)*XEL(3,I) 88 : 300 CONTINUE 89 : C INVERSION DE LA MATRICE D(3,3) 90 : DINV(1,1)= D(2,2)*D(3,3)-D(2,3)*D(3,2) 91 : DINV(2,2)= D(1,1)*D(3,3)-D(1,3)*D(3,1) 92 : DINV(3,3)= D(1,1)*D(2,2)-D(1,2)*D(2,1) 93 : DINV(1,2)=-D(1,2)*D(3,3)+D(3,2)*D(1,3) 94 : DINV(2,1)=-D(2,1)*D(3,3)+D(2,3)*D(3,1) 95 : DINV(1,3)= D(1,2)*D(2,3)-D(2,2)*D(1,3) 96 : DINV(3,1)= D(2,1)*D(3,2)-D(2,2)*D(3,1) 97 : DINV(2,3)=-D(1,1)*D(2,3)+D(2,1)*D(1,3) 98 : DINV(3,2)=-D(1,1)*D(3,2)+D(1,2)*D(3,1) 99 : DJAC=D(1,1)*DINV(1,1)+D(2,1)*DINV(1,2)+D(3,1)*DINV(1,3) 100 : XXXX = DJAC 101 : IF(DJAC.NE.0.D0) XXXX=1.D0/DJAC 102 : DO 350 IA=1,3 103 : DO 350 IB=1,3 104 : DINV(IA,IB)=DINV(IA,IB)*XXXX 105 : 350 CONTINUE 106 : DO 360 IA=1,NBNO 107 : DO 370 IB=1,3 108 : V(IB)=0.D0 109 : DO 370 IC=1,3 110 : V(IB)=V(IB)+DINV(IB,IC)*SHP(IC+1,IA) 111 : 370 CONTINUE 112 : SHP(2,IA)=V(1) 113 : SHP(3,IA)=V(2) 114 : SHP(4,IA)=V(3) 115 : 360 CONTINUE 116 : 666 CONTINUE 117 : RETURN 118 : END 119 :
© Cast3M 2003 - All rights reserved.
Disclaimer