Download umapol.eso
Back to the list
1 : C UMAPOL SOURCE CHAT 06/07/13 21:16:14 5488
2 : SUBROUTINE UMAPOL ( STRESS, STATEV, DDSDDE,
3 : & STRAN, DSTRAN, DTIME,
4 : & TEMP, DTEMP, PREDEF, DPRED,
5 : & NTENS, NSTATV,
6 : & PROPS, NPROPS, COORDS,
7 : & NOEL, NPT, KINC )
8 : IMPLICIT INTEGER(I-N)
9 : IMPLICIT REAL*8(A-H,O-Z)
10 : C-----------------------------------------------------------------------
11 : C
12 : C DESCRIPTION FONCTIONNELLE :
13 : C -------------------------
14 : C
15 : C Dependance du module Utilisateur UMAT
16 : C Appelee dans le cas suivant :
17 : C
18 : C ILOI = 21 : Modele 'FLUAGE' 'POLYNOMIAL' de CAST3M, INPLAS=21
19 : C RESTRICTION A LA FORMULATION MASSIVE
20 : C
21 : C Composantes de materiau
22 : C -----------------------
23 : C PROPS ( 1) : 'YOUN' PROPS ( 5) : 'SMAX'
24 : C PROPS ( 2) : 'NU ' PROPS ( 6) : 'AF0 '
25 : C PROPS ( 3) : 'RHO ' PROPS ( 7) : 'AF1 '
26 : C PROPS ( 4) : 'ALPH' PROPS ( 8) : 'AF2 '
27 : C PROPS ( 9) : 'AF3 '
28 : C PROPS (10) : 'AF4 '
29 : C PROPS (11) : 'AF5 '
30 : C PROPS (12) : 'AF6 '
31 : C Variables internes
32 : C ------------------
33 : C STATEV( 1) : 'PXX ' STATEV(17) : 'ENXX'
34 : C STATEV( 2) : 'PYY ' STATEV(18) : 'ENYY'
35 : C STATEV( 3) : 'PZZ ' STATEV(19) : 'ENZZ'
36 : C STATEV( 4) : 'PXY ' STATEV(20) : 'GNXY'
37 : C STATEV( 5) : 'PYZ ' STATEV(21) : 'GNXZ'
38 : C STATEV( 6) : 'PZX ' STATEV(22) : 'GNYZ'
39 : C STATEV( 7) : 'VNXX'
40 : C STATEV( 8) : 'VNYY'
41 : C STATEV( 9) : 'VNZZ'
42 : C STATEV(10) : 'VNXY'
43 : C STATEV(11) : 'VNYZ'
44 : C STATEV(12) : 'VNZX'
45 : C STATEV(13) : 'EPCE'
46 : C STATEV(14) : 'EPSH'
47 : C STATEV(15) : 'LL '
48 : C STATEV(16) : 'EPSE'
49 : C
50 : C-----------------------------------------------------------------------
51 : C Arguments de l'interface
52 : C
53 : INTEGER NTENS, NSTATV, NPROPS,
54 : & NOEL, NPT, KINC
55 : C
56 : REAL*8 STRESS(NTENS), STATEV(*),
57 : & DDSDDE(NTENS,NTENS),
58 : & STRAN(NTENS), DSTRAN(NTENS),
59 : & DTIME,
60 : & TEMP, DTEMP, PREDEF(*), DPRED(*),
61 : & PROPS(NPROPS),
62 : & COORDS(3)
63 : C
64 : C Variables locales
65 : C
66 : INTEGER NVARI
67 : PARAMETER (NVARI=16)
68 : C
69 : C Tableaux de travail
70 : C
71 : REAL*8 VAR(NVARI),VAR1(NVARI),VARP1(NVARI),VARP2(NVARI),
72 : & VAR12(NVARI),VARP3(NVARI),VARP4(NVARI),VAR13(NVARI),
73 : & SIG(8),SIG1(8),SIG12(8),SIG13(8),
74 : & DSPT(8),XX(8),
75 : & EVP1(8),EVP2(8),EVP3(8),EVP4(8),
76 : & EPSV(8),EPSV1(8),EPSV12(8),EPSV13(8)
77 : C
78 : REAL*8 EPIN0(8),DSIGT(8),VARF(NVARI),SIGF(8),EPINF(8)
79 : C
80 : C Variables scalaires
81 : C
82 : LOGICAL DTLIBR
83 : REAL*8 DT, DTLEFT, BORNE, RMAX, RMIN, DIV, FAC,
84 : & YOU, XNU, AUX0, AUX, AUX1, GEGE,
85 : & XMAX,
86 : CCCCC& TI0, TI1, TPOINT, DELTAT, TI12, FII0, FII1, FII12,
87 : & TI1, TI12, FII0, FII1, FII12,
88 : & TAU, ASIG, ERRABS, TAU2, RA, SQRA, RA1
89 : INTEGER MSOUPA, I, J, ITERO, NSSINC, NITERA,
90 : & INPLAS, MFR1, IFOURB
91 : C
92 : REAL*8 UN, DEUX, UNDEMI
93 : PARAMETER (UN=1.D0,DEUX=2.D0,UNDEMI=.5D0)
94 : C
95 : REAL*8 PRECIS
96 : C
97 : REAL*8 PROCON
98 : EXTERNAL PROCON
99 : C
100 : C------------------- Debut du code executable --------------------------
101 : C
102 : C=======================================================================
103 : C 1 - INITIALISATIONS
104 : C Pas de temps et parametres de pilotage des iterations internes
105 : C=======================================================================
106 : C
107 : PRECIS=1.d-8
108 : DTLIBR=.TRUE.
109 : C
110 : DT = DTIME
111 : IF (DT.LT.0.0) THEN
112 : KINC = -211
113 : RETURN
114 : ENDIF
115 : IF (DT.EQ.0.0) DT = 1.e-20
116 : C
117 : DTLEFT = DT
118 : BORNE = 2.0
119 : RMAX = 1.3
120 : RMIN = 0.7
121 : DIV = 7.0
122 : FAC = 3.0
123 : C
124 : MSOUPA = 200
125 : C
126 : INPLAS = 21
127 : MFR1 = 1
128 : IFOURB = 2
129 : FII0 = 0.0D0
130 : FII1 = 0.0D0
131 : FII12 = 0.0D0
132 : TI1 = 0.0D0
133 : TI12 = 0.0D0
134 : C
135 : C=======================================================================
136 : C 2 - Prediction elastique de l'increment de contraintes
137 : C La sequence ci-dessous se substitue a l'appel a CALSIG
138 : C Le calcul des coefficients de la matrice de Hooke reprend la
139 : C sequence d'operations de HOOKIS
140 : C=======================================================================
141 : C
142 : YOU = PROPS(1)
143 : XNU = PROPS(2)
144 : C
145 : AUX0 = UN/((UN+XNU)*(UN-DEUX*XNU))
146 : AUX = YOU*AUX0*(UN-XNU)
147 : AUX1 = YOU*AUX0*XNU
148 : GEGE = YOU*UNDEMI/(UN+XNU)
149 : C
150 : DO 10 J=1,NTENS
151 : DO 11 I=1,NTENS
152 : DDSDDE(I,J)=0.0D0
153 : 11 CONTINUE
154 : 10 CONTINUE
155 : C
156 : DDSDDE(1,1) = AUX
157 : DDSDDE(1,2) = AUX1
158 : DDSDDE(1,3) = AUX1
159 : C
160 : DDSDDE(2,1) = AUX1
161 : DDSDDE(2,2) = AUX
162 : DDSDDE(2,3) = AUX1
163 : C
164 : DDSDDE(3,1) = AUX1
165 : DDSDDE(3,2) = AUX1
166 : DDSDDE(3,3) = AUX
167 : C
168 : DDSDDE(4,4) = GEGE
169 : DDSDDE(5,5) = GEGE
170 : DDSDDE(6,6) = GEGE
171 : C
172 : DO 20 I=1,NTENS
173 : DSIGT(I)=0.0D0
174 : 20 CONTINUE
175 : DO 21 J=1,NTENS
176 : DO 22 I=1,NTENS
177 : DSIGT(I)=DSIGT(I)+DDSDDE(I,J)*DSTRAN(J)
178 : 22 CONTINUE
179 : 21 CONTINUE
180 : C
181 : C=======================================================================
182 : C 3 - Pre-traitement des deformations inelastiques au debut du pas
183 : C Passage gamma -> epsilon pour les termes extradiagonaux
184 : C=======================================================================
185 : C
186 : DO 34 J=1,3
187 : EPIN0(J)=STATEV(NVARI+J)
188 : 34 CONTINUE
189 : DO 35 J=4,6
190 : EPIN0(J)=0.5D0*STATEV(NVARI+J)
191 : 35 CONTINUE
192 : C
193 : C=======================================================================
194 : C 4 - TEST SUR XMAX MILL 8/3/91
195 : C=======================================================================
196 : C
197 : XMAX=PROPS(5)
198 : IF (XMAX.EQ.0.D0) XMAX=PROPS(1)*1.D-3
199 : C
200 : C=======================================================================
201 : C 5 - INITIALISATIONS AVANT ITERATIONS EN SOUS-INCREMENTS
202 : C=======================================================================
203 : C
204 : ITERO = 0
205 : 6543 CONTINUE
206 : C
207 : ITERO = 1 + ITERO
208 : IF ( ITERO.NE.1) THEN
209 : DTLIBR = .TRUE.
210 : PRECIS = PRECIS * 7.D0
211 : IF (ITERO.GT.3) THEN
212 : KINC = -212
213 : RETURN
214 : ENDIF
215 : ENDIF
216 : C
217 : DTLEFT = DT
218 : TAU = DTLEFT
219 : CCCCC TI0=TEMP
220 : CCCCC TI1=TEMP+DTEMP
221 : CCCCC TPOINT=DTEMP/DT
222 : ASIG = SQRT(PROCON(STRESS,STRESS,NTENS))
223 : ERRABS = PRECIS*ASIG
224 : IF (XMAX.GT.ASIG) ERRABS = PRECIS*XMAX
225 : C
226 : DO 40 I=1,NTENS
227 : SIG(I) = STRESS(I)
228 : EPSV(I) = EPIN0(I)
229 : DSPT(I) = DSIGT(I)/DT
230 : 40 CONTINUE
231 : C
232 : DO 50 I=1,NVARI
233 : VAR(I)=STATEV(I)
234 : 50 CONTINUE
235 : C
236 : C=======================================================================
237 : C 6 - ITERATIONS EN SSINCREMENTS /FIN SI DTLEFT = 0
238 : C=======================================================================
239 : C
240 : NSSINC = 0
241 : NITERA = 0
242 : 60 CONTINUE
243 : C
244 : NSSINC = NSSINC + 1
245 : IF (NSSINC.GT.MSOUPA) THEN
246 : DTLIBR=.FALSE.
247 : GOTO 6543
248 : ENDIF
249 : C
250 : CCCCC DELTAT=TPOINT*TAU
251 : CCCCC TI1=TI0+DELTAT
251 : CALL INCRE1(TAU,SIG,EPSV,VAR,EVP1,VARP1,PROPS,NTENS,NVARI,
252 : & INPLAS,NPROPS,MFR1,FII0,FII1,TI1)
253 : C
254 : C ---------------------------------------------------------------
255 : C DEBUT DES ITERATIONS SUR TAU OPTIMAL /FIN SI RA PETIT
256 : C
257 : NITERA = 0
258 : 70 CONTINUE
259 : C
260 : NITERA = NITERA + 1
261 : C
262 : C ____________________________________________________________
263 : TAU2=TAU*0.5D0
263 : CALL ADVAN1(TAU2,SIG,EPSV,VAR,SIG1,EPSV1,VAR1,DSPT,EVP1,
264 : & VARP1,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
265 : C
266 : C ____________________________________________________________
267 : CCCCC DELTAT=TPOINT*TAU
268 : CCCCC TI1=TI0+DELTAT
268 : CALL INCRE1(TAU,SIG1,EPSV1,VAR1,EVP2,VARP2,PROPS,NTENS,
269 : & NVARI,INPLAS,NPROPS,MFR1,FII0,FII1,TI1)
270 : DO 701 I=1,NTENS
271 : EVP2(I) = 0.5D0*(EVP1(I)+EVP2(I))
272 : 701 CONTINUE
273 : DO 702 I=1,2*NTENS+2
274 : VARP2(I)= 0.5D0*(VARP1(I) + VARP2(I))
275 : 702 CONTINUE
276 : DO 703 I=2*NTENS+4,NVARI
277 : VARP2(I)= 0.5D0*(VARP1(I) + VARP2(I))
278 : 703 CONTINUE
279 : C
280 : C ____________________________________________________________
281 : CCCCC DELTAT=TPOINT*TAU2
282 : CCCCC TI12=TI0+DELTAT
282 : CALL ADVAN1(TAU2,SIG,EPSV,VAR,SIG12,EPSV12,VAR12,DSPT,EVP2,
283 : & VARP2,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
283 : CALL INCRE1(TAU2,SIG12,EPSV12,VAR12,EVP3,VARP3,PROPS,NTENS,
284 : & NVARI,INPLAS,NPROPS,MFR1,FII0,FII12,TI12)
285 : C
286 : C ____________________________________________________________
287 : CCCCC TI1=TI12+DELTAT
287 : CALL ADVAN1(TAU2,SIG12,EPSV12,VAR12,SIG13,EPSV13,VAR13,DSPT,
288 : & EVP3,VARP3,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
288 : CALL INCRE1(TAU2,SIG13,EPSV13,VAR13,EVP4,VARP4,PROPS,NTENS,
289 : & NVARI,INPLAS,NPROPS,MFR1,FII12,FII1,TI1)
290 : DO 711 I=1,NTENS
291 : EVP4(I) = 0.5D0*(EVP3(I)+EVP4(I))
292 : 711 CONTINUE
293 : DO 712 I=1,2*NTENS+2
294 : VARP4(I)= 0.5D0*(VARP3(I) + VARP4(I))
295 : 712 CONTINUE
296 : DO 713 I=2*NTENS+4,NVARI
297 : VARP4(I)= 0.5D0*(VARP3(I) + VARP4(I))
298 : 713 CONTINUE
298 : CALL ADVAN1(TAU2,SIG12,EPSV12,VAR12,SIGF,EPINF,VARF,DSPT,
299 : & EVP4,VARP4,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
300 : C
301 : C ____________________________________________________________
302 : CCCCC DELTAT=TPOINT*TAU2
303 : CCCCC TI12=TI0+DELTAT
303 : CALL INCRE1(TAU2,SIGF,EPINF,VARF,EVP4,VARP4,PROPS,NTENS,
304 : & NVARI,INPLAS,NPROPS,MFR1,FII0,FII12,TI12)
305 : DO 721 I=1,NTENS
306 : EVP2(I) = (EVP1(I)+EVP4(I))/6.D0+EVP3(I)*2.D0/3.D0
307 : 721 CONTINUE
308 : DO 722 I=1,2*NTENS+2
309 : VARP2(I)=(VARP1(I) + VARP4(I))/6.D0+VARP3(I)*2.D0/3.D0
310 : 722 CONTINUE
311 : DO 723 I=2*NTENS+4,NVARI
312 : VARP2(I)= (VARP1(I) + VARP4(I))/6.D0+VARP3(I)*2.D0/3.D0
313 : 723 CONTINUE
313 : CALL ADVAN1(TAU,SIG,EPSV,VAR,SIG1,EPSV1,VAR1,DSPT,EVP2,
314 : & VARP2,PROPS,NTENS,NVARI,IFOURB,INPLAS,MFR1)
315 : C
316 : C ____________________________________________________________
317 : C CALCUL DU RAPPORT : ERREUR CALCULEE / ERREUR ADMISE
318 : C
319 : DO 730 I=1,NTENS
320 : XX(I) = SIGF(I)-SIG1(I)
321 : 730 CONTINUE
322 : RA = SQRT(PROCON(XX,XX,NTENS))/(ERRABS)
323 : SQRA = SQRT(RA)
324 : C
325 : C ____________________________________________________________
326 : C TEST DE FIN D'ITERATIONS / MISE A JOUR DE TAU
327 : C DIV =7 BORNE = 2
328 : C SI SQRA>7 TAU = TAU/7 ET NOUVEL ESSAI
329 : C SI 2RMAX on vise RA = 1
368 : C Fin des boucles en ss increments si tau = DTLEFT
369 : C
370 : IF ( TAU.LT.DTLEFT ) THEN
371 : DTLEFT = DTLEFT - TAU
372 : IF ( FAC*SQRA.LT.1.D0) THEN
373 : TAU=TAU*FAC
374 : ELSE IF ( (SQRA.LT.RMIN).OR.(SQRA.GT.RMAX) ) THEN
375 : TAU=TAU/SQRA
376 : ENDIF
377 : IF (TAU.GT.DTLEFT) TAU = DTLEFT
378 : CCCCC TI0=TI1
379 : CCCCC DELTAT=TPOINT*TAU
380 : CCCCC TI1=TI0+DELTAT
381 : GOTO 60
382 : ENDIF
383 : C
384 : IF (ABS(TAU-DTLEFT).GT.(TAU/1000.)) THEN
385 : KINC = -213
386 : RETURN
387 : ENDIF
388 : C
389 : C=======================================================================
390 : C 7 - STOCKAGE DES RESULTATS AVANT SORTIE
391 : C=======================================================================
392 : C
393 : C Contraintes
394 : C
395 : DO 100 J=1,NTENS
396 : STRESS(J)=SIGF(J)
397 : 100 CONTINUE
398 : C
399 : C Variables internes
400 : C
401 : DO 110 J=1,NVARI
402 : STATEV(J)=VARF(J)
403 : 110 CONTINUE
404 : C
405 : C Deformations inelastiques
406 : C Passage epsilon -> gamma pour les termes extradiagonaux
407 : C
408 : DO 120 J=1,3
409 : STATEV(NVARI+J)=EPINF(J)
410 : 120 CONTINUE
411 : DO 121 J=4,6
412 : STATEV(NVARI+J)=2.0D0*EPINF(J)
413 : 121 CONTINUE
414 : C
415 : RETURN
416 : END
417 :
418 :
419 :
420 :
421 :