1 : C JACOBIX SOURCE AS218369 10/01/29 21:15:07 6623 2 : C======================================================================= 3 : C= J A C O B I X = 4 : C= ----------- = 5 : C= Fonction : = 6 : C= ---------- = 7 : C= Calcul du jacobien (DJac) et des derivees des fonctions de forme = 8 : C= par rapport aux coordonnees REELLES x,y,z (SHPTOT(2 a 4,*)) en un = 9 : C= point donne d'un element (en general, un point de Gauss) = 10 : C= = 11 : C= Parametres : (E)=Entree (S)=Sortie = 12 : C= ------------ = 13 : C= XEL(3,NBNO) (E) Coordonnees GLOBALES des noeuds de l element = 14 : C= TABA0(IDIM,NBNO) (E) Valeurs des ddl de saut aux noeuds = 15 : C= SHP(6,NBNO) (E/S) Valeurs des fonctions de forme et de leurs = 16 : C= derivees au point considere (point de Gauss) = 17 : C= IDIM (E) Dimension du probleme traite (1 a 3) = 18 : C= Cas particulier de l'element JOI3 : IDIM = 86 = 19 : C= NBNO (E) Nombre de NOEUDS de l element fini = 20 : C= DJac (S) Jacobien calcule en ce point de Gauss = 21 : C= (egal a 0 si probleme en ce point) = 22 : C= = 23 : C= Remarque : = 24 : C= ---------- = 25 : C= Lors de l'entree dans le sous-programme, SHP(2 a 4,*) contient = 26 : C= les DERIVEES des fonctions de forme par rapport aux coordonnees = 27 : C= de REFERENCE Qsi,Eta,Dzeta. = 28 : C= En sortie du sous-programme, SHP(2 a 4,*) contient les DERIVEES = 29 : C= des fonctions de FORME par rapport aux coordonnees REELLES x,y,z. = 30 : C======================================================================= 31 : 32 : SUBROUTINE JACOBIX (XEL,TABA0,TABDH,SHP,IDIM,NBNO,DJac) 33 : 34 : IMPLICIT INTEGER(I-N) 35 : IMPLICIT REAL*8 (A-H,O-Z) 36 : 37 : DIMENSION XEL(3,*),SHP(6,*) 38 : DIMENSION TABA0(IDIM,*),TABDH(*) 39 : 40 : PARAMETER (XZero=0.d0,XUn=1.d0) 41 : 42 : dJInv=XZero 43 : 44 : C =========================== 45 : C 1 - Cas de la DIMENSION 1 46 : C =========================== 47 : IF (IDIM.EQ.1) THEN 48 : write(6,*) 'On ne traite pas la dimension 1 en XFEM' 49 : return 50 : * DJac=XZero 51 : * DO i=1,NBNO 52 : * DJac=DJac+SHP(2,i)*XEL(1,i) 53 : * ENDDO 54 : * IF (DJac.NE.XZero) dJInv=XUn/DJac 55 : * DO i=1,NBNO 56 : * SHP(2,i)=SHP(2,i)*dJInv 57 : * ENDDO 58 : C =========================== 59 : C 2 - Cas de la DIMENSION 2 60 : C =========================== 61 : ELSE IF (IDIM.EQ.2) THEN 62 : dXdQsi=XZero 63 : dXdEta=XZero 64 : dYdQsi=XZero 65 : dYdEta=XZero 66 : DO i=1,NBNO 67 : * dXdQsi=dXdQsi+SHP(2,i)*(XEL(1,i)) 68 : * dXdEta=dXdEta+SHP(3,i)*(XEL(1,i)) 69 : * dYdQsi=dYdQsi+SHP(2,i)*(XEL(2,i)) 70 : * dYdEta=dYdEta+SHP(3,i)*(XEL(2,i)) 71 : dXdQsi=dXdQsi+SHP(2,i)*(XEL(1,i)+TABDH(i)*TABA0(1,i)) 72 : dXdEta=dXdEta+SHP(3,i)*(XEL(1,i)+TABDH(i)*TABA0(1,i)) 73 : dYdQsi=dYdQsi+SHP(2,i)*(XEL(2,i)+TABDH(i)*TABA0(2,i)) 74 : dYdEta=dYdEta+SHP(3,i)*(XEL(2,i)+TABDH(i)*TABA0(2,i)) 75 : ENDDO 76 : DJac=dXdQsi*dYdEta-dXdEta*dYdQsi 77 : IF (DJac.NE.XZero) dJInv=XUn/DJac 78 : dQsidX= dYdEta*dJInv 79 : dQsidY=-dXdEta*dJInv 80 : dEtadX=-dYdQsi*dJInv 81 : dEtaDY= dXdQsi*dJInv 82 : DO i=1,NBNO 83 : V1=SHP(2,i)*dQsidX+SHP(3,i)*dEtadX 84 : SHP(3,i)=SHP(2,i)*dQsidY+SHP(3,i)*dEtadY 85 : SHP(2,i)=V1 86 : ENDDO 87 : C =========================== 88 : C 3 - Cas de la DIMENSION 3 89 : C =========================== 90 : ELSE IF (IDIM.EQ.3) THEN 91 : * 14/01/10: Il faut encore tester qu'on obtient de bons résultats !!! 92 : D11=XZero 93 : D21=XZero 94 : D31=XZero 95 : D12=XZero 96 : D22=XZero 97 : D32=XZero 98 : D13=XZero 99 : D23=XZero 100 : D33=XZero 101 : DO i=1,NBNO 102 : D11=D11+SHP(2,i)*(XEL(1,i)+TABDH(i)*TABA0(1,i)) 103 : D21=D21+SHP(3,i)*(XEL(1,i)+TABDH(i)*TABA0(1,i)) 104 : D31=D31+SHP(4,i)*(XEL(1,i)+TABDH(i)*TABA0(1,i)) 105 : D12=D12+SHP(2,i)*(XEL(2,i)+TABDH(i)*TABA0(2,i)) 106 : D22=D22+SHP(3,i)*(XEL(2,i)+TABDH(i)*TABA0(2,i)) 107 : D32=D32+SHP(4,i)*(XEL(2,i)+TABDH(i)*TABA0(2,i)) 108 : D13=D13+SHP(2,i)*(XEL(3,i)+TABDH(i)*TABA0(3,i)) 109 : D23=D23+SHP(3,i)*(XEL(3,i)+TABDH(i)*TABA0(3,i)) 110 : D33=D33+SHP(4,i)*(XEL(3,i)+TABDH(i)*TABA0(3,i)) 111 : ENDDO 112 : DInv11=D22*D33-D23*D32 113 : DInv12=D32*D13-D12*D33 114 : DInv13=D12*D23-D22*D13 115 : DJac=D11*DInv11+D21*DInv12+D31*DInv13 116 : IF (DJac.NE.XZero) dJInv=XUn/DJac 117 : DInv11=DInv11*dJInv 118 : DInv12=DInv12*dJInv 119 : DInv13=DInv13*dJInv 120 : DInv21=(D23*D31-D21*D33)*dJInv 121 : DInv22=(D11*D33-D13*D31)*dJInv 122 : DInv23=(D21*D13-D11*D23)*dJInv 123 : DInv31=(D21*D32-D22*D31)*dJInv 124 : DInv32=(D12*D31-D11*D32)*dJInv 125 : DInv33=(D11*D22-D12*D21)*dJInv 126 : DO i=1,NBNO 127 : V1=DInv11*SHP(2,i)+DInv12*SHP(3,i)+DInv13*SHP(4,i) 128 : V2=DInv21*SHP(2,i)+DInv22*SHP(3,i)+DInv23*SHP(4,i) 129 : V3=DInv31*SHP(2,i)+DInv32*SHP(3,i)+DInv33*SHP(4,i) 130 : SHP(2,i)=V1 131 : SHP(3,i)=V2 132 : SHP(4,i)=V3 133 : ENDDO 134 : ENDIF 135 : 136 : RETURN 137 : END 138 : 139 : 140 : 141 : 142 :
© Cast3M 2003 - All rights reserved.
Disclaimer