Télécharger ccoin0.eso

Retour à la liste

Numérotation des lignes :

ccoin0
  1. C CCOIN0 SOURCE SP204843 23/06/14 21:15:04 11675
  2.  
  3. SUBROUTINE CCOIN0(wrk52,wrk53,wrk54,wrk2,wrk3,
  4. & IB,IGAU,NBPGAU,LTRAC,IFOURB,iecou)
  5.  
  6. C-----------------------------------------------------------------------
  7. C ECOULEMENT PLASTIQUE POUR UN POINT
  8. C ALGORITHME ORTIZ ET SIMO
  9. C
  10. C EN ENTREE :
  11. C
  12. C SIG0 CONTRAINTES AU DEBUT DU PAS
  13. C DEPST INCREMENT DE DEFORMATIONS TOTALES
  14. C ( THERMIQUE ENLEVEE )
  15. C VAR0 VARIABLES INTERNES DEDUT DU PAS
  16. C VAREX0 VARIABLES EXTERNES DEBUT DU PAS
  17. C VAREXF VARIABLES EXTERNES FIN DU PAS
  18. C XMAT COEFFICIENTS DU MATERIAU
  19. C PRECIS PRECISION POUR ITERATIONS INTERNES
  20. C WORK DES CARACTERISTIQUES
  21. C TRAC COURBE DE TRACTION
  22. C MFR1 INDICE DE FORMULATION
  23. C NSTRSS NOMBRE DE CONTRAINTES CA2000
  24. C INPLAS NUMERO DU MODELE DE PLASTICITE
  25. C DDAUX = MATRICE DE HOOKE ELASTIQUE
  26. C CMATE = NOM DU MATERIAU
  27. C VALMAT= TABLEAU DE CARACTERISTIQUES DU MATERIAU
  28. C VALCAR= TABLEAU DE CARACTERISTIQUES GEOMETRIQUES
  29. C N2EL = NBRE D ELEMENTS DANS SEGMENT DE HOOKE
  30. C N2PTEL= NBRE DE POINTS DANS SEGMENT DE HOOKE
  31. C IFOU = OPTION DE CALCUL
  32. C IB = NUMERO DE L ELEMENT COURANT
  33. C IGAU = NUMERO DU POINT COURANT
  34. C EPAIST= EPAISSEUR
  35. C NBPGAU= NBRE DE POINTS DE GAUSS
  36. C MELE = NUMERO DE L ELEMENT FINI
  37. C NPINT = NBRE DE POINTS D INTEGRATION
  38. C NBGMAT= NBRE DE POINTS DANS SEGMENT DE CARACTERISTIQUES
  39. C NELMAT= NBRE D ELEMENTS DANS SEGMENT DE CARACTERISTIQUES
  40. C SECT = SECTION
  41. C LHOOK = TAILLE DE LA MATRICE DE HOOKE
  42. C TXR,XLOC,XGLOB,D1HOOK,ROTHOO,DDHOMU,CRIGI = TABLEAUX UTILISES
  43. C UTILISES POUR LE CALCUL DE LA MATRICE DE HOOKE
  44. C
  45. C EN SORTIE :
  46. C
  47. C SIGF CONTRAINTES A LA FIN DU PAS
  48. C VARF VARIABLES INTERNES FIN DU PAS
  49. C DEFP INCR. DE DEFORMATIONS PLASTIQUES
  50. C KERRE CODE D'ERREUR
  51. C = 0 SI TOUT OK
  52. C = 99 SI FORMULATION NON DISPONIBLE
  53. C EN PLASTICITE
  54. C = 1 SI DEPS EST NEGATIF
  55. C = 2 SI NOMBRE MAX D'ITERATIONS INTERNES DEPASSE
  56. C
  57. C IFOUR : OPTION DE CALCUL
  58. C
  59. C IFOUR = -3 DEFORMATION PLANE GENERALISEE
  60. C IFOUR = -2 CONTRAINTES PLANES
  61. C IFOUR = -1 DEFORMATIONS PLANES
  62. C IFOUR = 0 AXISYMETRIQUE
  63. C IFOUR = 1 SERIE DE FOURIER
  64. C IFOUR = 2 TRIDIMENSIONNEL
  65. C
  66. C MFR1 : NUMERO DE LA FORMULATION ELEMENT FINI
  67. C
  68. C MFR1 = 1 MASSIF
  69. C MFR1 = 3 COQUE
  70. C MFR1 = 5 COQUE EPAISSE
  71. C MFR1 = 7 POUTRE
  72. C MFR1 = 9 COQUE AVEC CISAILLEMENT TRANSVERSE
  73. Ckich MFR1 = 31 pondération réduite termes diagonaux matrice B,
  74. C dite MASSIF INCOMPRESSIBLE. Utilisation en contraintes planes a justifier
  75. c doublon historique MFR/MFR1
  76. C
  77. C INPLAS : NUMERO DU MODELE DE PLASTICITE
  78. C
  79. C INPLAS = 1 PARFAIT
  80. C INPLAS = 4 CINEMATIQUE
  81. C INPLAS = 5 ISOTROPE
  82. C INPLAS = 7 CHABOCHE1
  83. C INPLAS = 12 CHABOCHE2
  84. C
  85. C-----------------------------------------------------------------------
  86. C CONVENTION DE REMPLISSAGE DES MEMOIRES
  87. C-----------------------------------------------------------------------
  88. C
  89. C XMAT(1) MODULE D'YOUNG
  90. C XMAT(2) COEFFICIENT DE POISSON
  91. C
  92. C TRAC(1 A 2*LTRAC) COURBE DE TRACTION
  93. C WORK( " +1) ALFAH POUR COQUES PLASTICITE GLOBALE
  94. C WORK( " +..) DONNEES POUR CRITERE POUTRES GLOBALES
  95. C
  96. C MODELE ISOTROPE
  97. C VAR0(1) EPS*
  98. C
  99. C MODELE CINEMATIQUE LINEAIRE
  100. C VAR0(1) EPS*
  101. C VAR0(2) A VAR0(1+NSTRSS) DEPLACEMENT DE LA SPHERE
  102. C
  103. C MODELE CHABOCHE
  104. C XMAT(5) .... COEFFICIENTS A,C,...
  105. C VAR0(1) EPS*
  106. C VAR0(2) A VAR0(1+NSTRSS) DEPLACEMENT DE LA SPHERE 1
  107. C VAR0(2+NSTRSS) A VAR0(1+2*NSTRSS) " " " " 2
  108. C
  109. C-----------------------------------------------------------------------
  110. C 20/09/2017 : modif SG critere de convergence trop serre
  111. C TEST=PETI*APHI0 + utilisation CCREEL
  112. C voir aussi ecoin0.eso, syco12.eso
  113.  
  114. IMPLICIT INTEGER(I-N)
  115. IMPLICIT REAL*8(A-H,O-Z)
  116.  
  117. -INC PPARAM
  118. -INC CCOPTIO
  119. -INC CCREEL
  120. -INC DECHE
  121.  
  122. SEGMENT IECOU
  123. INTEGER icow1,icow2,icow3,icow4,icow5,icow6,icow7,
  124. 1 icow8,icow9,icow10,icow11,icow12,icow13,icow14,icow15,icow16,
  125. 2 icow17,icow18,icow19,icow20,icow21,icow22,icow23,icow24,
  126. 3 icow25,icow26,icow27,icow28,icow29,icow30,icow31,
  127. 4 icow32,icow33,NSTRSS,MFR1, NBGMAT,NELMAT,icow38,
  128. 5 icow39,icow40,icow41,icow42,icow43,icow44
  129. INTEGER icow45,icow46,icow47,icow48,icow49,icow50,
  130. . icow51,icow52,icow53,icow54,icow55,icow56
  131. . icow57,icow58
  132. ENDSEGMENT
  133. *
  134. SEGMENT WRK2
  135. REAL*8 TRAC(LTRAC)
  136. ENDSEGMENT
  137. *
  138. SEGMENT WRK3
  139. REAL*8 WORK(LW),WORK2(LW2)
  140. ENDSEGMENT
  141.  
  142. DIMENSION SIG(130),EPS(130)
  143. DIMENSION S(8),SX(8),DS(8),DSIG(8),SPHER(8),SPHER1(8),SPHER2(8)
  144. DIMENSION DSPHER1(8),DSPHER2(8),DEPSE(8),DEPSP(8),DDEPSE(8)
  145. DIMENSION F(8),W1(8),W2(8),SIGB(8),Z1(8),DIV(8),DDA(8,8)
  146. DIMENSION CRIGI(12)
  147.  
  148. * logical lvisu
  149. * lvisu = ib.eq.125.and.(igau.eq.1.or.igau.eq.5)
  150.  
  151. ncara = commat(/2)
  152. if(ib.eq.1.and.igau.eq.1) then
  153. do iaca = 1,ncara
  154. if(commat(iaca).eq.'LIMP') icow21 = iaca
  155. enddo
  156. endif
  157. if(icow21.gt.0) xlimp = valma0(icow21)
  158. do jj = 1,8
  159. sx(jj) = 0.d0
  160. enddo
  161. C---------COQUES AVEC CT------------------------------------------------
  162. C ON NE TRAVAILLE QUE SUR LES 6 PREMIERES COMPOSANTES
  163. IF(MFR1.EQ.9) THEN
  164. NCONT=6
  165. ELSE
  166. NCONT=NSTRSS
  167. ENDIF
  168. itracb=0
  169.  
  170. C-----------------------------------------------------------------------
  171. 2222 continue
  172. KERRE=0
  173.  
  174. DO I=1,NSTRSS
  175. S(I)=0.D0
  176. SPHER(I)=0.D0
  177. SPHER1(I)=0.D0
  178. SPHER2(I)=0.D0
  179. ENDDO
  180. YUNG=XMAT(1)
  181. XNU=XMAT(2)
  182.  
  183. C---------CARACTERISTIQUES GEOMETRIQUES---------------------------------
  184. C
  185. C COQUES
  186. C
  187. ALFAH=1.D0
  188. EP1 =1.D0
  189. IF(MFR1.EQ.3.OR.MFR1.EQ.5.OR.MFR1.EQ.9) THEN
  190. EP1=WORK(1)
  191. IF(MFR1.NE.5) ALFAH=WORK(2)**2
  192. ENDIF
  193.  
  194. C---------COEFFICIENTS CHABOCHE-----------------------------------------
  195.  
  196. IF(INPLAS.EQ.7) THEN
  197. A1=XMAT(5)
  198. C1=XMAT(6)
  199. R0=XMAT(7)
  200. PSI=XMAT(8)
  201. OME=XMAT(9)
  202. RM=XMAT(10)
  203. B=XMAT(11)
  204. A2=0.D0
  205. C2=0.D0
  206. ELSE IF(INPLAS.EQ.12) THEN
  207. A1=XMAT(5)
  208. C1=XMAT(6)
  209. A2=XMAT(7)
  210. C2=XMAT(8)
  211. R0=XMAT(9)
  212. PSI=XMAT(10)
  213. OME=XMAT(11)
  214. RM=XMAT(12)
  215. B=XMAT(13)
  216. ELSE IF(INPLAS.EQ.4) THEN
  217. SLI0 = XMAT(5)
  218. HECL = XMAT(6)
  219. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  220. c*dbg write(ioimp,*) 'Ecrouissage',INPLAS,SLI0,HECL
  221. c*dbg endif
  222. ELSE
  223. DO I=1,LTRAC
  224. SIG(I)=TRAC(2*I-1)
  225. EPS(I)=TRAC(2*I)
  226. ENDDO
  227. ENDIF
  228.  
  229. EPST=VAR0(1)
  230. C---------ECROUISSAGE CINEMATIQUE---------------------------------------
  231.  
  232. IF(INPLAS.EQ.4.OR.INPLAS.EQ.7.OR.INPLAS.EQ.12) THEN
  233. DO I=1,NSTRSS
  234. SPHER1(I)=VAR0(I+1)
  235. ENDDO
  236. IF(INPLAS.EQ.12) THEN
  237. DO I=1,NSTRSS
  238. SPHER2(I)=VAR0(NSTRSS+1+I)
  239. ENDDO
  240. ENDIF
  241. DO I=1,NSTRSS
  242. SPHER(I)=SPHER1(I)+SPHER2(I)
  243. ENDDO
  244. c*dbg if (ib.eq.1 .and.igau.eq.1) then
  245. c*dbg write(ioimp,*) ' SPHERi',NSTRSS
  246. c*dbg write(ioimp,*) ' ',(spher1(iof),iof=1,NSTRSS)
  247. c*dbg write(ioimp,*) ' ',(spher(iof),iof=1,NSTRSS)
  248. c*dbg endif
  249. C-----------------------------------------------------------------------
  250. C PREDICTEUR ELASTIQUE
  251. C S : CONTRAINTE
  252. C SPHER : VARIABLE D'ECROUISSAGE CINEMATIQUE
  253. C SX = S - X
  254. C-----------------------------------------------------------------------
  255.  
  256. * en elastique non lineaire on annule les contraintes initiales
  257. * mais on cumule les epsilons elastiques
  258. ELSE IF(INPLAS.EQ.87) THEN
  259. EPST=0.D0
  260. DO I=1,NSTRSS
  261. SIG0(I)=0.D0
  262. DEPST(I)=DEPST(I) + VAR0(I+1)
  263. ENDDO
  264. ENDIF
  265. CALL CALSIG(DEPST,DDAUX,NSTRSS,CMATE,VALMAT,VALCAR,
  266. & N2EL,N2PTEL,MFR1,IFOURB,IB,IGAU,EPAIST,
  267. & NBPGAU,MELE,NPINT,NBGMAT,NELMAT,SECT,LHOOK,TXR,
  268. & XLOC,XGLOB,D1HOOK,ROTHOO,DDHOMU,CRIGI,DSIGT,IRTD)
  269.  
  270. IF(IRTD.NE.1) THEN
  271. KERRE=69
  272. GOTO 1000
  273. ENDIF
  274.  
  275. IF ((mfr.eq.1.or.mfr.eq.31).and. ifourb.eq.-2) then
  276. * en cas de contraintes planes on repasse en 3D
  277. do ju=1,6
  278. do iu=1,6
  279. DDA(iu,ju)=0.d0
  280. enddo
  281. enddo
  282. cte_cp = xnu / (1.d0 - xnu)
  283. aux= YUNG / ((1.d0+xnu)*(1.d0-2.d0*xnu))
  284. aux1= aux * xnu
  285. aux2= aux * (1.d0-xnu)
  286. gege = Yung / (2.d0 *(1.d0 +xnu))
  287. DDA(1,1)=AUX2
  288. DDA(2,1)=AUX1
  289. DDA(1,2)=AUX1
  290. DDA(2,2)=AUX2
  291. DDA(3,3)=aux2
  292. DDA(1,3)=aux1
  293. DDA(2,3)=aux1
  294. DDA(3,1)=aux1
  295. DDA(3,2)=aux1
  296. DDA(4,4)=gege
  297. DDA(5,5)=gege
  298. DDA(6,6)=GEGE
  299.  
  300. ELSE IF ((MFR.EQ.3.AND.NPINT.EQ.0).OR.MFR.EQ.9) THEN
  301. AUX=YUNG/(1.D0-XNU*XNU)
  302. AUX1=AUX*XNU
  303. DO J=1,NCONT
  304. DO I=1,NCONT
  305. DDAUX(I,J)=0.D0
  306. ENDDO
  307. ENDDO
  308. C
  309. C CAS TRIDIMENSIONNEL ET FOURIER
  310. C
  311. IF(IFOURB.EQ.2.OR.IFOURB.EQ.1) THEN
  312.  
  313. GEGE=0.5D0*YUNG/(1.D0+XNU)
  314. DDAUX(1,1)=AUX
  315. DDAUX(2,1)=AUX1
  316. DDAUX(1,2)=AUX1
  317. DDAUX(2,2)=AUX
  318. DDAUX(3,3)=GEGE
  319. DDAUX(4,4)=AUX
  320. DDAUX(5,4)=AUX1
  321. DDAUX(4,5)=AUX1
  322. DDAUX(5,5)=AUX
  323. DDAUX(6,6)=GEGE
  324. C
  325. C CAS AXISYMETRIQUE ET DEFORMATIONS PLANES
  326. C
  327. ELSE IF(IFOURB.EQ.0.OR.IFOURB.EQ.-1.OR.IFOURB.EQ.-3) THEN
  328.  
  329. DDAUX(1,1)=AUX
  330. DDAUX(2,1)=AUX1
  331. DDAUX(1,2)=AUX1
  332. DDAUX(2,2)=AUX
  333. DDAUX(3,3)=AUX
  334. DDAUX(4,3)=AUX1
  335. DDAUX(3,4)=AUX1
  336. DDAUX(4,4)=AUX
  337. C
  338. C CAS CONTRAINTES PLANES
  339. C
  340. ELSE IF(IFOURB.EQ.-2) THEN
  341. DDAUX(1,1)=YUNG
  342. DDAUX(3,3)=YUNG
  343. ENDIF
  344. ENDIF
  345. *
  346. DO I=1,NSTRSS
  347. S(I)=SIG0(I)+DSIGT(I)
  348. SIGB(I)=S(I)
  349. SX(I)=S(I)-SPHER(I)
  350. ENDDO
  351. C---------CAS DES POUTRES-----------------------------------------------
  352.  
  353. IF(MFR1.EQ.7) THEN
  354. DIV(1)=1.D0/WORK(4)
  355. DIV(2)=1.D0
  356. DIV(3)=1.D0
  357. DIV(4)=WORK(10)/WORK(1)
  358. DIV(5)=WORK(11)/WORK(2)
  359. DIV(6)=WORK(12)/WORK(3)
  360. IF(DIV(4).EQ.0.D0) DIV(4)=1.D-10/SQRT(WORK(1)*WORK(4))
  361. IF(DIV(5).EQ.0.D0) DIV(5)=1.D-10/SQRT(WORK(2)*WORK(4))
  362. IF(DIV(6).EQ.0.D0) DIV(6)=1.D-10/SQRT(WORK(3)*WORK(4))
  363. DO I=1,NCONT
  364. S(I) = S(I) *DIV(I)
  365. SX(I) = SX(I)*DIV(I)
  366. SPHER(I) = SPHER(I) *DIV(I)
  367. SPHER1(I) = SPHER1(I)*DIV(I)
  368. SPHER2(I) = SPHER2(I)*DIV(I)
  369. END DO
  370. ENDIF
  371. c*dbg if (ib.eq.1 .AND.igau.eq.1) then
  372. c*dbg write(ioimp,*) ' SIGi',NSTRSS
  373. c*dbg write(ioimp,*) ' S',(s(iof),iof=1,NSTRSS)
  374. c*dbg write(ioimp,*) ' SX',(sx(iof),iof=1,NSTRSS)
  375. c*dbg endif
  376.  
  377. C-----------------------------------------------------------------------
  378. C CALCUL DE LA LIMITE ELASTIQUE SI
  379. C-----------------------------------------------------------------------
  380. IF (INPLAS.EQ.1) THEN
  381. SI=TRAC(1)
  382. C* Modele a ECROUISSAGE CINEMATIQUE LINEAIRE : la limite d'elasticite
  383. C* est la limite initiale donnee par SIGY
  384. ELSE IF (INPLAS.EQ.4) THEN
  385. SI = SLI0
  386. ELSE IF (INPLAS.EQ.5.OR.INPLAS.EQ.87) THEN
  387. CALL TRACTI(SI,EPST,SIG,EPS,LTRAC,2,IBI)
  388. IF(IBI.EQ.1) THEN
  389. KERRE=75
  390. GOTO 1000
  391. ENDIF
  392. C* Modele de CHABOCHE (prise en compte ecrouissage)
  393. ELSE IF (INPLAS.EQ.7 .OR. INPLAS.EQ.12) THEN
  394. RMmRR = (RM - R0) * EXP(-B*EPST)
  395. SI = RM - RMmRR
  396. ENDIF
  397. c*dbg if (ib.eq.1 .AND.igau.eq.1) then
  398. c*dbg write(ioimp,*) 'Limite d elasticite SI=',SI,TRAC(1),EPST
  399. c*dbg endif
  400.  
  401. **
  402. * kich : application pression limite trace sigma
  403. **
  404. yxsxii = 0.d0
  405. if (icow21.gt.0) then
  406. ytr = sx(1) + sx(2) + sx(3)
  407. ytr = ytr/3.D0
  408. if (ytr.gt.xlimp) then
  409. yxsxii = ytr - xlimp
  410. else if((ytr + xlimp).lt.0) then
  411. yxsxii = ytr + xlimp
  412. endif
  413. if (ysxii.ne.0) then
  414. dsigii = dsigt(1) + dsigt(2) + dsigt(3)
  415. if (dsigii.ne.0.D0) then
  416. do jj = 1,3
  417. s(jj) = s(jj) - (dsigt(jj)/dsigii*3.D0*yxsxii)
  418. sx(jj) = sx(jj) - (dsigt(jj)/dsigii*3.D0*yxsxii)
  419. enddo
  420. else
  421. do jj = 1,3
  422. s(jj) = s(jj) - yxsxii
  423. sx(jj) = sx(jj) - yxsxii
  424. enddo
  425. endif
  426. endif
  427. endif
  428.  
  429. C-----------------------------------------------------------------------
  430. C CALCUL DE LA CONTRAINTE EQUIVALENTE SEQ
  431. C-----------------------------------------------------------------------
  432. * attention en contraintes planes on se declare en 3D
  433. * rien besoin de faire dans vonmis0 car ifourb n'est pas utilise
  434. * et les contraintes sont dimensionnees a 6
  435. SEQ=VONMIS0(SX,NSTRSS,MFR1,IFOURB,EP1,ALFAH)
  436. c*dbg if (ib.eq.1 .AND.igau.eq.1) then
  437. c*dbg write(ioimp,*) 'VONMISES',SEQ,MFR1,IFOURB,EP1,ALFAH
  438. c*dbg endif
  439.  
  440. C-----------------------------------------------------------------------
  441. C LE CRITERE EST-IL VERIFIE
  442. C-----------------------------------------------------------------------
  443. IF (SEQ.LE.SI) THEN
  444. ITRY = 1
  445. ELSE
  446. PETI=1.1D0*0.5D0*PRECIS*SEQ
  447. CALL EPSPRE(SEQ,SI,PETI,ITRY)
  448. ENDIF
  449. IF (ITRY.EQ.1) THEN
  450. * rien a faire on n'a pas plastifie
  451. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  452. c*dbg write(ioimp,*) 'pas de plastification'
  453. c*dbg endif
  454. IF (MFR1.EQ.7) THEN
  455. DO I=1,NCONT
  456. S(I)=S(I)/DIV(I)
  457. ENDDO
  458. ENDIF
  459. DO I=1,NCONT
  460. SIGF(I)=S(I)
  461. DEFP(I)=0.D0
  462. ENDDO
  463. IF(MFR1.EQ.9) THEN
  464. DEFP(7)=0.D0
  465. DEFP(8)=0.D0
  466. SIGF(7)=S(7)
  467. SIGF(8)=S(8)
  468. ENDIF
  469. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  470. c*dbg write(ioimp,*) 'SIGF',(SIGF(iof),iof=1,NSTRSS)
  471. c*dbg endif
  472.  
  473. VARF(1)=VAR0(1)
  474. IF(INPLAS.EQ.4.OR.INPLAS.EQ.7) THEN
  475. DO I=1,NSTRSS
  476. VARF(I+1)=VAR0(I+1)
  477. ENDDO
  478. ELSE IF(INPLAS.EQ.12) THEN
  479. DO I=1,2*NSTRSS
  480. VARF(I+1)=VAR0(I+1)
  481. ENDDO
  482. ELSE IF(INPLAS.EQ.87) THEN
  483. DO I=1,NSTRSS
  484. VARF(I+1)=DEPST(I)
  485. ENDDO
  486. ENDIF
  487. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  488. c*dbg write(ioimp,*) 'VARF',EPST,'---',(VARF(iof),iof=2,1+NSTRSS)
  489. c*dbg endif
  490. RETURN
  491. ENDIF
  492.  
  493. C-----------------------------------------------------------------------
  494. C ON A PLASTIFIE
  495. C-----------------------------------------------------------------------
  496. NITER=0
  497. PHI=SEQ-SI
  498. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  499. c*dbg write(ioimp,*) 'niter,phi,si,seq,peti,precis=',
  500. c*dbg $ niter,phi,si,seq,peti,precis
  501. c*dbg endif
  502. cri0= si * 1.D-8
  503. PHI0=PHI
  504. SI0=SI
  505. RR=0.D0
  506.  
  507. DO I=1,NCONT
  508. DSIG(I) =0.D0
  509. DEPSP(I) =0.D0
  510. DSPHER1(I)=0.D0
  511. DSPHER2(I)=0.D0
  512. ENDDO
  513.  
  514. C-----------------------------------------------------------------------
  515. C DEBUT DE LA BOUCLE D'ITERATIONS INTERNES
  516. C-----------------------------------------------------------------------
  517. sx1in=sx(1)
  518. sx2in=sx(2)
  519. sx3in=sx(3)
  520. s1in=s(1)
  521. s2in=s(2)
  522.  
  523. iderin=0
  524. 10 CONTINUE
  525. NITER=NITER+1
  526.  
  527. phi= seq - si
  528.  
  529. C---------CALCUL DE W1=DF/D(SIGMA)--------------------------------------
  530.  
  531. C---------ELEMENTS MASSIFS----------------------------------------------
  532.  
  533. IF(MFR1.EQ.1.OR.MFR1.EQ.31) THEN
  534.  
  535. F(1)=(2.D0*SX(1)-SX(2)-SX(3))/3.D0
  536. F(2)=(2.D0*SX(2)-SX(1)-SX(3))/3.D0
  537. F(3)=(2.D0*SX(3)-SX(1)-SX(2))/3.D0
  538. DO I=4,NSTRSS
  539. F(I)=SX(I)
  540. ENDDO
  541. DO I=1,3
  542. W1(I)=1.5D0*F(I)/SEQ
  543. Z1(I)=W1(I)
  544. ENDDO
  545. DO I=4,NCONT
  546. W1(I)=3.D0*F(I)/SEQ
  547. Z1(I)=1.5D0*F(I)/SEQ
  548. ENDDO
  549.  
  550. C---------COQUES MINCES-------------------------------------------------
  551.  
  552. ELSE IF(MFR1.EQ.3.OR.MFR1.EQ.9) THEN
  553.  
  554. IF(IFOURB.GE.1) THEN
  555. W1(1)=(2.D0*SX(1)-SX(2))/(2.D0*SEQ*EP1)
  556. W1(2)=(2.D0*SX(2)-SX(1))/(2.D0*SEQ*EP1)
  557. W1(3)=3.D0*SX(3)/(SEQ*EP1)
  558. W1(4)=3.D0*WORK(2)*(2.D0*SX(4)-SX(5))/(SEQ*EP1*EP1)
  559. W1(5)=3.D0*WORK(2)*(2.D0*SX(5)-SX(4))/(SEQ*EP1*EP1)
  560. W1(6)=18.D0*WORK(2)*SX(6)/(SEQ*EP1*EP1)
  561. ELSE
  562. W1(1)=(2.D0*SX(1)-SX(2))/(2.D0*SEQ*EP1)
  563. W1(2)=(2.D0*SX(2)-SX(1))/(2.D0*SEQ*EP1)
  564. W1(3)=3.D0*WORK(2)*(2.D0*SX(3)-SX(4))/(SEQ*EP1*EP1)
  565. W1(4)=3.D0*WORK(2)*(2.D0*SX(4)-SX(3))/(SEQ*EP1*EP1)
  566. ENDIF
  567.  
  568.  
  569. C---------COQUES EPAISSES-----------------------------------------------
  570. ELSE IF(MFR1.EQ.5) THEN
  571. F(1)=(2.D0*SX(1)-SX(2))/3.D0
  572. F(2)=(2.D0*SX(2)-SX(1))/3.D0
  573. DO I=3,5
  574. F(I)=SX(I)
  575. ENDDO
  576. DO I=1,2
  577. W1(I)=1.5D0*F(I)/SEQ
  578. Z1(I)=W1(I)
  579. ENDDO
  580. DO I=3,5
  581. W1(I)=3.D0*F(I)/SEQ
  582. Z1(I)=1.5D0*F(I)/SEQ
  583. ENDDO
  584.  
  585. C---------POUTRES-------------------------------------------------------
  586.  
  587. ELSE IF(MFR1.EQ.7) THEN
  588.  
  589. DO J=1,NCONT
  590. DO I=1,NCONT
  591. DDA(I,J)=0.D0
  592. ENDDO
  593. ENDDO
  594. DDA(1,1)=YUNG
  595. DDA(4,4)=0.5D0*YUNG/(1.D0+XNU)
  596. DDA(5,5)=YUNG
  597. DDA(6,6)=YUNG
  598. W1(1)=SX(1)/SEQ
  599. W1(2)=0.D0
  600. W1(3)=0.D0
  601. W1(4)=SX(4)/SEQ
  602. W1(5)=SX(5)/SEQ
  603. W1(6)=SX(6)/SEQ
  604. ENDIF
  605.  
  606. DO I=1,NCONT
  607. W2(I)=0.D0
  608. ENDDO
  609.  
  610. IF(MFR1.EQ.7) THEN
  611. DO J=1,NCONT
  612. XFLO1=W1(J)
  613. DO I=1,NCONT
  614. W2(I)=W2(I)+XFLO1*DDA(I,J)
  615. ENDDO
  616. ENDDO
  617.  
  618. ELSE
  619. IF((mfr.eq.1.or.mfr.eq.31).and. ifourb.eq.-2) then
  620. DO J=1,NCONT
  621. XFLO1=W1(J)
  622. DO I=1,NCONT
  623. W2(I)=W2(I)+XFLO1*DDA(I,J)
  624. ENDDO
  625. ENDDO
  626. else
  627. DO J=1,NCONT
  628. XFLO1=W1(J)
  629. DO I=1,NCONT
  630. W2(I)=W2(I)+XFLO1*DDAUX(I,J)
  631. ENDDO
  632. ENDDO
  633. endif
  634. ENDIF
  635. COEF=0.D0
  636. DO I=1,NCONT
  637. COEF=COEF+W1(I)*W2(I)
  638. ENDDO
  639.  
  640. C-----------------------------------------------------------------------
  641. C PLASTICITE PARFAITE, ECROUISSAGE ISOTROPE ET CINEMATIQUE ZIEGLER
  642. C-----------------------------------------------------------------------
  643.  
  644. IF (INPLAS.EQ.1.OR.INPLAS.EQ.4.OR.INPLAS.EQ.5
  645. $ .OR.INPLAS.EQ.87) THEN
  646.  
  647. IF (INPLAS.EQ.4) THEN
  648. RP=0.D0
  649. C=HECL
  650. ELSE
  651. CALL TRACTI(PENTE,EPST,SIG,EPS,LTRAC,1,IBI)
  652.  
  653. IF (IBI.EQ.1) THEN
  654. KERRE=75
  655. GOTO 1000
  656. ENDIF
  657.  
  658. IF (INPLAS.EQ.1) THEN
  659. RP=0.D0
  660. C=0.D0
  661. ELSE IF(INPLAS.EQ.5.OR.INPLAS.EQ.87) THEN
  662. RP=PENTE
  663. C=0.D0
  664. ENDIF
  665. ENDIF
  666. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  667. c*dbg write(ioimp,*) 'RP,C=',RP,C,EPST
  668. c*dbg endif
  669.  
  670. DENOM=COEF+C+RP
  671. DELTA=PHI/DENOM
  672. DMU=C*DELTA/SEQ
  673.  
  674. DO I=1,NCONT
  675. DSIG(I)=-DELTA*W2(I)
  676. DSPHER1(I)=DMU*SX(I)
  677. ENDDO
  678.  
  679. * Cas des contraintes planes en massif
  680. if ((mfr.eq.1.or.mfr.eq.31).and.ifourb.eq.-2) then
  681.  
  682. bb= abs(dsig(3)+ s(3) )
  683. r_z = dsig(3) * cte_cp
  684. sx(3)= sx3in - dsig(3)
  685. sx(1)= sx1in - r_z
  686. sx(2)= sx2in - r_z
  687. SEQ=VONMIS0(SX,NSTRSS,MFR1,IFOURB,EP1,ALFAH)
  688. s(3)= - dsig(3)
  689. s(1)= s1in - r_z
  690. s(2)= s2in - r_z
  691. if (bb.gt.cri0) then
  692. if (iderin.eq.0) then
  693. niter=niter - 1
  694. endif
  695. iderin=iderin+1
  696. if (iderin.gt.50) then
  697. write(ioimp,*) ' probleme dans iterations internes'
  698. KERRE=2
  699. GO TO 1000
  700. endif
  701. go to 10
  702. endif
  703. DMU=C*DELTA/SEQ
  704. DO I=1,NCONT
  705. DSPHER1(I)=DMU*SX(I)
  706. ENDDO
  707. endif
  708.  
  709. iderin=0
  710. DP=DELTA
  711. DR=RP*DP
  712.  
  713. ELSE
  714.  
  715. C-----------------------------------------------------------------------
  716. C MODELE DE CHABOCHE
  717. C-----------------------------------------------------------------------
  718.  
  719. C---------UNIQUEMENT POUR LES ELEMENTS MASSIFS--------------------------
  720.  
  721. XPRO1=0.D0
  722. XPRO2=0.D0
  723. DO I=1,NCONT
  724. XPRO1=XPRO1+W1(I)*SPHER1(I)
  725. XPRO2=XPRO2+W1(I)*SPHER2(I)
  726. ENDDO
  727.  
  728. FIP=1.D0+(PSI-1.D0)*EXP(-OME*EPST)
  729.  
  730. DENOM=COEF+(A1*C1+A2*C2)*FIP-C1*XPRO1-C2*XPRO2+B*RMmRR
  731. DELTA=PHI/DENOM
  732.  
  733. DO I=1,NCONT
  734. DSIG(I)=-DELTA*W2(I)
  735. DSPHER1(I)=(2.D0*A1*FIP*Z1(I)/3.D0-SPHER1(I))*C1*DELTA
  736. DSPHER2(I)=(2.D0*A2*FIP*Z1(I)/3.D0-SPHER2(I))*C2*DELTA
  737. ENDDO
  738.  
  739. DR=B* RMmRR *DELTA
  740. DP=DELTA
  741. ENDIF
  742.  
  743. RR=RR+DR
  744. EPST=EPST+DP
  745. EPST=MAX(0.D0,EPST)
  746. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  747. c*dbg write(ioimp,*) 'MAJ=',EPST,DP,RR,DR
  748. c*dbg endif
  749.  
  750. IF(MFR1.EQ.3.OR.MFR1.EQ.9) THEN
  751. IF(IFOURB.GE.1) THEN
  752. DO I=1,3
  753. DSIG(I)=DSIG(I)*EP1
  754. ENDDO
  755. DO I=4,6
  756. DSIG(I)=EP1*EP1*DSIG(I)/(6.D0*WORK(2))
  757. ENDDO
  758. ELSE
  759. DSIG(1)=DSIG(1)*EP1
  760. DSIG(2)=DSIG(2)*EP1
  761. DSIG(3)=EP1*EP1*DSIG(3)/(6.D0*WORK(2))
  762. DSIG(4)=EP1*EP1*DSIG(4)/(6.D0*WORK(2))
  763. ENDIF
  764. ENDIF
  765. C mise a jour des contraintes
  766. DO I=1,NCONT
  767. S(I)=S(I)+DSIG(I)
  768. SPHER1(I)=SPHER1(I)+DSPHER1(I)
  769. SPHER2(I)=SPHER2(I)+DSPHER2(I)
  770. SPHER(I)=SPHER1(I)+SPHER2(I)
  771. SX(I)=S(I)-SPHER(I)
  772. ENDDO
  773. if(ifourb.eq.-2.and.(mfr.eq.1.or.mfr.eq.31)) then
  774. s(3)=0.d0
  775. endif
  776.  
  777. yxsxii = 0.D0
  778. if (icow21.gt.0) then
  779. * kich borne trace sigma
  780. ytr = sx(1) + sx(2) + sx(3)
  781. ytr = ytr/3.D0
  782. if (ytr.gt.xlimp) then
  783. yxsxii = ytr - xlimp
  784. else if((ytr + xlimp).lt.0) then
  785. yxsxii = ytr + xlimp
  786. endif
  787. dsigii = dsigt(1) + dsigt(2) + dsigt(3)
  788. if (dsigii.ne.0.D0) then
  789. do jj = 1,3
  790. s(jj) = s(jj) - (dsigt(jj)/dsigii*3.D0*yxsxii)
  791. sx(jj) = sx(jj) - (dsigt(jj)/dsigii*3.D0*yxsxii)
  792. enddo
  793. else
  794. do jj = 1,3
  795. s(jj) = s(jj) - yxsxii
  796. sx(jj) = sx(jj) - yxsxii
  797. enddo
  798. endif
  799. endif
  800.  
  801. SEQ=VONMIS0(SX,NSTRSS,MFR1,IFOURB,EP1,ALFAH)
  802.  
  803. C---------CONTRAINTES PLANES--------------------------------------------
  804.  
  805. IF(IFOURB.EQ.-2) THEN
  806.  
  807. IF(MFR1.EQ.1.OR.MFR1.EQ.31) THEN
  808. F(1)=(2.D0*SX(1)-SX(2)-SX(3))/3.D0
  809. F(2)=(2.D0*SX(2)-SX(1)-SX(3))/3.D0
  810. F(3)=(2.D0*SX(3)-SX(1)-SX(2))/3.D0
  811. DO I=4,NSTRSS
  812. F(I)=SX(I)
  813. ENDDO
  814. DO I=1,3
  815. W1(I)=1.5D0*F(I)/SEQ
  816. ENDDO
  817. DO I=4,NSTRSS
  818. W1(I)=3.D0*F(I)/SEQ
  819. ENDDO
  820.  
  821. ELSE IF(MFR1.EQ.3.OR.MFR1.EQ.9) THEN
  822. AUX=EP1*EP1*EP1*EP1
  823. W1(1)=(2.D0*SX(1)-SX(2))/(2.D0*SEQ*EP1*EP1)
  824. W1(2)=(2.D0*SX(2)-SX(1))/(2.D0*SEQ*EP1*EP1)
  825. W1(3)=18.D0*ALFAH*(2.D0*SX(3)-SX(4))/(SEQ*AUX)
  826. W1(4)=18.D0*ALFAH*(2.D0*SX(4)-SX(3))/(SEQ*AUX)
  827. ENDIF
  828.  
  829. DO I=1,NSTRSS
  830. DEPSP(I)=DEPSP(I)+DELTA*W1(I)
  831. ENDDO
  832. ENDIF
  833.  
  834. C-----------------------------------------------------------------------
  835. C TEST
  836. C CALCUL DE LA NOUVELLE VALEUR DE PHI
  837. C-----------------------------------------------------------------------
  838. IF (INPLAS.EQ.5.OR.INPLAS.EQ.87) THEN
  839. CALL TRACTI(SI,EPST,SIG,EPS,LTRAC,2,IBI)
  840. C* Modele de CHABOCHE (prise en compte ecrouissage)
  841. ELSE IF (INPLAS.EQ.7 .OR. INPLAS.EQ.12) THEN
  842. RMmRR = (RM - R0) * EXP(-B*EPST)
  843. SI = RM - RMmRR
  844. ELSE IF (INPLAS.EQ.4) THEN
  845. SI=SLI0
  846. ELSE
  847. SI=RR+SI0
  848. ENDIF
  849. PHI=SEQ-SI
  850.  
  851. PETI=1.D-8
  852. APHI=ABS(PHI)
  853. APHI0=ABS(PHI0)
  854. TEST=max(PETI*APHI0,XZPREC*100.D0*SEQ)
  855. *sg TEST=PETI*APHI0
  856. *sg write(ioimp,*) 'niter,phi,phi0,si,seq,rmmrr,test=',
  857. *sg $ niter,phi,phi0,si,seq,rmmrr,test
  858. IF(NITER.GT.50) THEN
  859. if(itracb.eq.0) then
  860. itracb=1
  861. go to 2222
  862. endif
  863. C write(ioimp,*) ' CCOIN0, APHI, TEST',APHI,TEST
  864. C write(ioimp,*) ' CCOIN0, KERRE = 2, NITER =',NITER
  865. KERRE=-2
  866. call soucis(268)
  867. C GO TO 1000
  868. ENDIF
  869. IF(APHI.LE.TEST.OR.KERRE.LT.0) THEN
  870. IF (KERRE.LT.0) KERRE = 0
  871.  
  872. IF(MFR1.EQ.7) THEN
  873. DO I=1,NCONT
  874. S(I)=S(I)/DIV(I)
  875. SPHER1(I)=SPHER1(I)/DIV(I)
  876. SPHER2(I)=SPHER2(I)/DIV(I)
  877. ENDDO
  878. ENDIF
  879.  
  880. C---------TOUTES FORMULATIONS SAUF CONTRAINTES PLANES-------------------
  881.  
  882. IF(IFOURB.NE.-2) THEN
  883. DO I=1,NCONT
  884. DS(I)=S(I)-SIGB(I)
  885. ENDDO
  886. CALL EPSIG0(DS,DDEPSE,MFR1,IFOURB,YUNG,XNU,WORK,NSTRSS)
  887. DO I=1,NCONT
  888. DEPSE(I)=DEPST(I)+DDEPSE(I)
  889. DEPSP(I)=DEPST(I)-DEPSE(I)
  890. ENDDO
  891. ENDIF
  892.  
  893. DO I=1,NSTRSS
  894. SIGF(I)=S(I)
  895. DEFP(I)=DEPSP(I)
  896. ENDDO
  897. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  898. c*dbg write(ioimp,*) 'SIGF',(SIGF(iof),iof=1,NSTRSS)
  899. c*dbg endif
  900.  
  901. C---------COQUES AVEC CISAILLEMENT TRANSVERSE---------------------------
  902.  
  903. IF(MFR1.EQ.9) THEN
  904. DEFP(7)=0.D0
  905. DEFP(8)=0.D0
  906. SIGF(7)=SIGB(7)
  907. SIGF(8)=SIGB(8)
  908. ENDIF
  909.  
  910. VARF(1)=EPST
  911. IF(INPLAS.EQ.4.OR.INPLAS.EQ.7.OR.INPLAS.EQ.12) THEN
  912. DO I=1,NSTRSS
  913. VARF(I+1)=SPHER1(I)
  914. ENDDO
  915. IF(INPLAS.EQ.12) THEN
  916. DO I=1,NSTRSS
  917. VARF(NSTRSS+1+I)=SPHER2(I)
  918. ENDDO
  919. ENDIF
  920. ELSE IF(INPLAS.EQ.87) THEN
  921. DO I=1,NSTRSS
  922. VARF(1+I)=DEPST(I)
  923. ENDDO
  924. ENDIF
  925. KERRE=0
  926. c*dbg if (ib.eq.1 .and. igau.eq.1) then
  927. c*dbg write(ioimp,*) 'VARF',EPST,'---',(VARF(iof),iof=2,1+NSTRSS)
  928. c*dbg endif
  929. RETURN
  930.  
  931. ELSE
  932. sx1in=sx(1)
  933. sx2in=sx(2)
  934. s1in=s(1)
  935. s2in=s(2)
  936. GOTO 10
  937. ENDIF
  938.  
  939. 1000 CONTINUE
  940.  
  941. C RETURN
  942. END
  943.  
  944.  
  945.  
  946.  
  947.  
  948.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales