Télécharger jacobp.eso

Retour à la liste

Numérotation des lignes :

  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.  
  120.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales