Download cmazar.eso
Back to the list
1 : C CMAZAR SOURCE AM 10/04/28 21:15:01 6672
2 : SUBROUTINE CMAZAR (WRK52,WRK53,WRK54,WRKK2,NSTRS1,NVARI,
3 : 1 ICARA,JDIM,IFOUR2)
4 : C MAZARS SOURCE AM 98/12/23 21:38:30 3409
5 : C
6 : C
7 : C variables en entree
8 : C
9 : C
10 : C WRK0 pointeur sur un segment deformation au pas precedent
11 : C
12 : C WRK1 pointeur sur un segment increment de deformation
13 : C
14 : C WRKK2 pointeur sur un segment variables internes au pas precedent
15 : C
16 : C WRK5 pointeur sur un segment de deformations inelastiques
17 : C
18 : C XMATER constantes du materiau
19 : C
20 : C NSTRS1 nombre de composantes dans les vecteurs des contraintes
21 : C et les vecteurs des deformations
22 : C
23 : C NVARI nombre de variables internes (doit etre egal a 2)
24 : C
25 : C NMATT nombre de constantes du materiau
26 : C
27 : C ISTEP flag utilise pour separer les etapes dans un calcul non local
28 : C ISTEP=0 -----> calcul local
29 : C ISTEP=1 -----> calcul non local etape 1 on calcule les seuils
30 : C ISTEP=2 -----> calcul non local etape 2 on continue le calcul
31 : C a partir des seuils moyennes
32 : C
33 : C Modif L.Bode - 14/10/92
34 : C Nouveaux parametres en entree
35 : C JDIM Dimension de travail
36 : C ( Coques JDIM =2 , Massifs JDIM = IDIM )
37 : C IFOUR2 Type de formulation
38 : C ( Coques IFOUR2 = -2 => contraintes planes ,
39 : C Massifs IFOUR2 = IFOUR)
40 : C
41 : C variables en sortie
42 : C
43 : C VARINF variables internes finales
44 : C
45 : C SIGMAF contraintes finales
46 : C
47 : C C. LA BORDERIE MARS 1992
48 : C declaration des variables
49 : C
50 : C
51 : IMPLICIT INTEGER(I-N)
52 : IMPLICIT REAL*8(A-H,O-Z)
53 : -INC CCOPTIO
54 : -INC DECHE
55 : *
56 : SEGMENT WRKK2
57 : REAL*8 EPSILI(NSTRS1)
58 : ENDSEGMENT
59 : C Modif L.Bode - 14/10/92
60 : SEGMENT WRK3
61 : REAL*8 EPSILO(NSTRS1)
62 : ENDSEGMENT
63 : C Definition dynamique de EPSILO
64 : C Fin modif L.Bode
65 : INTEGER NSTRS1,NVARI
66 : REAL*8 EPS33(3,3),EPSIPP(3),EPSILT(3),VALP33(3,3)
67 : REAL*8 SIGP(3),SIGPT(3),SIGPC(3),TRSIGT,TRSIGC
68 : REAL*8 YOUN,XNU,EPSD0,ACOM,BCOM,ATRA,BTRA,BETA
69 : INTEGER ISTRS,JSTRS,KCAS,IRTD
70 : REAL*8 XZERO,UN,DEUX,XPETIT
71 : REAL*8 DINI,D,DT,DC,EPSTIL,EPSTIM,ALFAT,ALFAC,GAMMA
72 : PARAMETER (XZERO=0.D0 , UN=1.D0 , DEUX=2.D0, XPETIT=1.D-12)
73 :
74 : C
75 : C
76 : C recuperation des variables initiales dans les tableaux
77 : C
78 : C
79 : N=NSTRS1
80 : CMATE = 'ISOTROPE'
81 : YOUN = XMAT(1)
82 : XNU = XMAT(2)
83 : EPSD0= XMAT(5)
84 : ACOM = XMAT(6)
85 : BCOM = XMAT(7)
86 : ATRA = XMAT(8)
87 : BTRA = XMAT(9)
88 : BETA = XMAT(10)
89 : DINI = VAR0(2)
90 :
91 : C
92 : C calcul des seuils
93 : C
94 : C
95 : C calcul de la deformation totale
96 : C
97 : SEGINI WRK3
98 : DO 100 ISTRS=1,NSTRS1
99 : EPSILO(ISTRS)=EPSILI(ISTRS)+DEPST(ISTRS)
100 : 100 CONTINUE
101 : C
102 : C calcul des deformations principales
103 : C
104 : C
105 : C on reecrit les deformations sous forme matricielle
106 : C
107 : C Modif L.Bode - 14/10/92
108 : C Rajout de IFOUR2 en argument de ENDOCA
109 : * print*,'on appelle ENDOCB'
109 : CALL ENDOCB (EPSILO,EPS33,2,IFOUR2)
110 : * print*,'apres endocb'
111 : C Fin modif L.Bode
112 : C
113 : C et on diagonalise
114 : C
115 : C Modif L.Bode - 14/10/92
116 : C Pour les elts Coques, on travaille en contraintes planes => JDIM =2
117 : C Pour les elts Massifs JDIM =IDIM
118 : * print*,'avant JACOB3'
118 : CALL JACOB3 (EPS33,JDIM,EPSIPP,VALP33)
119 : * print*,'apres JACOB3'
120 : C Fin modif L.Bode
121 : IF (ISTEP .EQ. 0 .OR. ISTEP.EQ.2) THEN
122 : C
123 : C on calcule la matrice de hooke et les contraintes ppales
124 : C
125 : CMATE = 'ISOTROPE'
126 : KCAS=1
127 : C Modif L.Bode - 14/10/92
128 : C IFOUR --> IFOUR2 dans appel DOHMAS
129 : * print*,'avant dohmas'
129 : CALL DOHMAS(XMAT,CMATE,IFOUR2,NSTRS1,KCAS,DDHOOK,IRTD)
130 : * print*,'apres dohmas'
131 : C Fin modif L.Bode
132 : DO 200 ISTRS=1,3
133 : SIGP(ISTRS)= XZERO
134 : DO 210 JSTRS=1,3
135 : SIGP(ISTRS)=SIGP(ISTRS)+DDHOOK(ISTRS,JSTRS)*EPSIPP(JSTRS)
136 : 210 CONTINUE
137 : 200 CONTINUE
138 : END IF
139 : C
140 : C on complete la deformation dans le cas des contraintes planes
141 : C
142 : C Modif L.Bode - 14/10/92
143 : C IFOUR remplace par IFOUR2
144 : IF (IFOUR2.EQ. -2) THEN
145 : EPSIPP(3)= -(EPSIPP(1) + EPSIPP(2))*XNU / (UN-XNU)
146 : END IF
147 : C Fin modif L.Bode
148 : C
149 : C on calcule le epsilontild
150 : C
151 : EPSTIL=MAX( XZERO , EPSIPP(1) )**2 +
152 : 1 MAX( XZERO , EPSIPP(2) )**2 +
153 : 2 MAX( XZERO , EPSIPP(3) )**2
154 : EPSTIL=SQRT (EPSTIL)
155 :
156 : IF (ISTEP .EQ. 0) THEN
157 : EPSTIM=EPSTIL
158 : VARF(1)=EPSTIL
159 : ELSE IF (ISTEP .EQ. 1) THEN
160 : VARF(2)=DINI
161 : VARF(1)=EPSTIL
162 : ELSE IF (ISTEP .EQ. 2) THEN
163 : EPSTIM=VAR0(1)
164 : VARF(1)=EPSTIM
165 : ELSE
166 : PRINT*,'DANS MAZARS ISTEP = 0,1,2 ET PAS ',ISTEP
167 : END IF
168 : IF ( (ISTEP .EQ. 0) .OR. (ISTEP .EQ. 2)) THEN
169 : C
170 : C on calcule l'endommagement et les contraintes
171 : C
172 : C
173 : C on calcule ALFAT ALFAC DT et DC puis D
174 : C dans le cas ou le seuil initial est depasse
175 : C
176 : IF ( EPSTIM .GT. EPSD0) THEN
177 : C
178 : C calcul de l'endommagement
179 : C
180 : C
181 : C on calcule le signe des contraintes elastiques
182 : C
183 : DO 300 ISTRS=1,3
184 : IF (SIGP(ISTRS).LT. XZERO) THEN
185 : SIGPC(ISTRS)=SIGP(ISTRS)
186 : SIGPT(ISTRS)=XZERO
187 : ELSE
188 : SIGPT(ISTRS)=SIGP(ISTRS)
189 : SIGPC(ISTRS)=XZERO
190 : END IF
191 : 300 CONTINUE
192 : TRSIGT=SIGPT(1)+SIGPT(2)+SIGPT(3)
193 : TRSIGC=SIGPC(1)+SIGPC(2)+SIGPC(3)
194 : C
195 : C on calcule les deformations dues aux contraintes positives
196 : C
197 : DO 400 ISTRS=1,3
198 : EPSILT(ISTRS)=(SIGPT(ISTRS)*(UN+XNU)-TRSIGT*XNU)/YOUN
199 : 400 CONTINUE
200 : C
201 : C on en deduit ALFAT et ALFAC
202 : C
203 : ALFAT = MAX(XZERO,EPSIPP(1))*EPSILT(1) +
204 : 1 MAX(XZERO,EPSIPP(2))*EPSILT(2) +
205 : 2 MAX(XZERO,EPSIPP(3))*EPSILT(3)
206 : ALFAT = ALFAT/(EPSTIL*EPSTIL)
207 : ALFAC = UN - ALFAT
208 : C
209 : C modification pour la bi ou tricompression
210 : C
211 : IF (TRSIGC.LT. -XPETIT .AND. TRSIGT.LT.XPETIT) THEN
212 : GAMMA=SIGPC(1)*SIGPC(1)+SIGPC(2)*SIGPC(2)+
213 : 1 SIGPC(3)*SIGPC(3)
214 : GAMMA=-SQRT(GAMMA)/TRSIGC
215 : EPSTIM=EPSTIM*GAMMA
216 : END IF
217 : C
218 : C amelioration de la reponse en cisaillement pour beta > 1.
219 : C
220 : IF (BETA .GT. UN) THEN
221 : IF ( ALFAT .GT. XPETIT ) THEN
222 : ALFAT=ALFAT**BETA
223 : END IF
224 : IF ( ALFAC .GT. XPETIT ) THEN
225 : ALFAC=ALFAC**BETA
226 : END IF
227 : END IF
228 :
229 : C
230 : C on calcule DT et DC puis D
231 : C dans le cas ou le seuil initial est depasse
232 : C
233 : C on est oblige de verifier car on a pu multiplier par gamma
234 : C
235 : C IF (EPSTIM .GT. EPSD0) THEN
236 : C DT=UN - EPSD0*(UN-ATRA)/EPSTIM -
237 : C 1 ATRA*EXP(-BTRA*(EPSTIM-EPSD0))
238 : C DC=UN - EPSD0*(UN-ACOM)/EPSTIM -
239 : C 1 ACOM*EXP(-BCOM*(EPSTIM-EPSD0))
240 : C ELSE
241 : C DT=XZERO
242 : C DC=XZERO
243 : C END IF
244 : CLB 2010/01/20
245 : C trois lois d'évolution:
246 : C ATRA > 0 : loi de mazars classique
247 : C -10 BTRA est alors
249 : C la déformation pour laquelle la contrainte s'annule C
250 : C
251 : IF (EPSTIM .GT. EPSD0) THEN
252 : IF (ATRA .GT. 0.D0) THEN
253 : DT=UN - EPSD0*(UN-ATRA)/EPSTIM -
254 : 1 ATRA*EXP(BTRA*(EPSD0-EPSTIM))
255 : ELSE IF (ATRA . GT. -10.) THEN
256 : DT=UN - epsd0/epstim*EXP(BTRA*(EPSD0-EPSTIM))
257 : ELSE
258 : IF(EPSTIM .LT. BTRA) THEN
259 : DT=UN - EPSD0*(BTRA - EPSTIM)/EPSTIM/(BTRA - EPSD0)
260 : ELSE
261 : DT=1.D0
262 : ENDIF
263 : END IF
264 : DC=UN - EPSD0*(UN-ACOM)/EPSTIM -
265 : 1 ACOM*EXP(-BCOM*(EPSTIM-EPSD0))
266 : ELSE
267 :
268 : DT=XZERO
269 : DC=XZERO
270 : END IF
271 : D = ALFAT*DT + ALFAC*DC
272 : C
273 : C on borne la valeur de D a 0.99999..
274 : C
275 : D=MIN ( D , UN-1.D-20 )
276 : ELSE
277 : D=XZERO
278 : END IF
279 : C
280 : C on teste la croissance de D
281 : C
282 : D=MAX ( D , DINI )
283 : C
284 : C on le stocke dans les variables internes finales
285 : C
286 : VARF(2)= D
287 : C
288 : C on calcule les contraintes finales
289 : C
289 : CALL MATVE1 (DDHOOK,EPSILO,NSTRS1,NSTRS1,SIGF,2)
290 : DO 500 ISTRS=1,NSTRS1
291 : SIGF(ISTRS)=SIGF(ISTRS)*(UN-D)
292 : 500 CONTINUE
293 : C
294 : C et les deformations inelastiques finales
295 : C
296 : DO 600 ISTRS=1,NSTRS1
297 : EPINF(ISTRS)=EPSILO(ISTRS)*D
298 : 600 CONTINUE
299 :
300 : END IF
301 : SEGSUP WRK3
302 : RETURN
303 : END
304 :
305 :
306 :
307 :
308 :
309 :