Esope source files

   1 : C VISPL1    SOURCE    CHAT      05/01/13    04:07:26     5004           
   2 : C VISPL1    SOURCE    AM        00/12/21    21:26:00     4062           
   3 :       SUBROUTINE VISPL1(SIGP,DSIGT,NSTRS,SIGF,SIGRV,DSTRN,
   4 :      &         DEPS1,DEPS2,BETJEF,VISCO,NECH0,NECH1)
   5 : C
   6 :       IMPLICIT INTEGER(I-N)
   7 :       IMPLICIT REAL*8(A-H,O-Z)
   8 :       DIMENSION SIGP(4),DSIGT(4),DSTRN(4),D(4,4)
   9 :       DIMENSION DSIG(4),DSIG2(4),DSIG3(4),DP(4,4)
  10 :       DIMENSION SIGF(4),SIGRV(4),DF(4)
  11 : C
  12 : *
  13 :       SEGMENT BETJEF
  14 :        REAL*8 AA,BETA,RB,ALFA,YOUN,XNU,GFC,GFT,CAR,ETA,TDEF,
  15 :      &          TCON,DPS1,DPS2,TETA,PDT,TP0
  16 :        INTEGER ICT,ICC,IMOD,IVIS,ITER,
  17 :      &                ISIM,IBB1,IGAU1,IZON
  18 :       ENDSEGMENT
  19 :       SEGMENT VISCO
  20 :        REAL*8 DPSITV1,DPSITV2,SIGV1,SIGV2,ENDV
  21 :       ENDSEGMENT
  22 :       SEGMENT NECH0
  23 :        REAL*8 DT,DC,ALFG,S0,ENDO
  24 :       ENDSEGMENT
  25 :       SEGMENT NECH1
  26 :        REAL*8 ENDL
  27 :       ENDSEGMENT
  28 : C
  29 : C
  30 : C ********************************************************************
  31 : C **              COMPORTEMENT VISCOPLASTIQUE_ENDOMMAGEABLE         **
  32 : C **                    DE TYPE DUVAUT LIONS                        **
  33 : C ********************************************************************
  34 : C
  35 :       IPAS=0
  36 :       DAM = 0.D0
  36 :       CALL ZERO(DSIG2,4,1)
  36 :       CALL ZERO(DSIG3,4,1)
  37 :       DPSTV1 = DPSITV1
  38 :       DPSTV2 = DPSITV2
  39 :       ENDV0 = ENDV
  40 :       ENDO0 = ENDO
  41 : C
  42 : C ---------------- calcul des termes viscoplastiques -----------------
  43 : C
  44 :       GAMA = 0.5
  45 :       A = (ETA)/(ETA+GAMA*PDT)
  46 :       B = (GAMA*PDT)/ETA
  47 :       C = PDT
  48 :       G=EXP(-PDT/ETA)
  49 :       IF (C.EQ.0.D0) THEN
  50 :       WRITE(*,*)'l increment de temps est de longueur nulle'
  51 :       stop
  52 :       ENDIF
  53 : C
  54 : C -------------------- creation matrice de young -------------------
  55 : C 
  56 :       IF (IMOD.EQ.1.OR.IMOD.EQ.3) THEN
  56 :         CALL CREMAT(D,YOUN,XNU,3,-2) 
  57 :       ENDIF 
  58 :       IF (IMOD.EQ.2.OR.IMOD.EQ.4) THEN
  58 :         CALL CREMAT(D,YOUN,XNU,4,-1) 
  59 :       ENDIF
  60 : C
  61 : C~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  62 : C        MISE A JOUR DES VARIABLES D'ECROUISSAGES
  63 : C~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  64 : C
  65 :            DPSTV1=DPSTV1*G+(1.D0-G)*DEPS1
  65 :            CALL UNICO1(DPSTV1,PAEC,1,SIGV1,BETJEF,NECH0,NECH1) 
  66 :            DPSTV2=DPSTV2*G+(1.D0-G)*DEPS2
  66 :            CALL UNICO1(DPSTV2,PAEC,2,SIGV2,BETJEF,NECH0,NECH1)
  67 : C
  68 : C~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  69 : C       MISE A JOUR DE LA VARIABLE D'ENDOMMAGEMENT VISQUEUSE
  70 : C~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  71 : C
  71 :         CALL DAMAG1(DPSTV1,DPSTV2,ENDC,
  72 :      &          ENDT,IFIS,IPLA,3,BETJEF,NECH0,NECH1)
  73 :        ENDV = 1.D0-((1.D0-ENDC)*(1.D0-ENDT))
  74 : C
  75 : C ------------ mise a jour contrainte ----------------------
  76 : C ------ methode explicite
  77 : CC      DD = (ENDV - ENDV0) / PDT
  78 : CC      FF = 1. - ENDV0
  79 : CC      FK = -((1./ETA)+(DD/FF))
  80 : CC      DO 15 I=1,NSTRS
  81 : CC      DF(I) = FK*SIGRV(I)+FF*(DSIGT(I)/PDT)+(SIGP(I)/ETA)
  82 : CC   15 CONTINUE
  83 : C ------------ resolution numerique de l equa diff ---------------
  83 : CC      CALL ZERO(SIGF,4,1)
  84 : CC      DO 36 I=1,NSTRS
  85 : CC      SIGF(I)= SIGRV(I) + PDT * DF(I)         
  86 : CC   36 CONTINUE
  87 : C
  88 :       DD = (ENDV - ENDV0)/(1. - ENDV) 
  89 :       FK = (PDT/ETA)
  90 :       CC = 1./(1 + FK + DD)
  91 : C
  92 : C ------------ resolution numerique de l equa diff ---------------
  92 :       CALL ZERO(SIGF,4,1)
  93 :       DO 37 I=1,NSTRS
  94 :       SIGF(I)= CC *(SIGRV(I) + (1.-ENDV)*DSIGT(I) + FK*SIGP(I))         
  95 :    37 CONTINUE
  96 : C
  97 : C
  98 : C -------  taux de deformation en norme --------
  99 : C
 100 :       TDEF=0.D0
 101 :       D1=DSTRN(1)
 102 :       D2=DSTRN(2)
 103 :       D3=DSTRN(3)
 104 :       D4=DSTRN(4)
 105 :       TDEF=SQRT(D1*D1+D2*D2+D3*D3+D4*D4)/PDT
 106 : C
 107 : C -------  taux de contrainte en norme --------
 108 : C
 109 :       TCON=0.D0
 110 :       DO 62 I=1,NSTRS
 111 :       DSIG(I) = (SIGF(I) - SIGRV(I))
 112 :    62 CONTINUE
 113 :       S1=DSIG(1)
 114 :       S2=DSIG(2)
 115 :       S3=DSIG(3)
 116 :       S4=DSIG(4)
 117 :       TCON=SQRT(S1*S1+S2*S2+S3*S3+S4*S4)/PDT
 118 : C
 118 :       CALL INDIC1(DPSTV1,DPSTV2,IFIS,IPLA,3,BETJEF,NECH0,NECH1)
 119 : C
 120 :       DPSITV1 = DPSTV1
 121 :       DPSITV2 = DPSTV2
 122 :       RETURN
 123 :       END
 124 : 
 125 :  
 126 :  
 127 :