1 : C PROJVM SOURCE CHAT 05/01/13 02:32:35 5004
2 : SUBROUTINE PROJVM(SIG,EPS,EPST,EPSTAR,DEPS,PENTE,SN,
3 : . PANTIN,SELAS,ITER,SO,SSTAR,SI,NK,IBOU,COEF,JKI,
4 : . KERRE,ecou,necou)
5 : IMPLICIT INTEGER(I-N)
6 : IMPLICIT REAL*8(A-H,O-Z)
7 : C
8 : C CALCULE L ECOULEMENT EN CAS DE CRITERE DE VON MISES
9 : C
10 : DIMENSION SIG(*),EPS(*),COEF(*)
11 : C
12 : SEGMENT ECOU
13 : *** COMMON/ECOU/TEST,ALFAH,
14 : REAL*8 TEST, ALFAH,
15 : 1 HPAS,TEMPS,COVNMS(6),VECPRO(9),VALPRO(6),
16 : 2 CVNMSD(12),STOT(6),SIGEL(6),DSIGP(6),SIGT(6),W1(6),W2(6),
17 : 1 DALPHA(6),DSIGO(6),E(12),XINV(3),
18 : 2 SIPLAD(6),DSIGP0(6),TETM,TETI
19 : ENDSEGMENT
20 : C COMMON/ECOU/TEST,ALFAH,
21 : C 1 HPAS,TEMPS,COVNMS(6),VECPRO(9),VALPRO(6),
22 : C 2 CVNMSD(12),STOT(6),SIGEL(6),DSIGP(6),SIGT(6),W1(6),W2(6),
23 : C 1 DALPHA(6),DSIGO(6),E(12),XINV(3),
24 : C 2 SIPLAD(6),DSIGP0(6),TETM,TETI
25 : C
26 : SEGMENT NECOU
27 : * COMMON/NECOU/NCOURB,IPLAST,IT,IMAPLA,ISOTRO,
28 : INTEGER NCOURB,IPLAST,IT,IMAPLA,ISOTRO,
29 : . ITYP,IFOUR,IFLUAG,
30 : . ICINE,ITHER,IFLUPL,ICYCL,IBI,
31 : . JFLUAG,KFLUAG,LFLUAG,
32 : . IRELAX,JNTRIN,MFLUAG,JSOUFL,JGRDEF
33 : ENDSEGMENT
34 : C COMMON/NECOU/NCOURB,IPLAST,IT,IMAPLA,ISOTRO,
35 : C . ITYP,IFOUR,IFLUAG,
36 : C . ICINE,ITHER,IFLUPL,ICYCL,IBI,
37 : C . JFLUAG,KFLUAG,LFLUAG,
38 : C . IRELAX,JNTRIN,MFLUAG,JSOUFL,JGRDEF
39 : C
40 : DATA ITMAX/200/
41 : C
42 : JEBOUC=0
43 : PENTEJ=PENTE
44 : DEPSJ=DEPS
45 : 2020 JEBOUC=JEBOUC+1
46 : PENTE=PENTEJ
47 : DEPS=DEPSJ
48 : IF(JEBOUC.GT.1) WRITE(6,77881) PENTE,DEPS
49 : 77881 FORMAT('0 VALEURS INITIALES PENTE=',1PE12.5,2X,
50 : . 'DEPS=',1PE12.5/)
51 : C
52 : C
53 : C CAS PLASTIQUE
54 : PENTE1=PENTE
55 : EPST=EPSTAR+DEPS
56 : C CONTRAINTE SUR LA COURBE DE TRACTION CORRESPONDANTE
57 : IF(ICINE.EQ.1) GO TO 81
58 : C CAS ISOTROPE
58 : CALL TRACTI(SN,EPST,SIG,EPS,NCOURB,2,IBI)
59 : IF(IBI.EQ.1) THEN
60 : KERRE=75
61 : GO TO 57
62 : ENDIF
63 : *
64 : C CAS MIXTE CALCUL DU NOUVEAU RAYON DE LA SPHERE
65 : IF(ICINE.EQ.2) SN=SN-PANTIN*EPST
66 : GO TO 82
67 : C CAS CINEMATIQUE
68 : 81 SN=SELAS
69 : 82 CONTINUE
70 : C CONTRAINTE MOYENNE
71 : SM=0.5*(SN+SELAS)
72 : 2222 CONTINUE
73 : C
74 : C ITERATIONS INTERNES CALCUL DE DELTA EPSILON
75 : C
76 : ITER=0
77 : DEPS0=DEPS
78 : DEPSI=DEPS0
79 : SO=SN
80 : DEPS00=DEPS0
81 : JYCOMP=0
82 : UNSSM=1./SM
83 : SIGT(1)=SIGEL(1)+DSIGP(1)
84 : SI=SSTAR
85 : IF(JEBOUC.GT.1) WRITE(6,77882) SN,SM,SI
86 : 77882 FORMAT('0 SN =',1PE12.5,2X,'SM=',1PE12.5,2X,
87 : . 'SI=',1PE12.5/)
88 : 555 ITER=ITER+1
89 : IF(JEBOUC.GT.1) WRITE(6,77883) ITER
90 : 77883 FORMAT('0 -----> ITER=',I3/)
91 : DEPSI1=DEPSI
92 : SI1=SI
93 : STEST=TEST*SN
94 : Z=0.
95 : DO 52 IB=NK,IBOU
96 : Y=UNSSM*E(IB)
97 : X=Y*DEPS
98 : IF(X.EQ.0.) GO TO 51
99 : EXPMX=EXP(-X)
100 : UNSX=1./X
101 : UNSX2=UNSX*UNSX
102 : X2=X*X
103 : IF (ABS(X).LT.0.01) GO TO 50
104 : SIGT(IB)=(SIGEL(IB)-DSIGP(IB)*UNSX)*EXPMX+DSIGP(IB)*UNSX
105 : Z=Z+COEF(JKI+IB)*SIGT(IB)*(EXPMX*(DSIGP(IB)*UNSX2-SIGEL(IB)+
106 : 1DSIGP(IB)*UNSX)-DSIGP(IB)*UNSX2)*E(IB)
107 : GO TO 52
108 : 50 CONTINUE
109 : SIGT(IB)=SIGEL(IB)*(1.-X+0.5*X2)+DSIGP(IB)*(1-0.5*X)
110 : Z=Z+COEF(JKI+IB)*SIGT(IB)*E(IB)*(DSIGP(IB)*(X*.333333-0.5)-
111 : 1 SIGEL(IB)*(1-X+0.5*X2))
112 : GO TO 52
113 : C UNE VALEUR PROPRE EST NULLE
114 : 51 SIGT(IB)=SIGEL(IB)+DSIGP(IB)
115 : 52 CONTINUE
116 : C NOUVEL DEPS
117 : SI=VNMISD(SIGT,ITYP,ALFAH,CVNMSD)
118 : X=PENTE-Z*UNSSM/SI
119 : ZOB1 = UNSSM/SI
120 : ZOB2 = ZOB1 * Z
121 : IF(JEBOUC.GT.1) WRITE(6,77888) PENTE,Z,UNSSM
122 : 77888 FORMAT('0 PENTE =',1PE12.5,2X,'Z=',1PE12.5,2X,'UNSSM=',1PE12.5/)
123 : DEPSI=(SI-SN)/X
124 : DEPS=DEPS+DEPSI
125 : EPST=DEPS+EPSTAR
126 : IF(JEBOUC.GT.1) WRITE(6,77884) SI,SN,X,DEPSI,DEPS
127 : 77884 FORMAT('0 SI =',1PE12.5,2X,'SN=',1PE12.5,2X,'X=',1PE12.5,2X,
128 : . 'DEPSI=',1PE12.5,2X,'DEPS=',1PE12.5/)
129 : IF(ICINE.EQ.1) GO TO 580
130 : C ON CALCULE LA CONTRAINTE SUR LA COURBE DE TRACTION A
131 : C LA FIN DU PAS
132 : EPSTP=EPST
133 : C ON CONTROLE SI LE DEPS EST POSITIF
134 : C SI DEPS EST NEGATIF ON PREND EPSTP =EPST A L ITERATION PRECE
135 : IF(DEPS.LT.0.) EPSTP=EPSTP-DEPSI
135 : CALL TRACTI(SN,EPSTP,SIG,EPS,NCOURB,2,IBI)
136 : IF(IBI.EQ.1) THEN
137 : KERRE=75
138 : GO TO 57
139 : ENDIF
140 : C ON CALCULE LA PENTE A LA COURBE DE TRACTION A LA FIN DU PAS
140 : CALL TRACTI(PENTE1,EPSTP,SIG,EPS,NCOURB,1,IBI)
141 : IF(ICINE.EQ.2) SN=SN-EPST*PANTIN
142 : 580 CONTINUE
143 : IF(ITER.GE.ITMAX) GO TO 56
144 : C ON N A PAS DEPASSE LE NOMBRE MAX D ITERATIONS INTERNES
145 : IF(DEPS) 54,500,500
146 : C DEPS EST POSITIF
147 : 500 CONTINUE
148 : PENTE=PENTE1
149 : C ON PREND COMME NOUVELLE PENTE LA PENTE A LA FIN DU PAS
150 : IF(ABS(SI-SN)-STEST) 57,57,555
151 : C A T ON CONVERGE
152 : 56 SSTST=ABS(SI-SN)/SN
153 : KERRE=2
154 : write(6,fmt='('' kerre projvm'',i6)')kerre
155 : GO TO 57
156 : 54 CONTINUE
157 : C DEPS EST NEGATIF
158 : C LA PENTE A LA FIN DU PAS EST-ELLE SUPERIEURE A LA PENTE
159 : C AU DEBUT DANS LES CAS ISOTROPE ?
160 : IF(ICINE.EQ.1) GO TO 543
161 : *
162 : * MILL 26/3/93
163 : IF(PENTE.EQ.0.D0.AND.DEPS.LE.0.D0) GOTO 540
164 : *
165 : X=ABS(PENTE)
166 : IF(X.LT.ABS(PENTE1)) X=ABS(PENTE1)
167 : X=ABS((PENTE-PENTE1)/(X+1.D-10))
168 : IF(X.LT..2D00) GO TO 543
169 : C LA PENTE VARIE DE PLUS DE 20%
170 : C ON REGARDE SI LA COURBE EST CONCAVE
171 : IF(PENTE1.LT.PENTE) GO TO 543
172 : EPST=EPST-DEPS
173 : DEPS=DEPS-DEPSI
174 : C ON PREND PENTE =.5*(PENTE AU DEBUT + PENTE A LA FIN )
175 : PENTE=0.5D00*(PENTE+PENTE1)
176 : C ON CALCULE LE NOUVEAU DEPSI
177 : X=PENTE-Z*UNSSM/SI
178 : DEPSI=(SI-SN)/X
179 : DEPS=DEPS+DEPSI
180 : IF(DEPS.GT.0.D00) GO TO 542
181 : C DEPS EST ENCORE NEGATIF
182 : GO TO 540
183 : 543 CONTINUE
184 : C LE POINT SI EST AU DESSUS DE LA SURFACE
185 : C ON RECOMMENCE AVEC LE PREMIER DEPS DIVISE PAR 2
186 : C LE POINT SI EST EN DESSOUS ON ASSIMILE LA TGT A LA CORDE
187 : EPST=EPST-DEPS
188 : DEPS=DEPS-DEPSI
189 : IF(SI.GT.SN) GO TO 540
190 : X=PENTE-(SI-SI1)/DEPSI1
191 : GO TO 541
192 : 540 CONTINUE
193 : IF(JYCOMP.GT.0) DEPS00=0.5*DEPS00
194 : JYCOMP=JYCOMP+1
195 : IF(JEBOUC.GT.1) WRITE(6,77887) JYCOMP
196 : 77887 FORMAT('0 JYCOMP = ',I4/)
197 : DEPS=0.5*DEPS00
198 : EPSOM=EPSTAR+DEPS
199 : DEPSI=DEPS
199 : CALL TRACTI(SN,EPSOM,SIG,EPS,NCOURB,2,IBI)
200 : IF(IBI.EQ.1) THEN
201 : KERRE=75
202 : GO TO 57
203 : ENDIF
204 : SM=0.5D00*(SN+SELAS)
205 : UNSSM=1.D00/SM
206 : GO TO 542
207 : 541 CONTINUE
208 : DEPSI=(SI-SN)/X
209 : DEPS=DEPS+DEPSI
210 : IF(DEPS.LE.0.D00) GO TO 540
211 : 542 CONTINUE
212 : EPST=EPSTAR+DEPS
213 : GO TO 555
214 : 57 CONTINUE
215 : IF( DEPS.LT.0.D0) THEN
216 : IF(JEBOUC.EQ.1) GO TO 2020
217 : ENDIF
218 : RETURN
219 : END
220 :
221 :
222 :
223 :
224 :
© Cast3M 2003 - All rights reserved.
Disclaimer