C GRAIN     SOURCE    CHAT      05/01/13    00:20:57     5004
C GRAIN     SOURCE    LLB
      SUBROUTINE GRAIN(TAU,EVP1,EVP2,SIG1,SIG2,VAR1,VAR2,XMAT,
     &     NSTRS,NVARI,KERRE)
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C       CALCUL DE LA TAILLE DE GRAIN A L'AIDE D'UN SHEMA DE TYPE
C    RUNGE KUTTA 2.1. LE MODELE UTILISE EST CELUI DE BESSON.
C ________________________________________________________________________
C
C   EN ENTREE :
C
C    DE1     : LE TAUX DE DEFAUTS AU TEMPS i
C    GR1     : LA TAILLE DE GRAIN AU TEMPS i
C    TAU     : LE PAS DE TEMPS
C    EVP1    : LA VITESSE DE DEFORMATION VISCOPLASTIQUE AU TEMPS i
C    EVP2    : LA VITESSE DE DEFORAMTIOM VISCOPLATIQUE AU TEMPS i+TAU/2
C    SIG1    : LA CONTRAINTE AU TEMPS i
C    SIG2    : LA CONTRAINTE AU TEMPS i+TAU
C    XMAT    : LE TABLEAU DES CARACTERISTIQUES MATERIELLES AU TEMPS i+TAU
C    NSTRS   : LA DIMENSION DES VECTEURS EVP1...SIG1...
C
C   EN SORTIE :
C
C    DE2     : LE TAUX DE DEFAUT AU TEMPS i+TAU
C    GR2     : LA TAILLE DE GRAIN AU TEMPS i+TAU
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
-INC PPARAM
-INC CCOPTIO
      DIMENSION EVP1(*),EVP2(*),SIG1(*),VAR1(*),SIG2(*),VAR2(*)
      DIMENSION XMAT(*),SIGS2(6)
C
C--------------RAPPEL DES DONNES MATERIAU CALCULEE EN DEBUT DE PAS DE TEMPS
C--------------ATTENTION, DANS L'ETAT ACTUEL, CE SONT LES SEULES DONNEES
C--------------DISPONIBLES. DANS L'AVENIR, IL FAUDRA REACTUALISER LES K0, K1
C--------------K2, K3 EN COURS DE PAS DE TEMPS.
C
      GR=VAR1(NVARI-2)
      DE=VAR1(NVARI-1)
      G0 = XMAT(20)
      D0 = XMAT(21)
      RM = XMAT(22)
      RK0 = XMAT(23)
      RK1 = XMAT(24)
      RK2 = XMAT(25)
      RK3 = XMAT(26)
C
C--------------ADIMENSIONNALISATION DES VARIABLES INTERNES
C
      DE1=DE/D0
      GR1=GR/G0
C
C--------------CALCUL DE L'ENERGIE DE DEFORMATION AU TEMPS i
C
      ENERI=PROCON(EVP1,SIG1,NSTRS)
C
C--------------CALCUL DE L'ENERGIE DE DEFORMATION AU TEMPS i+TAU/2
C
      DO 10 I=1,NSTRS
       SIGS2(I)= 0.5*(SIG1(I)+SIG2(I))
10    CONTINUE
      ENERIS2=PROCON(EVP2,SIGS2,NSTRS)
C
C--------------CALCUL DE D AU TEMPS i+TAU/2
C
      D1S2=DE1+((TAU/2.D0)*(((RK1/D0)*ENERI)-(RK2*DE1)))
C
C--------------CALCUL G AU TEMPS i+TAU/2
C
      G1S2=GR1+(0.5D0*TAU*(((RK0/((G0*GR1)**RM))+(RK3*D0*DE1))/G0))
C
C--------------CALCUL DE G AU TEMPS i+1
C
      GR2=GR1+TAU*(((RK0/((G0*G1S2)**RM))+RK3*D0*D1S2)/G0)
C
C--------------CALCUL DE D AU TEMPS i+1
C
      DE2=DE1+(TAU*(((RK1/D0)*ENERIS2)-(RK2*D1S2)))
C
C   REMPLISSAGE DU TABLEAU DE VARIABLES INTERNES
C
      VAR2(NVARI-2)=GR2*G0
      VAR2(NVARI-1)=DE2*D0
C
      RETURN
      END