1 : C MICROP SOURCE CHAT 05/01/13 01:45:08 5004 2 : SUBROUTINE MICROP (WRK0,WRK1,NSTRS,NVARI,NMATT, 3 : & ICARA,KERRE,MFR1,IFOURB) 4 : * 5 : * modele d'endommagement microplan couple a la plasticite 6 : * C. La Borderie + S. Fichant Oct. 95 7 : * routines utilisees: 8 : * micro1: plasticite nadai 9 : * IDECAL=3 DANS LE CAS ISO IDECAL=8 DANS LE CAS ANISO 10 : * jacob3: diagonalisation: 11 : * attention jacob3 modifie la matrice a diagonaliser!! 12 : * prodt et prodt2 13 : * attention prodt2 ne fonctionne qu'avec la matrice des V. P. !! 14 : * 15 : IMPLICIT INTEGER(I-N) 16 : C IMPLICIT REAL*8(A-H,O-Z) 17 : INTEGER NSTRS,NVARI,NMATT,ICARA,KERRE,MFR1,IFOURB 18 : INTEGER ISTRS,I,J 19 : REAL*8 YOUNG,XNU,EPSD0,BT,LAMB,DEUXMU,ALFA 20 : REAL*8 DOM(6),SIGAN(6),TRSIG,DEF33(3,3),EPSIPP(3),VECP(3,3) 21 : REAL*8 VECPT(3,3),DOM33(3,3),DEFRPDE(3,3),DOMRPDE(3,3) 22 : REAL*8 D1,D2,D3,DOM3(3) 23 : REAL*8 SIGPP(3),SIGPM(3),SIG33(3,3),SIG33P(3,3),SIG33M(3,3) 24 : REAL*8 S33PRD(3,3),S33MRD(3,3) 25 : REAL*8 LAMBDAP(6), LAMBDAM(6) 26 : REAL*8 DSIGT(6) 27 : LOGICAL COMP 28 : INTEGER IDECAL 29 : * 30 : -INC CCOPTIO 31 : * 32 : SEGMENT WRK0 33 : REAL*8 XMAT(NMATT) 34 : ENDSEGMENT 35 : * 36 : SEGMENT WRK1 37 : REAL*8 DDHOOK(LHOOK,LHOOK),SIG0(NSTRS),DEPST(NSTRS) 38 : REAL*8 SIGF(NSTRS),VAR0(NVARI),VARF(NVARI) 39 : REAL*8 DEFP(NSTRS),XCAR(ICARA) 40 : ENDSEGMENT 41 : SEGMENT WRKK1 42 : REAL*8 DEFELA(NSTRS) 43 : REAL*8 DDEFEL(NSTRS) 44 : ENDSEGMENT 45 : * 46 : * on recupere les variable materielles 47 : * 48 : YOUNG=XMAT(1) 49 : XNU=XMAT(2) 50 : EPSD0=XMAT(5) 51 : BT=XMAT(6) 52 : ALFA=XMAT(9) 53 : DEUXMU=YOUNG/(1.D0+XNU) 54 : LAMB=XNU*DEUXMU/(1.D0-2.D0*XNU) 55 : * print*,'dans microp' 56 : * 57 : * recuperation des variables internes d'endommagement 58 : * 59 : DO ISTRS=1,NSTRS 60 : DOM(ISTRS)=VAR0(2+ISTRS) 61 : END DO 62 : SEGINI WRKK1 63 : * 64 : * on ecoule plastiquement sur la contrainte effective 65 : * 66 : IDECAL=8 67 : * print*,'appel a micro1' 67 : CALL MICRO1(WRK0,NMATT,WRK1,NSTRS,DSIGT, 68 : 1 NVARI,IDECAL,KERRE,.false.) 69 : IF (KERRE .NE. 0) THEN 70 : print*,'on n''a pas converge dans micro1' 70 : CALL MICRO1(WRK0,NMATT,WRK1,NSTRS,DSIGT, 71 : 1 NVARI,IDECAL,KERRE,.true.) 72 : SEGSUP WRKK1 73 : RETURN 74 : ENDIF 75 : * 76 : * on ecoule en endommagement sur les deformations elastiques 77 : * 78 : * print*,'apres micro1' 79 : * 80 : * calcul de l'increment deformations elastiques DDEFEL 81 : * 1) on calcule l'increment de deformations totales avec 82 : * l'increment de contraintes elastique DSIGT 83 : * 2) on retranche l'increment de deformations plastiques DEFP 84 : * 85 : TRSIG=0.D0 86 : DO ISTRS=1,3 87 : TRSIG=TRSIG+DSIGT(ISTRS) 88 : END DO 89 : DO ISTRS=1,3 90 : DDEFEL(ISTRS)=( (1.D0+XNU)*DSIGT(ISTRS)-XNU*TRSIG)/YOUNG 91 : 1 - DEFP(ISTRS) 92 : END DO 93 : DO ISTRS=4,NSTRS 94 : DDEFEL(ISTRS)= (1.D0+XNU)*DSIGT(ISTRS)/YOUNG 95 : 1 - 0.5d0*DEFP(ISTRS) 96 : END DO 97 : * 98 : * on diagonalise l'increment de deformations elastiques: 99 : * 1) on met sous forme 3x3 avec endoca 100 : * 2) on diagonalise avec jacob3 101 : * 102 : * print*,'incr de def el dans rpg' 103 : * print*,ddefel 104 : * print*,'avant endoca ddefl' 104 : CALL ENDOCA (DDEFEL,DEF33,1) 105 : * print*,'apres endoca ddefl' 105 : CALL JACOB3 (DEF33,IDIM,EPSIPP,VECP) 106 : * print*,'apres jacob3 def33' 107 : * print*,'incr deformations principales' 108 : * print*,epsipp 109 : * 110 : * calcul des deformations elastiques 111 : * 112 : TRSIG=0.D0 113 : * print*,'-----sigf------' 114 : * print*,sigf(1),sigf(2),sigf(3) 115 : DO ISTRS=1,3 116 : TRSIG=TRSIG+SIGF(ISTRS) 117 : END DO 118 : DO ISTRS=1,3 119 : DEFELA(ISTRS)=( (1.D0+XNU)*SIGF(ISTRS)-XNU*TRSIG)/YOUNG 120 : END DO 121 : * print*,'-----defela------' 122 : * print*,defela(1),defela(2),defela(3) 123 : DO ISTRS=4,NSTRS 124 : DEFELA(ISTRS)= (1.D0+XNU)*SIGF(ISTRS)/YOUNG 125 : END DO 126 : * print*,'deformations elastiques dans rpg' 127 : * print*,defela 128 : * 129 : * on met les deformations sous forme de matrice 3x3 130 : * puis on ecrit la matrice dans le repere de depsilon:DEFRPDE 131 : * 131 : CALL ENDOCA (DEFELA,DEF33,1) 132 : * print*,'deformations elastiques dans rpg(3x3)' 133 : * print*,def33 133 : CALL PRODT (DEFRPDE,DEF33,VECP,3,3) 134 : * print*,'def elast dans rpddeps' 135 : * print*,defrpde 136 : * 137 : * on calcule l'endommagement resultant 138 : * 139 : IF (DEFRPDE(1,1) .GT. EPSD0) THEN 140 : D1=1.D0-EPSD0/DEFRPDE(1,1)*DEXP(BT*( EPSD0 - DEFRPDE(1,1))) 141 : ELSE 142 : D1=0.D0 143 : END IF 144 : IF (DEFRPDE(2,2) .GT. EPSD0) THEN 145 : D2=1.D0-EPSD0/DEFRPDE(2,2)*DEXP(BT*( EPSD0 - DEFRPDE(2,2))) 146 : ELSE 147 : D2=0.D0 148 : END IF 149 : IF (DEFRPDE(3,3) .GT. EPSD0) THEN 150 : D3=1.D0-EPSD0/DEFRPDE(3,3)*DEXP(BT*( EPSD0 - DEFRPDE(3,3))) 151 : ELSE 152 : D3=0.D0 153 : END IF 154 : * print*,defrpde(1,1),epsd0,D1 155 : * print*,defrpde(2,2),epsd0,D2 156 : * print*,defrpde(3,3),epsd0,D3 157 : * 158 : * on met l'endommagement initial dans le meme repere 159 : * 159 : CALL ENDOCA(DOM,DOM33,1) 160 : * print*,'prodt DOMRPDE' 160 : CALL PRODT (DOMRPDE,DOM33,VECP,3,3) 161 : * print*,'endommagement initial dans rpddeps' 162 : * print*,domrpde 163 : * 164 : * et on en deduit l'increment d'endommagement dans RPDE 165 : * 166 : IF(d1.gt.domrpde(1,1))then 167 : domrpde(1,1)=d1 168 : endif 169 : IF(d2.gt.domrpde(2,2))then 170 : domrpde(2,2)=d2 171 : endif 172 : IF(d3.gt.domrpde(3,3))then 173 : domrpde(3,3)=d3 174 : endif 175 : * print*,'endommagement final dans rpddeps' 176 : * print*,domrpde 177 : * 178 : * on remet D dans RPG 179 : * 180 : DO I=1,3 181 : DO J=1,3 182 : VECPT(I,J)=VECP(J,I) 183 : END DO 184 : END DO 184 : call prodt(DOM33,DOMRPDE,VECPT,3,3) 185 : * print*,'endommagement final dans rpg' 186 : * print*,dom33 187 : 188 : * 189 : * on separe les contraintes effectives en + et - dans rpsigma 190 : * 190 : CALL ENDOCA (SIGF,SIG33,1) 190 : CALL JACOB3 (SIG33,3,SIGPP,VECP) 191 : * print*,'contraintes ppales' 192 : * print*,sigpp 193 : DO I=1,3 194 : IF (SIGPP(I) .LT. 0.D0)THEN 195 : SIGPM(I)=SIGPP(I) 196 : SIGPP(I)=0.D0 197 : ELSE 198 : SIGPM(I)=0.D0 199 : END IF 200 : END DO 200 : CALL PRODT2(SIG33P,SIGPP,VECP,3) 200 : CALL PRODT2(SIG33M,SIGPM,VECP,3) 201 : * print*,'contraintes dans rpg' 202 : * print*,sig33p 203 : * print*,sig33m 204 : * 205 : * on met le tout dans le repere ppal d'endo 206 : * attention jacob3 modifie la matrice fournie 207 : * --> on passe une copie 208 : DO I=1,3 209 : DO J=1,3 210 : SIG33(I,J)=DOM33(I,J) 211 : END DO 212 : END DO 212 : CALL JACOB3(SIG33,3,DOM3,VECP) 213 : * print*,'endom dans rpd' 214 : * print*,dom3 215 : ************** 216 : * A REVOIR EN ATENDANT MIEUX ON BORNE LES VALEURS PROPRES DE D ENTRE 0. ET 1.* 217 : ************** 218 : DO I=1,3 219 : DOM3(I)=DMAX1 (DOM3(I),0.D0) 220 : DOM3(I)=DMIN1 (DOM3(I),1.D0-1.d-6) 221 : END DO 221 : CALL PRODT2 (DOM33,DOM3,VECP,3) 222 : ************** 223 : * FIN A REVOIR 224 : ************** 225 : * print*,'dom3=',dom3 226 : * print*,'dom33=',dom33 226 : CALL PRODT (S33PRD,SIG33P,VECP,3,3) 226 : CALL PRODT (S33MRD,SIG33M,VECP,3,3) 227 : 228 : * print*,'s33prd',s33prd(3,3) 229 : * print*,'s33mrd',s33mrd(3,3) 229 : * print*,'callambdap' 230 : COMP=.FALSE. 230 : CALL CLMBDA(DOM3(1),DOM3(2),DOM3(3),LAMBDAP,YOUNG,XNU,ALFA, 231 : 1 COMP) 232 : COMP=.TRUE. 232 : CALL CLMBDA(DOM3(1),DOM3(2),DOM3(3),LAMBDAM,YOUNG,XNU,ALFA, 233 : 1 COMP) 233 : * print*,'apres callambdap' 234 : SIG33P(1,1)=LAMBDAP(1)*S33PRD(1,1)+LAMBDAM(1)*S33MRD(1,1) 235 : SIG33P(1,2)=LAMBDAP(4)*S33PRD(1,2)+LAMBDAM(4)*S33MRD(1,2) 236 : SIG33P(1,3)=LAMBDAP(5)*S33PRD(1,3)+LAMBDAM(5)*S33MRD(1,3) 237 : SIG33P(2,1)=SIG33P(1,2) 238 : SIG33P(2,2)=LAMBDAP(2)*S33PRD(2,2)+LAMBDAM(2)*S33MRD(2,2) 239 : SIG33P(2,3)=LAMBDAP(6)*S33PRD(2,3)+LAMBDAM(6)*S33MRD(2,3) 240 : SIG33P(3,1)=SIG33P(1,3) 241 : SIG33P(3,2)=SIG33P(2,3) 242 : SIG33P(3,3)=LAMBDAP(3)*S33PRD(3,3)+LAMBDAM(3)*S33MRD(3,3) 243 : * print*,'contraintes dans RPD' 244 : * print*,sig33p 245 : * print*,'sig33p(1,2)',sig33p(1,2) 246 : * print*,'dsigt',dsigt(1),dsigt(2),dsigt(3) 247 : * print*,'dsigt',dsigt(4),dsigt(5),dsigt(6) 248 : * print*,'sigf',sigf(1),sigf(2),sigf(3) 249 : * print*,'sigf',sigf(4),sigf(5),sigf(6) 250 : * 251 : * on remet le tout dans le repere global 252 : * 253 : * print*,'matrice des vect proprs' 254 : * print*,vecp 255 : DO I=1,3 256 : DO J=1,3 257 : VECPT(I,J)=VECP(J,I) 258 : END DO 259 : END DO 259 : CALL PRODT (SIG33,SIG33P,VECPT,3,3) 260 : * print*,'contraintes dans rpg' 261 : * print*,sig33 262 : * 263 : * on rend les contraintes et les variables internes finales 264 : * 265 : SIGAN(1)=SIGF(1)-SIG33(1,1) 266 : SIGF(1)=SIG33(1,1) 267 : VARF(3)=DMAX1 (DOM33(1,1),0.d0) 268 : SIGAN(2)=SIGF(2)-SIG33(2,2) 269 : SIGF(2)=SIG33(2,2) 270 : VARF(4)=DMAX1 (DOM33(2,2),0.d0) 271 : SIGAN(3)=SIGF(3)-SIG33(3,3) 272 : SIGF(3)=SIG33(3,3) 273 : VARF(5)=DMAX1 (DOM33(3,3),0.d0) 274 : SIGAN(4)=SIGF(4)-SIG33(1,2) 275 : SIGF(4)=SIG33(1,2) 276 : VARF(6)=DMAX1 (DOM33(1,2),0.d0) 277 : VARF(7)=DMAX1 (DOM33(1,3),0.d0) 278 : VARF(8)=DMAX1 (DOM33(2,3),0.d0) 279 : IF(IFOUR.GE.1.OR.IFOUR.LE.-3) THEN 280 : SIGAN(5)=SIGF(5)-SIG33(1,3) 281 : SIGF(5)=SIG33(1,3) 282 : SIGAN(6)=SIGF(6)-SIG33(2,3) 283 : SIGF(6)=SIG33(2,3) 284 : ELSE 285 : SIGAN(5)=0.D0 286 : SIGAN(6)=0.D0 287 : END IF 288 : DO ISTRS=1,6 289 : VARF(ISTRS+8)=SIGAN(ISTRS) 290 : END DO 291 : * print*,'sigf',sigf(1),sigf(2),sigf(3) 292 : * print*,'sigf',sigf(4),sigf(5),sigf(6) 293 : * print*,'sigan',sigan(1),sigan(2),sigan(3) 294 : * print*,'sigan',sigan(4),sigan(5),sigan(6) 295 : SEGSUP WRKK1 296 : RETURN 297 : END 298 : 299 :
© Cast3M 2003 - All rights reserved.
Disclaimer