Télécharger gdp.eso

Retour à la liste

Numérotation des lignes :

gdp
  1. C GDP SOURCE CHAT 05/01/13 00:16:43 5004
  2. SUBROUTINE GDP(TENS,DTENS,RDP,ADP,GAMDP)
  3. IMPLICIT INTEGER(I-N)
  4. IMPLICIT REAL*8(A-H,O-Z)
  5. DIMENSION TENS(3),DTENS(3)
  6. C
  7. C----------------------------------------------
  8. C CETTE ROUTINE PERMET DE CALCULER LE
  9. C GAMMA POUR ARRIVER SUR LA LIMITE DU
  10. C CRITERE DRUCKER PRAGER
  11. C----------------------------------------------
  12. C
  13. GAMDP=0.
  14. A=-ADP*(DTENS(1)+DTENS(2))*ADP*(DTENS(1)+DTENS(2))
  15. 1 +DTENS(1)*DTENS(1)+DTENS(2)*DTENS(2)+3.D0*DTENS(3)*DTENS(3)
  16. 1 -DTENS(1)*DTENS(2)
  17. B=2.D0*ADP*(DTENS(1)+DTENS(2))*(RDP-ADP*(TENS(1)+TENS(2)))
  18. 1 +2.D0*DTENS(1)*TENS(1)+2.D0*DTENS(2)*TENS(2)+
  19. 1 6.D0*DTENS(3)*TENS(3)-DTENS(1)*TENS(2)-DTENS(2)*TENS(1)
  20. C=RDP-ADP*(TENS(1)+TENS(2))
  21. C=C*C
  22. C=-C+TENS(1)*TENS(1)+TENS(2)*TENS(2)+3.D0*TENS(3)*TENS(3)
  23. 1 -TENS(2)*TENS(1)
  24. X1=0.D0
  25. X2=0.D0
  26. IF(A.EQ.0.D0) GOTO 1000
  27. B=B/A
  28. C=C/A
  29. A=1.D0
  30. DIS=B*B-4.D0*A*C
  31. ADIS=ABS(DIS)
  32. IF(DIS.GT.0.D0) THEN
  33. DIS=SQRT(ADIS)
  34. X2=(-B+DIS)/2.D0
  35. X1=(-B-DIS)/2.D0
  36. GOTO 2000
  37. ENDIF
  38. IF(ADIS.LE.1.E-10) THEN
  39. X1=-B/2.D0
  40. X2=X1
  41. GOTO 2000
  42. ENDIF
  43. GAMDP=0.D0
  44. RETURN
  45. 1000 IF(B.EQ.0.D0) THEN
  46. X1=0.D0
  47. X2=0.D0
  48. GOTO 2000
  49. ENDIF
  50. X1=-C/B
  51. X2=X1
  52. 2000 CONTINUE
  53. S1=RDP-ADP*(TENS(1)+TENS(2)+X1*(DTENS(1)+DTENS(2)))
  54. S2=RDP-ADP*(TENS(1)+TENS(2)+X2*(DTENS(1)+DTENS(2)))
  55. IF(RDP.NE.0.D0) THEN
  56. S1=S1/RDP
  57. S2=S2/RDP
  58. ENDIF
  59. IF(S1.GE.-1.E-10.AND.S2.GE.-1.E-10) THEN
  60. GAMDP=MAX(X1,X2)
  61. GOTO 3000
  62. ENDIF
  63. IF(S1.GE.-1.E-10.AND.S2.LT.-1.E-10) THEN
  64. GAMDP=X1
  65. GOTO 3000
  66. ENDIF
  67. IF(S1.LT.-1.E-10.AND.S2.GE.-1.E-10) THEN
  68. GAMDP=X2
  69. GOTO 3000
  70. ENDIF
  71. 3000 IF(GAMDP.LT.0.D0) GAMDP=0.D0
  72. IF(GAMDP.GT.1.D0) GAMDP=1.D0
  73. RETURN
  74. END
  75.  
  76.  

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