Download ceraca.eso
Back to the list
1 : C CERACA SOURCE CHAT 05/01/12 21:52:12 5004
2 : SUBROUTINE CERACA(WRK0,WRK1,WRK5,WTRAV,INPLAS,MFR,DT,NSTRS,NVARI,
3 : 1 PRECIS,MSOUPA,JECHER,DTT,NSSINC,INV,KERRE,
4 : 2 ICARA,IFOURB,CMATE,N2EL,N2PTEL,IB,IGAU,EPAIST,
5 : 4 NBPGAU,MELE,NPINT,NBGMAT,NELMAT,SECT,LHOOK,
6 : 5 CRIGI)
7 : C
8 : C---------------------------------------------------------------------
9 : C Objet: Calculer au cours d'un pas de temps DT, l'evolution des
10 : C variables internes a l'aide d'un schema Runge-Kutta 1.2
11 : C ---------------------------------------------------------------------
12 : C
13 : C---------------------------------------------------------------------
14 : C Entree: INPLAS type de materiau
15 : C MFR indice de la formulation mecanique(seulement massif ou coque
16 : C pour les materiaux endommageables)
17 : C DEPST(NSTRS) increment des deformations totales
18 : C SIG0(NSTRS) contraintes initiales
19 : C EPIN0(NSTRS) deformations viscoplastiques initiales
20 : C VAR0(NVARI) variables internes initiales
21 : C NVARI nombre de variables internes
22 : C XMAT(NCOMAT) materiau
23 : C XCAR(ICARA) caracteristiques geometriques
24 : C PRECIS precision relative sur SIGMA
25 : C MSOUPA nombre maximal de sous pas autorises
26 : C JECHER = 0 avancer
27 : C = 1 rechercher sortie avec DTT
28 : C IFOURB = -2 EN CONTR.PLANES
29 : C -1 EN DEFORM. PLANES
30 : C 0 EN AXISYMETRIE
31 : C 1 EN SERIE DE FOURIER
32 : C 2 EN TRIDIM
33 : C CMATE = NOM DU MATERIAU
34 : C VALMAT= TABLEAU DE CARACTERISTIQUES DU MATERIAU
35 : C VALCAR= TABLEAU DE CARACTERISTIQUES GEOMETRIQUES
36 : C N2EL = NBRE D ELEMENTS DANS SEGMENT DE HOOKE
37 : C N2PTEL= NBRE DE POINTS DANS SEGMENT DE HOOKE
38 : C IB = NUMERO DE L ELEMENT COURANT
39 : C IGAU = NUMERO DU POINT COURANT
40 : C EPAIST= EPAISSEUR
41 : C NBPGAU= NBRE DE POINTS DE GAUSS
42 : C MELE = NUMERO DE L ELEMENT FINI
43 : C NPINT = NBRE DE POINTS D INTEGRATION
44 : C NBGMAT= NBRE DE POINTS DANS SEGMENT DE CARACTERISTIQUES
45 : C NELMAT= NBRE D ELEMENTS DANS SEGMENT DE CARACTERISTIQUES
46 : C SECT = SECTION
47 : C LHOOK = TAILLE DE LA MATRICE DE HOOKE
48 : C DT pas de temps
49 : C-----------------------------------------------------------------------
50 : C
51 : C-----------------------------------------------------------------------
52 : C Sortie: SIGF(NSTRS) contraintes finales
53 : C EPINF(NSTRS) deformations viscoplastiques finales
54 : C VARF(NVARI) variables internes finales
55 : C DTT sous-increment de temps optimal (si JECHER=1)
56 : C NSSINC nombre de sous-increments si JECHER=0
57 : C INV = 1 si inversion
58 : C 0 sinon
59 : C KERRE = 0 si tout OK
60 : C <> 0 si entrees incoherentes
61 : C-----------------------------------------------------------------------
62 : C
63 : IMPLICIT INTEGER(I-N)
64 : IMPLICIT REAL*8(A-H,O-Z)
65 : -INC CCOPTIO
66 : *
67 : SEGMENT WRK0
68 : REAL*8 XMAT(NCXMAT)
69 : ENDSEGMENT
70 : *
71 : SEGMENT WRK1
72 : REAL*8 DDHOOK(LHOOK,LHOOK),SIG0(NSTRS),DEPST(NSTRS)
73 : REAL*8 SIGF(NSTRS),VAR0(NVARI),VARF(NVARI)
74 : REAL*8 DEFP(NSTRS),XCAR(ICARA)
75 : ENDSEGMENT
76 : *
77 : SEGMENT WRK5
78 : REAL*8 EPIN0(NSTRS),EPINF(NSTRS),EPST0(NSTRS)
79 : ENDSEGMENT
80 : *
81 : SEGMENT WTRAV
82 : REAL*8 DDAUX(LHOOK,LHOOK),VALMAT(NUMAT)
83 : REAL*8 VALCAR(NUCAR),DSIGT(NSTRS)
84 : REAL*8 TXR(IDIM,IDIM),DDHOMU(LHOOK,LHOOK)
85 : REAL*8 XLOC(IDIM,IDIM),XGLOB(IDIM,IDIM)
86 : REAL*8 D1HOOK(LHOOK,LHOOK),ROTHOO(LHOOK,LHOOK)
87 : ENDSEGMENT
88 :
89 : *
90 : DIMENSION SIG(6),EPSV(6),VAR(100)
91 : DIMENSION DSPT(6),XX(8),SIG1(6),EPSV1(6)
92 : DIMENSION VAR1(100),EVP1(6),VARP1(100)
93 : DIMENSION EVP2(6),VARP2(100)
94 : DIMENSION CRIGI(12)
95 : *
96 : CHARACTER*8 CMATE
97 : C**** debut ajout Eloi
98 : logical dtlibr
99 : dtlibr=.TRUE.
100 : C**** fin ajout Eloi
101 :
102 :
103 : C ==================================================
104 : C RECHERCHE DU NUMERO DE LA VARIABLE INTERNE EGALE A
105 : C 1 SI ON A ENDOMAGEMENT GENERALISE
106 : IBID2 = NVARI -1
107 : C*************************************
108 : C On cherche le numero de la propriété du matériau correspondant à ENDG
109 : IF (MFR.EQ.1.AND.IFOMOD.EQ.2) THEN
110 : IBID1 = 20
111 : ELSE
112 : IBID1 = 15
113 : ENDIF
114 : C Test pour savoir si on a dépassé la limite de déformation de
115 : C fluage
116 : IF (VAR0(IBID2).EQ.1.D0) THEN
117 : DO 1032 IJ=1,NVARI
118 : VARF(IJ) = VAR0(IJ)
119 : 1032 CONTINUE
120 : DO 1022 K=1,NSTRS
121 : EPINF(K)=EPIN0(K)
122 : SIGF(K) = 0.D0
123 : 1022 CONTINUE
124 : GO TO 998
125 : ENDIF
126 : DO 1010 I =1,NSTRS
127 : TOTO = DABS(EPIN0(I))
128 : IF(TOTO.GE.XMAT(IBID1)) THEN
129 : DO 1021 K=1,NSTRS
130 : EPINF(K)=EPIN0(K)
131 : SIGF(K) = 0.D0
132 : 1021 CONTINUE
133 : DO 1030 IJ=1,(NVARI-1)
134 : VARF(IJ) = VAR0(IJ)
135 : 1030 CONTINUE
136 : VARF(IBID2) = 1.D0
137 : GO TO 998
138 : ENDIF
139 : 1010 CONTINUE
140 : C*************************************
141 :
142 :
143 : C =========================================================
144 : KERRE = 0
145 : IF(MFR.NE.1.AND.MFR.NE.3.OR.IFOURB.EQ.1)THEN
146 : KERRE = 99
147 : RETURN
148 : ENDIF
149 :
150 : IF(MFR.EQ.3) THEN
151 : THICK = XCAR(1)
152 : ALFA = XCAR(2)
153 : ENDIF
154 : BORNE = 2.0
155 : RMAX = 1.3
156 : RMIN = 0.7
157 : DIV = 7.0
158 : FAC = 3.0
159 : C
160 : C CALCUL DES INCREMENTS DE CONTRAINTES ELASTIQUES
161 : C
161 : CALL CALSIG(DEPST,DDAUX,NSTRS,CMATE,VALMAT,
162 : 1 VALCAR,N2EL,N2PTEL,MFR,IFOURB,IB,IGAU,EPAIST,
163 : 2 NBPGAU,MELE,NPINT,NBGMAT,NELMAT,SECT,LHOOK,TXR,XLOC,
164 : 3 XGLOB,D1HOOK,ROTHOO,DDHOMU,CRIGI,DSIGT,IRTD)
165 : C PRINT *,'DSIGT =',DSIGT
166 : *
167 : IF(IRTD.NE.1) THEN
168 : KERRE=69
169 : GOTO 998
170 : ENDIF
171 : *
172 : IF(MFR.EQ.3) THEN
173 : DO 10 I=1,NSTRS/2
174 : SIG0( I) = SIG0( I)/THICK
175 : SIG0(NSTRS/2+I) = SIG0(NSTRS/2+I)*6.0D0/THICK/THICK
176 : DSIGT( I) = DSIGT( I)/THICK
177 : DSIGT(NSTRS/2+I)= DSIGT(NSTRS/2+I)*6.0D0/THICK/THICK
178 : 10 CONTINUE
179 : IF(IFOURB.EQ.-2) THEN
180 : SIG0 (2) =0.0
181 : SIG0 (4) =0.0
182 : DSIGT(2) =0.0
183 : DSIGT(4) =0.0
184 : ENDIF
185 : ENDIF
186 : C
187 : C------------------------------------------
188 : C CONTROLE DE LA COHERENCE DES ENTREES
189 : C Supprimée dans ce cas
190 : C------------------------------------------
191 : IF (DT.LE.0.0) KERRE = 414
192 : MOTERR(1:8) = ' CONST '
193 : C IF (INCOMAT.NE. 15)KERRE = 146
194 : C IF (IFOURB.EQ.-2.AND.NCOMAT.NE. 9)KERRE = 146
195 : IF(IFOURB.NE.2.AND.IFOURB.NE.-2.AND.IFOURB.NE.0
196 : & .AND.IFOURB.NE.-1) THEN
197 : KERRE = 194
198 : MOTERR(1:8) = 'FLUAGE'
199 : ENDIF
200 : IF(MFR.EQ.1.AND.IFOMOD.EQ.2) THEN
201 : IBID = 15
202 : ELSE
203 : IBID = 10
204 : ENDIF
205 : XMAX=XMAT(IBID)
206 : *
207 : * TEST SUR XMAX MILL 8/3/91
208 : *
209 : IF (XMAX.EQ.0.D0) THEN
210 : XMAX=XMAT(1)*1.D-3
211 : ENDIF
212 : *
213 : 30 CONTINUE
214 : C
215 : C-----------------------------
216 : IF (KERRE.NE.0) THEN
217 : GOTO 999
218 : ENDIF
219 : C
220 : C===========================================================
221 : C A PARTIR DE MAINTENANT, LES DEFORMATIONS
222 : C DE CISAILLEMENT NE SONT PLUS
223 : C DEFINIES PAR DES GAMA.
224 : C ON DIVISE DONC LES TERMES DE CISAILLEMENT PAR 2.
225 : C SEULES LES FORMULATIONS SUIVANTES SONT ACCEPTEES PAR CERACARO:
226 : C MFR=1 (MASSIF)
227 : C MFR=3 (COQUES MINCES)
228 : C
229 : C
230 : C Cas de la formulation massive
231 : C Les termes de cisaillement apparaissent
232 : C au delà de la troisieme composante
233 : C
234 : IF (MFR.EQ.1) THEN
235 : DO 11 I=1,NSTRS
236 : A=1.D0
237 : IF (I.GT.3) A=2.D0
238 : EPIN0(I)=EPIN0(I)/A
239 : 11 CONTINUE
240 :
241 : C Cas des coques minces
242 : C Les termes de cisaillement apparaissent
243 : C pour la troisieme et la sixieme composante
244 : C uniquement dans les cas de calculs
245 : C tridimensionnels ou d'analyse de Fourier
246 : C
247 : ELSE IF (MFR.EQ.3) THEN
248 : IF (IFOURB.EQ.2) THEN
249 : DO 13 I=1,NSTRS
250 : A=1.D0
251 : IF (I.EQ.3) A=2.D0
252 : IF (I.EQ.6) A=2.D0
253 : EPIN0(I)=EPIN0(I)/A
254 : 13 CONTINUE
255 : ENDIF
256 : ENDIF
257 :
258 : C
259 : C===========================================================
260 : C
261 : C ----------------
262 : C INITIALISATION
263 : C ----------------
264 : C**** debut ajout Eloi
265 : ITERO = 0
266 : 1200 CONTINUE
267 :
268 : itero = 1 + itero
269 : if (itero.ne.1) THEN
270 : c write(6,*) 'itero ib igau', itero,ib,igau
271 : dtlibr = .true.
272 : precis = precis * 7.
273 : c write(6,*) ' precision modifiée ', precis
274 : if (itero.gt.3) then
275 : **** kerre = 460
276 : kerre = 268
277 : return
278 : endif
279 : endif
280 : C**** fin ajout Eloi
281 : DTLEFT = DT
282 : TAU = DTLEFT
283 : ASIG = SQRT(PROCON(SIG0,SIG0,NSTRS))
284 : ERRABS = PRECIS*ASIG
285 : IF (XMAX.GT.ASIG) ERRABS = PRECIS*XMAX
286 : DO 40 I=1,NSTRS
287 : SIG(I) = SIG0(I)
288 : EPSV(I) = EPIN0(I)
289 : DSPT(I) = DSIGT(I)/DT
290 : 40 CONTINUE
291 : C
292 :
293 : DO 48 I=1,NVARI
294 : 48 VAR(I)=VAR0(I)
295 : C
296 : C ---------------------------------------------------------------------
297 : INV = 0
298 : NSSINC = 0
299 : C ---------------------------------------------------------------------
300 : C DEBUT DES ITERATIONS EN SSINCREMENTS /FIN SI DTLEFT = 0
301 : C ---------------------------------------------------------------------
302 : 70 NSSINC = NSSINC + 1
303 : IF (NSSINC.GT.MSOUPA) THEN
304 : C**** debut ajout Eloi
305 : DTLIBR=.FALSE.
306 : GOTO 1200
307 : C**** fin ajout Eloi
308 : C*** KERRE = 460
309 : C*** GOTO 999
310 : ENDIF
311 : C
312 : C---------------------------------------------------------------------
313 : C START OF CALCULATIONS
314 : C_____________________________________________________________________
315 :
316 : IF(MFR.NE.3) THEN
317 :
317 : CALL INCR11(TAU,SIG,EPSV,VAR,EVP1,VARP1,XMAT,NSTRS,NVARI,MFR)
318 : ELSE IF(MFR.EQ.3) THEN
318 : CALL INCR33(TAU,SIG,EPSV,VAR,XMAT,EVP1,VARP1,ALFA,NSTRS,NVARI)
319 : ENDIF
320 : C
321 : NITERA = 0
322 :
323 : C --------------------------------------------------------------------
324 : C DEBUT DES ITERATIONS SUR TAU OPTIMAL /FIN SI RA PETIT
325 : C --------------------------------------------------------------------
326 : 80 NITERA = NITERA + 1
327 : IF(MFR.EQ.3) GO TO 150
328 :
328 : CALL ADVA11(TAU,SIG,EPSV,VAR,SIG1,EPSV1,VAR1,DSPT,EVP1,VARP1,
329 : & XMAT,NSTRS,NVARI,MFR)
330 :
330 : CALL INCR11(TAU,SIG1,EPSV1,VAR1,EVP2,VARP2,XMAT,NSTRS,NVARI,
331 : & MFR)
332 : DO 110 I=1,NSTRS
333 : 110 EVP2(I) = 0.5*(EVP1(I)+EVP2(I))
334 : CCC Eloi : on peut effectuer la boucle suivante sur
335 : CCC les variables internes propres a Ottosen
336 : CCC qui n'ont pas ete modifiees
337 : DO 120 I=1,NVARI
338 : VARP2(I)= 0.5*(VARP1(I) + VARP2(I))
339 : 120 CONTINUE
339 : CALL ADVA11(TAU,SIG,EPSV,VAR,SIGF,EPINF,VARF,DSPT,EVP2,VARP2,
340 : & XMAT,NSTRS,NVARI,MFR)
341 :
342 : GO TO 250
343 : C _____________________________________________________________________
344 : 150 CONTINUE
345 : C----------------------------------------------------------------------
346 : C CALCULATIONS FOR GENERALISED STRESS/STRAIN FORMULATIONS
347 : C----------------------------------------------------------------------
347 : CALL ADVA22(TAU,SIG,EPSV,VAR,SIG1,EPSV1,VAR1,DSPT,EVP1,VARP1,
348 : & XMAT,NSTRS,NVARI)
348 : CALL INCR33(TAU,SIG1,EPSV1,VAR1,XMAT,EVP2,VARP2,ALFA,NSTRS,
349 : & NVARI)
350 : DO 160 I=1,NSTRS,1
351 : 160 EVP2(I) = 0.5*(EVP1(I)+EVP2(I))
352 : CCC Eloi : on peut effectuer la boucle suivante sur
353 : CCC les variables internes propres a Ottosen
354 : CCC qui n'ont pas ete modifiees
355 : DO 170 I=1,NVARI
356 : VARP2(I)= 0.5*(VARP1(I) + VARP2(I))
357 : 170 CONTINUE
357 : CALL ADVA22(TAU,SIG,EPSV,VAR,SIGF,EPINF,VARF,DSPT,
358 : & EVP2,VARP2,XMAT,NSTRS,NVARI)
359 : 250 CONTINUE
360 :
361 : C ---------------------------------------------------------------------
362 : C CALCUL DU RAPPORT : ERREUR CALCULEE / ERREUR ADMISE
363 : C ---------------------------------------------------------------------
364 : DO 260 I=1,NSTRS
365 : 260 XX(I) = SIGF(I)-SIG1(I)
366 : RA = SQRT(PROCON(XX,XX,NSTRS))/(ERRABS)
367 : SQRA = SQRT(RA)
368 :
369 : C ---------------------------------------------------------------------
370 : C TEST DE FIN D'ITERATIONS / MISE A JOUR DE TAU /OPTION JECHER
371 : C DIV =7 BORNE = 2
372 : C SI SQRA>7 TAU = TAU/7 ET NOUVEL ESSAI
373 : C SI 2RMAX on vise RA = 1
451 : C fin des boucles en ss increments si tau = dtleft
452 : C --------------------------------------------------------------------
453 : C
454 : IF ( TAU.LT.DTLEFT ) THEN
455 : DTLEFT = DTLEFT - TAU
456 : IF ( FAC*SQRA.LT.1.D0) THEN
457 : TAU=TAU*FAC
458 : ELSEIF ( (SQRA.LT.RMIN).OR.(SQRA.GT.RMAX) ) THEN
459 : TAU=TAU/SQRA
460 : ENDIF
461 : IF (TAU.GT.DTLEFT) TAU = DTLEFT
462 : GOTO 70
463 : ENDIF
464 : C
465 : IF (ABS(TAU-DTLEFT).GT.(TAU/1000.)) THEN
466 : WRITE ( IOIMP,* ) ' PROBLEME TAU > DTLEFT '
467 : KERRE = 223
468 : ENDIF
469 : C-----------------------------------------------------------------------
470 : 999 CONTINUE
471 : IF(MFR.EQ.3) THEN
472 : DO 1000 I=1,NSTRS/2
473 : SIGF( I) =SIGF( I)*THICK
474 : SIGF(NSTRS/2+I) =SIGF(NSTRS/2+ I)*THICK*THICK/6.0
475 : * DSIGT( I)=DSIGT( I)*THICK
476 : * DSIGT(NSTRS/2+I)=DSIGT(NSTRS/2+I)*THICK*THICK/6.0
477 : 1000 CONTINUE
478 : ENDIF
479 :
480 : C
481 : C===========================================================
482 : C RETOUR A LA DEFINITION NORMALE DES DEFORMATIONS
483 : C A SAVOIR: LES DEFORMATIONS DE CISAILLEMENT SONT
484 : C DEFINIES PAR DES GAMA.
485 : C SEULES LES FORMULATIONS SUIVANTES SONT ACCEPTEES PAR Ceracaro:
486 : C MFR=1 (MASSIF)
487 : C MFR=3 (COQUES MINCES)
488 : C
489 : C Cas de la formulation massive
490 : C Les termes de cisaillement apparaissent
491 : C au delà de la troisieme composante
492 : C
493 : IF (MFR.EQ.1) THEN
494 : DO 14 I=1,NSTRS
495 : A=1.D0
496 : IF (I.GT.3) A=2.D0
497 : EPIN0(I)=EPIN0(I)*A
498 : EPINF(I)=EPINF(I)*A
499 : 14 CONTINUE
500 : C
501 : C Cas des coques minces
502 : C Les termes de cisaillement apparaissent
503 : C pour la troisieme et la sixieme composante
504 : C uniquement dans les cas de calculs
505 : C tridimensionnels ou d'analyse de Fourier
506 : C
507 : ELSE IF (MFR.EQ.3) THEN
508 : IF (IFOURB.EQ.2) THEN
509 : DO 16 I=1,NSTRS
510 : A=1.D0
511 : IF (I.EQ.3) A=2.D0
512 : IF (I.EQ.6) A=2.D0
513 : EPIN0(I)=EPIN0(I)*A
514 : EPINF(I)=EPINF(I)*A
515 : 16 CONTINUE
516 : ENDIF
517 : ENDIF
518 : C
519 : C===========================================================
520 : C
521 : 998 RETURN
522 : END
523 :
524 :
525 :
526 :
527 :
528 :
529 :
530 :
531 :
532 :
533 :
534 :
535 :
536 :