Télécharger la plaquette Cast3M


fr en

Télécharger ziro.eso

Retour à la liste

Numérotation des lignes : 

ziro

  1. C ZIRO      SOURCE    CHAT      05/01/13    04:23:01     5004
  2. c******************************************************************
  3. c       sous routine utilisee dans l'ecoulement plastique siuvant *
  4. c       la loi de gurson                                          *
  5. c******************************************************************
  6.         subroutine ziro(dphi,wrkgur)
  7. c       this subroutine solve for phi1 the solution of
  8. c       phi=exp(3.|P|./2/sigbar)
  9. c
  10. c---    variables
  11.       IMPLICIT INTEGER(I-N)
  12.         IMPLICIT REAL*8 (A-H,O-Z)
  13.         segment wrkgur
  14.          real*8 sigbar, sy0,phi0,rho0,g,b,h
  15.          real*8 epn,phin,sqrtj2,rho,sig(6)
  16.          real*8 e(7),dt
  17.          real*8 conv,tol1,tol2
  18.         endsegment
  19. *        common/prop/sigbar,sy0,phi0,rho0,g,b,h
  20. *        common/state/epn,phin,sqrtj2,rho,sig(6)
  21. *       common /toler/ conv,tol1,tol2
  22. c
  23. c***    init data
  24.         phimax=1.d0-rho/rho0*(1.d0-phi0)
  25.         phi1=phimax
  26.         p=(1.d0-phin)*b*(1.d0-(1.d0-phin)*rho0/((1.d0-phi0)*rho))
  27.         f=log(phin)+1.5d0*b/sigbar*(1-phin)*(phin-phimax)/(1.d0-phimax)
  28. c
  29. c***    material in compression
  30.         if ( p .gt. 0.d0) then
  31. c
  32. c---    test if there is a root inside the intervall
  33.           if ( f .le. 0.d0) then
  34.             dphi=0.d0
  35.             return
  36.           endif
  37. c
  38. c---    test if the answer will be smaller than phi0
  39.           f=log(phi0)+1.5d0*b/sigbar*(1-phi0)
  40.      &          *(phi0-phimax)/(1.d0-phimax)
  41.           if (f .gt. 0.d0) then
  42.             dphi=phi0*0.99d0-phin
  43.             return
  44.           else
  45.             if (phimax .lt. phi0) phi1=phi0
  46.           endif
  47. c
  48. c---    iterate the tangent approximation
  49.           do 100 iter=1,50
  50.           f=log(phi1)+1.5d0*b/sigbar*(1-phi1)*(phi1-phimax)
  51.      &            /(1.d0-phimax)
  52. c
  53. c---    check for convergence
  54.           syp=sqrt(abs(1.d0+phi1*phi1-2.d0*phi1*cosh(-1.5d0*b*(1-phi1)*
  55.      &         (phi1-phimax)/(1.d0-phimax)/sigbar)))
  56. c
  57.           if ( syp .le. tol1) then
  58.             dphi=phi1-phin
  59.             return
  60.           endif
  61. c---
  62.           fprime=1.d0/phi1+1.5d0*b/sigbar*(1.d0+phimax-2.d0*phi1)/
  63.      &          (1.d0-phimax)
  64.           phi1=phi1-f/fprime
  65. c
  66.   100     continue
  67. c
  68. c***    material in tension
  69.         else
  70.           a1=1+phimax
  71. c
  72. c---    loop over the process
  73.           do 200 iter=1,90
  74. c
  75.           a2=2.d0/3.d0*sigbar/b*(1.d0-phimax)*log(phi1)+phimax
  76.           phi1=(a1-sqrt(a1*a1-4.d0*a2))*0.5d0
  77. c
  78. c---    test if the search is outside the interval
  79.           if ( phi1 .le. phin) then
  80.             dphi=0.d0
  81.             return
  82.           endif
  83. c
  84. c---    check for convergence
  85.           syp=sqrt(abs(1.d0+phi1*phi1-2.d0*phi1*cosh(-1.5d0*b*(1-phi1)*
  86.      &         (phi1-phimax)/(1.d0-phimax)/sigbar)))
  87. c
  88.           if (syp .le. tol1) then
  89.             dphi=phi1-phin
  90.             return
  91.           endif
  92. c
  93.   200     continue
  94. c
  95.         endif
  96.         dphi=phi1-phin
  97. c       the subroutine didnot converge return anyway
  98.         return
  99. c
  100.         end
  101.  
  102.  
  103.  
  104.  
  105.  

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