Télécharger betocy.eso

Retour à la liste

Numérotation des lignes :

  1. C BETOCY SOURCE PV 11/03/07 21:15:06 6885
  2. SUBROUTINE BETOCY(WRK0,WRK1,WRK2,NCURVT,NCURVC,KERRE)
  3. IMPLICIT INTEGER(I-N)
  4. IMPLICIT REAL*8(A-H,O-Z)
  5. C
  6. C MODELE DE MACONNERIE ENDOMMAGEABLE EN CYCLIQUE
  7. C
  8. C=======================================================================
  9. C CETTE ROUTINE EST APPELE DANS ECOUL2
  10. C
  11. C ENTREES:
  12. C -------
  13. C SIG0 = CONTRAINTES AU DEBUT DU PAS D'INTEGRATION
  14. C DSIGT = INCREMENT DE CONT. CALCULE ELASTIQUEMENT A PARTIR DE
  15. C L'INCREMENT DES DEFORMATIONS TOTALES
  16. C VAR0 = VARIABLES INTERNES AU DEBUT DU PAS D'INTEGRATION
  17. C
  18. C XMAT = CARACTERISTIQUES MECANIQUES DU MATERIAU
  19. C
  20. C SORTIES:
  21. C -------
  22. C SIGF = CONTR. A LA FIN DU PAS D'INTEGRATION
  23. C VARF = VARIABLES INTERNES A LA FIN DU PAS D'INTEGRATION
  24. C DEFP = INCREMENT DES DEFORM. PLASTIQUES A LA FIN DU PAS
  25. C D'INTEGRATION
  26. C KERRE = INDICE QUI REGIT LES ERREURS
  27. C
  28. C=========================================================================
  29. C VARIABLES NECESSAIRES
  30. C=========================================================================
  31. -INC CCOPTIO
  32. C
  33. C Segment pour le Broyden
  34. C
  35. SEGMENT QUASIN
  36. REAL*8 XH0(MN,MN)
  37. REAL*8 WUP(MN,ITM),YUP(MN,ITM)
  38. INTEGER IT
  39. REAL*8 R(MN)
  40. REAL*8 D(MN),DES(MN)
  41. ENDSEGMENT
  42. C Segment des materiaux
  43. SEGMENT WRK0
  44. REAL*8 XMAT(NCXMAT)
  45. ENDSEGMENT
  46. C Segment des courbes d'ecrouissage
  47. SEGMENT WRK2
  48. REAL*8 TRAC(LTRAC)
  49. ENDSEGMENT
  50. C Segment des entrees/sortie de ecoul2
  51. C Les contraintes sont sous la forme [SXX,SYY,0,SXY]
  52. C Vecteurs de contrainte totale
  53. C 4 composantes pour les contraintes SIGT
  54. SEGMENT WRK1
  55. REAL*8 DDHOOK(LHOOK,LHOOK),SIG0(NSTRS),DEPST(NSTRS)
  56. REAL*8 SIGF(NSTRS),VAR0(NVARI),VARF(NVARI)
  57. REAL*8 DEFP(NSTRS),XCAR(ICARA)
  58. ENDSEGMENT
  59. DIMENSION DSIGT(4)
  60. DIMENSION SIGT(4), SIGI(4), DSIGI(4), DSIGIDI(4)
  61. C Matrice de compliance (C-1)
  62. DIMENSION CINV(3)
  63. C Variables d'ecrouissage isotrope
  64. DIMENSION XKISO(2)
  65. C Travaux plastiques cycliques et Deformations plastiques equivalentes
  66. C associes au mecanisme de compression de la grande surf.
  67. REAL*8 WORK
  68. DIMENSION EPSI(2),EPSI0(2),EPSII(2)
  69. C Contraintes principales et variables d'ecrouissage cinematique
  70. DIMENSION SIGTP(2), SIGHP(2),SIGHP0(2)
  71. C Cos et Sin de l'angle definissant la premiere direction principale
  72. REAL*8 CPHI, CPHI0, CPHI00, SPHI, SPHI0, SPHI00
  73. C Indicateur du contact
  74. DIMENSION KONTAC(4),KONTA0(4)
  75. C Vecteurs contrainte transitoire
  76. DIMENSION SITPRO(4)
  77. C Jacobien initial
  78. DIMENSION XJ00(3,3)
  79. C Variable phi apparaissant au denominateur (calcul de la normale et de A)
  80. REAL*8 PHI
  81. C Matrices contenant A-1
  82. C AINV --> Inverse de A pour la grande surface
  83. C [ A B 0] A=AINV(1)
  84. C AINVGRANDESURFACE = [-B A 0] B=AINV(2)
  85. C [ 0 0 C] C=AINV(3)
  86. DIMENSION AINV(3)
  87. C Vecteurs normaux aux surface
  88. C XNOR1 --> 1eme mecanisme
  89. C XNOR2 --> 2eme mecanisme
  90. DIMENSION XNOR1(4), XNOR2(4)
  91. C
  92. C Multiplicateurs plastiques 2 mecanismes + 2 surfaces dans GAMBDA
  93. C Valeurs des criteres dans FFF
  94. C Valeurs des limites dans RESIST
  95. C Ordre pour les 3 tableaux:
  96. C 1 --> Traction petite surface direction 1
  97. C 2 --> Traction petite surface direction 2
  98. C 3 --> Compression petite surface direction 1
  99. C 4 --> Compression petite surface direction 2
  100. C 5 --> Traction grande surface
  101. C 6 --> Compression grande surface
  102. DIMENSION GAMBDA(6), FFF(6), RESIST(6)
  103. C Valeur necessaire pour le calcul des residus de la grande surface
  104. DIMENSION XKISO0(2),XKISOI(2)
  105. C
  106. C Separation de l'increment en deux ou trois parties dans le cas du contact
  107. C
  108. REAL*8 SPLITT
  109. REAL*8 XFFFS
  110. C
  111. C Pour differencier les differents cas
  112. C ICAS --> NUMERO DU MULTIPLICATEUR PLASTIQUE A METTRE A JOUR
  113. C IHARD --> NUMERO DE LA VARIABLE D'ECROUISSAGE CINEMATIQUE
  114. C ISGNTC --> Signe indiquant si on est en traction ou compression
  115. C ISGN2 --> Signe indiquant si on prend la plus grande ou la plus
  116. C petite contrainte principale
  117. C NCORNER --> Revient-on sur la pointe ?
  118. INTEGER ICAS, IHARD, ISGNTC, ISGN2
  119. INTEGER ICASB, IHARDB, ISGNTCB, ISGN2B
  120. INTEGER NCORNER
  121. C
  122. C Variable intermediaire valant la resistance de la grande surface
  123. C
  124. DIMENSION FG(2)
  125. C
  126. C Critere de convergence
  127. C
  128. DIMENSION EPSIL(6)
  129. CHARACTER*8 CMATE
  130. C========================================================================
  131. C Parametres utiles
  132. C========================================================================
  133. PARAMETER (XZER=0.D0,UN=1.D0,DEUX=2.D0,UNDEMI=.5D0)
  134. PARAMETER (UNRADE=0.707106781D0)
  135. PARAMETER (EPSILO=1.D-8,EPSILO2=-1.D-10)
  136. C========================================================================
  137. C
  138. C CALCUL DE DSIGT
  139. C
  140. C calcul de la matrice elastique
  141. C
  142. CMATE = 'ISOTROPE'
  143. KCAS=1
  144. CALL DOHMAS(XMAT,CMATE,IFOUR,4,KCAS,DDHOOK,IRTD)
  145. DO IE1=1,4
  146. DSIGT(IE1)=XZER
  147. DO IE2=1,4
  148. DSIGT(IE1)=DSIGT(IE1)+DDHOOK(IE1,IE2)*DEPST(IE2)
  149. ENDDO
  150. ENDDO
  151. DSIGT(3)=XZER
  152. C=========================================================================
  153. C LECTURE DES PARAMETRES DU MATERIAU
  154. C
  155. C=========================================================================
  156. C Module d'Young YOUN, Poisson XNU, Shear GG
  157. C Ecrouissage cinematique XHH1 ,
  158. YOUN = XMAT(1)
  159. XNU = XMAT(2)
  160. GGG = YOUN/(2*(UN+XNU))
  161. XHH1 = XMAT(5)
  162. CPI = YOUN/(UN-XNU)
  163. C
  164. RESIST(1) = XMAT(6)
  165. RESIST(2) = XMAT(6)
  166. RESIST(3) = XMAT(7)
  167. RESIST(4) = XMAT(7)
  168. RESIST(5) = XMAT(8)
  169. RESIST(6) = XMAT(9)
  170. C Deformation plastique equivalente et travail plastique de ref.
  171. C pour la dissipation cyclique
  172. EPS0 = XMAT(10)
  173. WOR0 = XMAT(11)
  174. DEGRAD = EPS0/WOR0
  175. XLCAT = XMAT(14)
  176. XLCAC = XMAT(15)
  177. C
  178. C
  179. C Evolution des variables d'ecrouissage isotrope
  180. C
  181. C Voila les positions des points d'entree des deux courbes dans TRAC
  182. IPT = 1
  183. IPC = IPT+2*NCURVT
  184. C NCURVT : nombre de point de la courbe de traction
  185. C NCURVC : nombre de point de la courbe de compression
  186. C===========================================================
  187. C Fin lecture des caracteristiques
  188. C On passe aux contraintes et variables d'etat
  189. C===========================================================
  190. C
  191. C CONTRAINTES ELASTIQUES
  192. C
  193. SIGI(1) = SIG0(1)
  194. SIGI(2) = SIG0(2)
  195. SIGI(4) = SIG0(4)
  196. C DELTA SIGMA TRIAL
  197. DSIGI(1) = DSIGT(1)
  198. DSIGI(2) = DSIGT(2)
  199. DSIGI(4) = DSIGT(4)
  200. C
  201. C ECROUISSAGE CINEMATIQUE
  202. C
  203. SIGHP0(1) = VAR0(1)
  204. SIGHP0(2) = VAR0(2)
  205. C
  206. C ECROUISSAGE ISOTROPE
  207. C
  208. EPSI0(1) = VAR0(3)
  209. EPSI0(2) = VAR0(4)
  210. C
  211. CPHI00 = VAR0(5) + 1.D0
  212. CPHI0 = CPHI00
  213. SPHI00 = VAR0(6)
  214. SPHI0 = SPHI00
  215. C
  216. KONTA0(1) = nint(VAR0( 7))
  217. KONTA0(2) = nint(VAR0( 8))
  218. KONTA0(3) = nint(VAR0( 9))
  219. KONTA0(4) = nint(VAR0(10))
  220. C
  221. C Variables de contact
  222. C
  223. SPLITT = 0.D0
  224. NCON = 0
  225. IDICHO = 0
  226. NDICHO = 0
  227. NPASS = 0
  228. DSIGIDI(1) = XZER
  229. DSIGIDI(2) = XZER
  230. DSIGIDI(4) = XZER
  231. C===============================================
  232. C Fin de la lecture des contraintes et des variables d'etat
  233. C================================================
  234. C
  235. C Premier calcul des residus puis DRIVER
  236. C================================================
  237. C
  238. C Entree quand on revient des routines, dernier test avant de finir
  239. C
  240. 10 CONTINUE
  241. SIGT(1) = SIGI(1) + DSIGI(1)
  242. SIGT(2) = SIGI(2) + DSIGI(2)
  243. SIGT(4) = SIGI(4) + DSIGI(4)
  244. CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  245. EPSI(1) = EPSI0(1)
  246. EPSI(2) = EPSI0(2)
  247. CALL VALEUR(EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  248. &XKISO(1),DUMM,KERRE)
  249. CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  250. &XKISO(2),DUMM,KERRE)
  251. SIGHP(1) = SIGHP0(1)
  252. SIGHP(2) = SIGHP0(2)
  253. NPASS2 = 0
  254. KONTAC(1) = KONTA0(1)
  255. KONTAC(2) = KONTA0(2)
  256. KONTAC(3) = KONTA0(3)
  257. KONTAC(4) = KONTA0(4)
  258. C
  259. C Calcul des residus initiaux
  260. C
  261. NMEC = 0
  262. FFF(1) = XZER
  263. FFF(2) = XZER
  264. FFF(3) = XZER
  265. FFF(4) = XZER
  266. FFF(5) = XZER
  267. FFF(6) = XZER
  268. FINT1 = XZER
  269. FINT2 = XZER
  270. GAMBDA(1) = XZER
  271. GAMBDA(2) = XZER
  272. GAMBDA(3) = XZER
  273. GAMBDA(4) = XZER
  274. GAMBDA(5) = XZER
  275. GAMBDA(6) = XZER
  276. ICAS = 0
  277. ICASB = 0
  278. C
  279. C Critere de convergence
  280. C
  281. EPSIL(1) = ABS(RESIST(1)*EPSILO)
  282. EPSIL(2) = ABS(RESIST(2)*EPSILO)
  283. EPSIL(3) = ABS(RESIST(3)*EPSILO)
  284. EPSIL(4) = ABS(RESIST(4)*EPSILO)
  285. EPSIL(5) = ABS(RESIST(5)*EPSILO)
  286. EPSIL(6) = ABS(RESIST(6)*EPSILO)
  287. C
  288. C
  289. IF (KONTAC(1).EQ.0) THEN
  290. CALL CRITPE (SIGTP(1), SIGHP0(1), RESIST(1), 1, FFF(1))
  291. FINT1 = 0.
  292. IF (FFF(1).GE.EPSIL(1)) THEN
  293. NMEC = NMEC + 1
  294. IF (NMEC.EQ.1) THEN
  295. ICAS = 1
  296. ISGNTC = 1
  297. IHARD = 1
  298. ELSE
  299. ICASB = 1
  300. ISGNTCB = 1
  301. IHARDB = 1
  302. ENDIF
  303. ENDIF
  304. ELSE
  305. FG(1) = XKISO(1)*RESIST(5)
  306. CALL CRITPE (SIGTP(1), XZER, FG(1), 1, FINT1)
  307. ENDIF
  308. C
  309. IF (KONTAC(2).EQ.0) THEN
  310. CALL CRITPE (SIGTP(2), SIGHP0(2), RESIST(2), 1, FFF(2))
  311. FINT2 = 0.
  312. IF (FFF(2).GE.EPSIL(2)) THEN
  313. NMEC = NMEC + 1
  314. IF (NMEC.EQ.1) THEN
  315. ICAS = 2
  316. ISGNTC = 1
  317. IHARD = 2
  318. ELSE
  319. ICASB = 2
  320. ISGNTCB = 1
  321. IHARDB = 2
  322. ENDIF
  323. ENDIF
  324. ELSE
  325. FG(1) = XKISO(1)*RESIST(5)
  326. CALL CRITPE (SIGTP(2), XZER, FG(1), 1, FINT2)
  327. ENDIF
  328. C
  329. IF ((FINT1.GE.EPSIL(5)).OR.(FINT2.GE.EPSIL(5))) THEN
  330. NMEC = NMEC + 1
  331. IF (FINT1.GE.FINT2) THEN
  332. FFF(5) = FINT1
  333. IF (NMEC.EQ.1) THEN
  334. IHARD = 1
  335. ICAS = 5
  336. ISGNTC = 1
  337. ELSE
  338. IHARDB = 1
  339. ICASB = 5
  340. ISGNTCB = 1
  341. ENDIF
  342. ELSE
  343. FFF(5) = FINT2
  344. IF (NMEC.EQ.1) THEN
  345. IHARD = 2
  346. ICAS = 5
  347. ISGNTC = 1
  348. ELSE
  349. IHARDB = 2
  350. ICASB = 5
  351. ISGNTCB = 1
  352. ENDIF
  353. ENDIF
  354. ENDIF
  355. C
  356. IF (KONTAC(3).EQ.0) THEN
  357. CALL CRITPE (SIGTP(1), SIGHP0(1), RESIST(3), -1, FFF(3))
  358. FINT1 = 0.
  359. IF (FFF(3).GE.EPSIL(3)) THEN
  360. NMEC = NMEC + 1
  361. IF (NMEC.EQ.1) THEN
  362. ICAS = 3
  363. ISGNTC = -1
  364. IHARD = 1
  365. ELSE
  366. ICASB = 3
  367. ISGNTCB = -1
  368. IHARDB = 1
  369. ENDIF
  370. ENDIF
  371. ELSE
  372. FG(1) = XKISO(2)*RESIST(6)
  373. CALL CRITPE (SIGTP(1), XZER, FG(1), -1, FINT1)
  374. ENDIF
  375. C
  376. IF (KONTAC(4).EQ.0) THEN
  377. CALL CRITPE (SIGTP(2), SIGHP0(2), RESIST(4), -1, FFF(4))
  378. FINT2 = 0.
  379. IF (FFF(4).GE.EPSIL(4)) THEN
  380. NMEC = NMEC + 1
  381. IF (NMEC.EQ.1) THEN
  382. ICAS = 4
  383. ISGNTC = -1
  384. IHARD = 2
  385. ELSE
  386. ICASB = 4
  387. ISGNTCB = -1
  388. IHARDB = 2
  389. ENDIF
  390. ENDIF
  391. ELSE
  392. FG(1) = XKISO(2)*RESIST(6)
  393. CALL CRITPE (SIGTP(2), XZER, FG(1), -1, FINT2)
  394. ENDIF
  395. C
  396. IF ((FINT1.GE.EPSIL(6)).OR.(FINT2.GE.EPSIL(6))) THEN
  397. NMEC = NMEC + 1
  398. IF (FINT1.GE.FINT2) THEN
  399. FFF(6) = FINT1
  400. IF (NMEC.EQ.1) THEN
  401. IHARD = 1
  402. ICAS = 6
  403. ISGNTC = -1
  404. ELSE
  405. IHARDB = 1
  406. ICASB = 6
  407. ISGNTCB = -1
  408. ENDIF
  409. ELSE
  410. FFF(6) = FINT2
  411. IF (NMEC.EQ.1) THEN
  412. IHARD = 2
  413. ICAS = 6
  414. ISGNTC = -1
  415. ELSE
  416. IHARDB = 2
  417. ICASB = 6
  418. ISGNTCB = -1
  419. ENDIF
  420. ENDIF
  421. ENDIF
  422. C
  423. C
  424. C Le tableau FFF contient les residus
  425. C
  426. C================================================
  427. C DRIVER
  428. C================================================
  429. IF (NMEC.EQ.0) GOTO 9999
  430. C
  431. IF (NMEC.EQ.1) GOTO 1000
  432. C
  433. IF (NMEC.EQ.2) GOTO 4000
  434. C
  435. GOTO 9999
  436. C==================================================
  437. C CAS 1:
  438. C Un seul mecanisme est active
  439. C==================================================
  440. 1000 CONTINUE
  441. C
  442. CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  443. C
  444. C Mis a jour de isgn2
  445. C
  446. IF (IHARD.EQ.1) THEN
  447. IF (SIGTP(1).GE.SIGTP(2)) THEN
  448. ISGN2 = 1
  449. ELSE
  450. ISGN2 = -1
  451. ENDIF
  452. HARDC = 2
  453. ELSE
  454. IF (IHARD.EQ.2) THEN
  455. IF (SIGTP(2).GE.SIGTP(1)) THEN
  456. ISGN2 = 1
  457. ELSE
  458. ISGN2 = -1
  459. ENDIF
  460. ENDIF
  461. HARDC = 1
  462. ENDIF
  463. C
  464. IF (IHARD.EQ.1) THEN
  465. IHARDC = 2
  466. ELSE
  467. IHARDC = 1
  468. ENDIF
  469. C
  470. NCORNER = 0
  471. C
  472. 1010 CONTINUE
  473. C
  474. C Definition du segment
  475. C
  476. ITM = 30
  477. MN = 1
  478. SEGINI,QUASIN
  479. IT = -1
  480. C
  481. C Calcul de xhh2 et xkiso
  482. C
  483. IF (ICAS.GE.5) THEN
  484. IF (ICAS.EQ.5) THEN
  485. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  486. &XKISO(1),XHH2,KERRE)
  487. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  488. &XKISO(2),DUMM,KERRE)
  489. ELSE
  490. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  491. &XKISO(1),DUMM,KERRE)
  492. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  493. &XKISO(2),XHH2,KERRE)
  494. ENDIF
  495. EPSI(1) = EPSI0(1)
  496. XKISO0(1) = XKISO(1)
  497. EPSI(2) = EPSI0(2)
  498. XKISO0(2) = XKISO(2)
  499. XHH2I = XHH2
  500. EPSII(1) = EPSI0(1)
  501. XKISOI(1) = XKISO(1)
  502. EPSII(2) = EPSI0(2)
  503. XKISOI(2) = XKISO(2)
  504. ENDIF
  505. C
  506. DO 1011,IE1=1,6
  507. FFF(IE1) = XZER
  508. 1011 CONTINUE
  509. C
  510. IF (ICAS.LE.4) THEN
  511. CALL CRITPE (SIGTP(IHARD),SIGHP0(IHARD),RESIST(ICAS),
  512. & ISGNTC,FFF(ICAS))
  513. ELSE
  514. FG(1) = XKISOI(ICAS-4)*RESIST(ICAS)
  515. CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,FINT1)
  516. CALL CRITPE (SIGTP(IHARDC),XZER,FG(1),ISGNTC,FINT2)
  517. FFF(ICAS) = MAX(FINT1,FINT2)
  518. ENDIF
  519. C
  520. C Calcul de la normale et du jacobien
  521. C
  522. IF (ICAS.LE.4) THEN
  523. CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,1,XNOR1)
  524. XJ00(1,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU)-XHH1
  525. ELSE
  526. CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,2,XNOR1)
  527. XJ00(1,1) = -XXCCYY(XNOR1,XNOR1,YOUN,XNU)
  528. & -ISGNTC*RESIST(ICAS)*XHH2
  529. ENDIF
  530. C
  531. R(1) = + FFF(ICAS)
  532. GAMBDA(1)=XZER
  533. GAMBDA(2)=XZER
  534. GAMBDA(3)=XZER
  535. GAMBDA(4)=XZER
  536. GAMBDA(5)=XZER
  537. GAMBDA(6)=XZER
  538. C
  539. XJ0INV = UN/XJ00(1,1)
  540. C
  541. CALL ZERO(D,1,MN)
  542. XH0(1,1) = XJ0INV
  543. C
  544. C On rentre dans les iterations internes
  545. C
  546. DO I=0,ITM
  547. C
  548. C Appel de Broyden
  549. C
  550. CALL BROYDE (QUASIN)
  551. C
  552. GAMBDA(ICAS) = D(1)
  553. C
  554. C Calcul du phi et des termes de la matrice A et A-1
  555. C
  556. SITPRO(1) = SIGI(1) + DSIGI(1)
  557. & - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  558. & + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  559. SITPRO(2) = SIGI(2) + DSIGI(2)
  560. & - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  561. & + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  562. SITPRO(4) = SIGI(4) + DSIGI(4)
  563. C
  564. C Mis a jour des variables d'ecrouissage
  565. C
  566. SIGHP(1) = SIGHP0(1) + XHH1*(GAMBDA(1)-GAMBDA(3))
  567. SIGHP(2) = SIGHP0(2) + XHH1*(GAMBDA(2)-GAMBDA(4))
  568. EPSI(1) = EPSI0(1) + GAMBDA(5)
  569. IF (ICAS.LE.4) THEN
  570. WORK = GAMBDA(ICAS)*RESIST(ICAS)*ISGNTC
  571. ELSE
  572. WORK = XZER
  573. ENDIF
  574. EPSI(2) = EPSI0(2) + GAMBDA(6) + DEGRAD*WORK
  575. C
  576. IF (ICAS.LE.4) THEN
  577. CALL VALEUR(EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  578. &XKISO(1),DUMM,KERRE)
  579. CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  580. &XKISO(2),DUMM,KERRE)
  581. ELSE
  582. IF (ICAS.EQ.5) THEN
  583. CALL VALEUR(EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  584. &XKISO(1),XHH2,KERRE)
  585. CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  586. &XKISO(2),DUMM,KERRE)
  587. ELSE
  588. CALL VALEUR(EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  589. &XKISO(1),DUMM,KERRE)
  590. CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  591. &XKISO(2),XHH2,KERRE)
  592. ENDIF
  593. ENDIF
  594. C
  595. C Calcul des phi
  596. C
  597. IF (ICAS.LE.4) THEN
  598. PHI = ISGNTC*(-UNDEMI*(SITPRO(1) + SITPRO(2)) + SIGHP(IHARD)
  599. & + RESIST(ICAS))
  600. ELSE
  601. PHI = ISGNTC*(-UNDEMI*(SITPRO(1)+SITPRO(2))
  602. & +(XKISOI(ICAS - 4)+(EPSI(ICAS-4) - EPSII(ICAS-4))*XHH2I)
  603. & *RESIST(ICAS))
  604. ENDIF
  605. C
  606. C Calcul de AINV
  607. C
  608. SOMLAM = GAMBDA(1) + GAMBDA(2) + GAMBDA(3) + GAMBDA(4)
  609. & + GAMBDA(5) + GAMBDA(6)
  610. AINV(1) = (DEUX*PHI+GGG*SOMLAM)/(DEUX*PHI+DEUX*GGG*SOMLAM)
  611. AINV(2) = GGG*SOMLAM/(DEUX*PHI + DEUX*GGG*SOMLAM)
  612. AINV(3) = PHI/(PHI + (GGG*SOMLAM))
  613. C
  614. C Calcul des contraintes a l'aide de AINV
  615. C
  616. SIGT(1) = AINV(1)*SITPRO(1) + AINV(2)*SITPRO(2)
  617. SIGT(2) = AINV(2)*SITPRO(1) + AINV(1)*SITPRO(2)
  618. SIGT(4) = AINV(3)*SITPRO(4)
  619. C
  620. C Nouvelles contraintes principales
  621. C
  622. CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  623. C
  624. C
  625. IF (ICAS.LE.4) THEN
  626. CALL CRITPE (SIGTP(IHARD),SIGHP(IHARD),RESIST(ICAS),
  627. & ISGNTC,FFF(ICAS))
  628. R(1)= + FFF(ICAS)
  629. ELSE
  630. FG(1) = (XKISOI(ICAS-4)+(EPSI(ICAS-4)-EPSII(ICAS-4))*XHH2I)
  631. & *RESIST(ICAS)
  632. C
  633. CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,FINT1)
  634. CALL CRITPE (SIGTP(IHARDC),XZER,FG(1),ISGNTC,FINT2)
  635. FFF(ICAS) = MAX (FINT1,FINT2)
  636. C
  637. R(1) = + FFF(ICAS)
  638. & - ISGNTC*RESIST(ICAS)
  639. & * (XKISO(ICAS-4) - XKISOI(ICAS-4)
  640. & -(EPSI(ICAS-4) - EPSII(ICAS-4))*XHH2I)
  641. IF (ABS(EPSI(ICAS-4) - EPSII(ICAS-4)).LT.EPSILO2) THEN
  642. XHH2I=XHH2
  643. ELSE
  644. XHH2I = (XKISO(ICAS-4) - XKISOI(ICAS-4))/
  645. & (EPSI(ICAS-4) - EPSII(ICAS-4))
  646. ENDIF
  647. XKISOI(ICAS-4) = XKISO(ICAS-4)
  648. EPSII(ICAS-4) = EPSI(ICAS-4)
  649. ENDIF
  650. C
  651. C Test de convergence
  652. C
  653. C IF (ABS(FFF(ICAS)).LT.EPSIL(ICAS)) GOTO 1020
  654. IF ((ABS(FFF(ICAS)).LT.EPSIL(ICAS))
  655. & .AND.(ABS(R(1)).LT.EPSIL(ICAS))) GOTO 1020
  656. ENDDO
  657. C
  658. C Pas de convergence
  659. C
  660. SEGSUP,QUASIN
  661. GOTO 6000
  662. C
  663. 1020 CONTINUE
  664. SEGSUP,QUASIN
  665. C
  666. IF (GAMBDA(ICAS).LT.EPSILO2) THEN
  667. WRITE (IOIMP,*)'BETOCY: Multiplicateurs plastiques < 0'
  668. KERRE = 2
  669. RETURN
  670. ELSE
  671. CONTINUE
  672. ENDIF
  673. C
  674. C On verifie qu'on ne viole pas les autres criteres
  675. C
  676. IF (IHARD.EQ.2) THEN
  677. IHARDC=1
  678. ELSE
  679. IHARDC=2
  680. ENDIF
  681. C
  682. IF (ICAS.LE.2) THEN
  683. FG(1) = XKISO(1)*RESIST(5)
  684. FG(2) = XKISO(2)*RESIST(6)
  685. ELSE
  686. IF (ICAS.LE.4) THEN
  687. FG(1) = XKISO(2)*RESIST(6)
  688. FG(2) = XKISO(1)*RESIST(5)
  689. ELSE
  690. IF (ICAS.EQ.5) THEN
  691. FG(1) = XKISO(1)*RESIST(5)
  692. FG(2) = XKISO(2)*RESIST(6)
  693. ELSE
  694. FG(1) = XKISO(2)*RESIST(6)
  695. FG(2) = XKISO(1)*RESIST(5)
  696. ENDIF
  697. ENDIF
  698. ENDIF
  699. C
  700. C Test GS
  701. C
  702. IF (ICAS.GE.5) THEN
  703. C Regles de suivi si on est sur la grande surface
  704. C 1- Direction de l'ecoulement
  705. SIGHP(IHARD) = FG(1) -RESIST(2*ICAS-10+IHARD)
  706. C 2- Direction perpendiculaire a l'ecoulement
  707. C (si on se trouve dans le coin)
  708. XLIMIT = FG(1) -RESIST(2*ICAS-10+IHARDC)
  709. IF (ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  710. SIGHP(IHARDC) = XLIMIT
  711. KONTAC(2*ICAS-10+IHARDC)=1
  712. ENDIF
  713. ENDIF
  714. C
  715. C
  716. C On teste les autres mecanismes
  717. C
  718. C Test des contacts lateraux
  719. C Meme signe signtc
  720. IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  721. ICASB = ICAS + 1
  722. ELSE
  723. IF ((ICAS.EQ.2).OR.(ICAS.EQ.4)) THEN
  724. ICASB = ICAS -1
  725. ELSE
  726. ICASB = 2*ICAS-10+IHARDC
  727. ENDIF
  728. ENDIF
  729. IF (KONTAC(ICASB).EQ.0) THEN
  730. CALL CRITPE(SIGTP(IHARDC), SIGHP(IHARDC), RESIST(ICASB),
  731. &ISGNTC,XFFFS)
  732. IF (((XFFFS.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  733. &((XFFFS.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  734. IHARDB = IHARDC
  735. ISGNTCB = ISGNTC
  736. GOTO 1100
  737. ELSE
  738. IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  739. &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  740. KONTAC(ICASB)=1
  741. ENDIF
  742. ENDIF
  743. ELSE
  744. CALL CRITPE(SIGTP(IHARDC), XZER, FG(1), ISGNTC,XFFFS)
  745. IF (ICASB.LE.2) THEN
  746. ICASB = 5
  747. ELSE
  748. ICASB = 6
  749. ENDIF
  750. IF (XFFFS.GE.EPSIL(ICASB)) THEN
  751. IHARDB = IHARDC
  752. ISGNTCB = ISGNTC
  753. GOTO 1100
  754. ENDIF
  755. ENDIF
  756. C Autre signe signtc
  757. IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  758. ICASB = ICAS + 1
  759. ELSE
  760. IF ((ICAS.EQ.2).OR.(ICAS.EQ.4)) THEN
  761. ICASB = ICAS - 1
  762. ELSE
  763. ICASB = 2*ICAS-10+IHARDC
  764. ENDIF
  765. ENDIF
  766. IF (ICASB.LE.2) THEN
  767. ICASB = ICASB +2
  768. ELSE
  769. ICASB = ICASB -2
  770. ENDIF
  771. IF (KONTAC(ICASB).EQ.0) THEN
  772. CALL CRITPE(SIGTP(IHARDC), SIGHP(IHARDC), RESIST(ICASB),
  773. &-ISGNTC,XFFFS)
  774. IF (((XFFFS.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  775. &((XFFFS.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  776. IHARDB = IHARDC
  777. ISGNTCB = -ISGNTC
  778. GOTO 1100
  779. ELSE
  780. IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  781. &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  782. KONTAC(ICASB)=1
  783. ENDIF
  784. ENDIF
  785. ELSE
  786. CALL CRITPE(SIGTP(IHARDC), XZER, FG(2), -ISGNTC,XFFFS)
  787. IF (ICASB.LE.2) THEN
  788. ICASB = 5
  789. ELSE
  790. ICASB = 6
  791. ENDIF
  792. IF (XFFFS.GE.EPSIL(ICASB)) THEN
  793. IHARDB = IHARDC
  794. ISGNTCB = -ISGNTC
  795. GOTO 1100
  796. ENDIF
  797. ENDIF
  798. C Test du contact frontal
  799. IF (ICAS.LE.4) THEN
  800. CALL CRITPE(SIGTP(IHARD), XZER, FG(1), ISGNTC,XFFFS)
  801. IF (((XFFFS.GE.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  802. &((XFFFS.GE.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  803. GOTO 1040
  804. ELSE
  805. IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  806. &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  807. KONTAC(ICAS) = 1
  808. ENDIF
  809. ENDIF
  810. ENDIF
  811. C
  812. C Il faut faire suivre la petite surface si la grande se retrecit trop vite
  813. C
  814. IF (ICAS.LE.4) THEN
  815. IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  816. ICASB = ICAS +1
  817. ELSE
  818. ICASB = ICAS -1
  819. ENDIF
  820. XLIMIT = FG(1) -RESIST(ICASB)
  821. IF (ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  822. SIGHP(IHARDC) = XLIMIT
  823. KONTAC(ICASB) = 1
  824. ENDIF
  825. IF (ICASB.LE.2) THEN
  826. ICASB = ICASB +2
  827. ELSE
  828. ICASB = ICASB -2
  829. ENDIF
  830. XLIMIT = FG(2) -RESIST(ICASB)
  831. IF (-ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  832. SIGHP(IHARDC) = XLIMIT
  833. KONTAC(ICASB) = 1
  834. ENDIF
  835. ENDIF
  836. C
  837. C Mis a jour des indicateurs de contact
  838. C
  839. IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).GT.
  840. &(RESIST(1)-RESIST(3))) THEN
  841. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) KONTAC(ICAS+2)=0
  842. IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) KONTAC(ICAS-2)=0
  843. ELSE
  844. KONTAC(1) = 1
  845. KONTAC(2) = 1
  846. KONTAC(3) = 1
  847. KONTAC(4) = 1
  848. ENDIF
  849. C
  850. ICASB = ICAS
  851. GOTO 9999
  852. C
  853. C Cas ou on trouve le contact frontal
  854. C
  855. 1040 CONTINUE
  856. SIGHP(1) = SIGHP0(1)
  857. SIGHP(2) = SIGHP0(2)
  858. SIGT(1) = SIGI(1) + DSIGI(1)
  859. SIGT(2) = SIGI(2) + DSIGI(2)
  860. SIGT(4) = SIGI(4) + DSIGI(4)
  861. GOTO 3000
  862. C
  863. C Cas ou on viole un critere lateral
  864. C
  865. 1100 CONTINUE
  866. SIGHP(1) = SIGHP0(1)
  867. SIGHP(2) = SIGHP0(2)
  868. SIGT(1) = SIGI(1) + DSIGI(1)
  869. SIGT(2) = SIGI(2) + DSIGI(2)
  870. SIGT(4) = SIGI(4) + DSIGI(4)
  871. IF (ICAS.EQ.ICASB) THEN
  872. IF (ICAS.GE.5) THEN
  873. GOTO 1000
  874. ELSE
  875. KERRE = 2
  876. RETURN
  877. ENDIF
  878. ELSE
  879. GOTO 4000
  880. ENDIF
  881. C================================================================
  882. C CAS 2
  883. C Un seul mecanisme de chaque surface
  884. C mais entree en contact des 2 mecanismes
  885. C================================================================
  886. 3000 CONTINUE
  887. C
  888. CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  889. C
  890. C
  891. C Calcul de la normale a l'estimation elastique et du jacobien initial
  892. C
  893. IF (IHARD.EQ.1) THEN
  894. IF (SIGTP(1).GE.SIGTP(2)) THEN
  895. ISGN2 = 1
  896. ELSE
  897. ISGN2 = -1
  898. ENDIF
  899. ELSE
  900. IF (IHARD.EQ.2) THEN
  901. IF (SIGTP(2).GE.SIGTP(1)) THEN
  902. ISGN2 = 1
  903. ELSE
  904. ISGN2 = -1
  905. ENDIF
  906. ENDIF
  907. ENDIF
  908. C
  909. 3010 CONTINUE
  910. C
  911. C Definition du segment
  912. C
  913. ITM = 30
  914. MN = 2
  915. SEGINI,QUASIN
  916. IT = -1
  917. CALL ZERO(D,1,MN)
  918. C
  919. C Calcul de xhh2 et xkiso
  920. C
  921. IF (ICAS.LE.2) THEN
  922. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  923. &XKISO(1),XHH2,KERRE)
  924. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  925. &XKISO(2),DUMM,KERRE)
  926. ELSE
  927. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  928. &XKISO(1),DUMM,KERRE)
  929. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  930. &XKISO(2),XHH2,KERRE)
  931. ENDIF
  932. EPSI(1) = EPSI0(1)
  933. XKISO0(1) = XKISO(1)
  934. EPSI(2) = EPSI0(2)
  935. XKISO0(2) = XKISO(2)
  936. XHH2I = XHH2
  937. EPSII(1) = EPSI0(1)
  938. XKISOI(1) = XKISO(1)
  939. EPSII(2) = EPSI0(2)
  940. XKISOI(2) = XKISO(2)
  941. C
  942. DO 3011,IE1=1,6
  943. FFF(IE1) = XZER
  944. 3011 CONTINUE
  945. C
  946. IF (ISGNTC.GT.0) THEN
  947. FG(1) = XKISO(1)*RESIST(5)
  948. ELSE
  949. FG(1) = XKISO(2)*RESIST(6)
  950. ENDIF
  951. CALL CRITPE (SIGTP(IHARD),SIGHP0(IHARD)
  952. & ,RESIST(ICAS),ISGNTC,FFF(ICAS))
  953. CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC, XFFFS)
  954. C
  955. CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,1,XNOR1)
  956. XJ00(1,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU)-XHH1
  957. XJ00(2,2) = -XXYY(XNOR1,DSIGI)
  958. XJ00(1,2) = -XXYY(XNOR1,DSIGI)
  959. XJ00(2,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU)
  960. IF (ISGNTC.LT.0) THEN
  961. XJ00(2,1)=XJ00(2,1)+RESIST(6)*XHH2*DEGRAD*ISGNTC*RESIST(ICAS)
  962. ENDIF
  963. R(1) = FFF(ICAS)
  964. R(2) = XFFFS
  965. GAMBDA(1)=XZER
  966. GAMBDA(2)=XZER
  967. GAMBDA(3)=XZER
  968. GAMBDA(4)=XZER
  969. GAMBDA(5)=XZER
  970. GAMBDA(6)=XZER
  971. C
  972. C Inverse du jacobien
  973. C
  974. XH0(1,1) = XJ00(1,1)
  975. XH0(1,2) = XJ00(1,2)
  976. XH0(2,1) = XJ00(2,1)
  977. XH0(2,2) = XJ00(2,2)
  978. CALL INVALM (XH0, MN , MN , IRD, 1.D-12)
  979. IF (IRD.NE.0) THEN
  980. SEGSUP,QUASIN
  981. GOTO 6000
  982. ENDIF
  983. C
  984. C On rentre dans les iterations internes
  985. C
  986. DO I=0,ITM
  987. C
  988. C Appel de Broyden
  989. C
  990. CALL BROYDE (QUASIN)
  991. C
  992. GAMBDA(ICAS) = D(1)
  993. SPLITT = D(2)
  994. C
  995. C Calcul du phi et des termes de la matrice A et A-1
  996. C
  997. SITPRO(1) = SIGI(1) + ((UN-SPLITT)*DSIGI(1))
  998. & - UNDEMI*(GAMBDA(1) + GAMBDA(2))*CPI
  999. & + UNDEMI*(GAMBDA(3) + GAMBDA(4))*CPI
  1000. SITPRO(2) = SIGI(2) + ((UN-SPLITT)*DSIGI(2))
  1001. & - UNDEMI*(GAMBDA(1) + GAMBDA(2))*CPI
  1002. & + UNDEMI*(GAMBDA(3) + GAMBDA(4))*CPI
  1003. SITPRO(4) = SIGI(4) + (UN-SPLITT)*DSIGI(4)
  1004. C
  1005. SIGHP(1) = SIGHP0(1) + XHH1*(GAMBDA(1)-GAMBDA(3))
  1006. SIGHP(2) = SIGHP0(2) + XHH1*(GAMBDA(2)-GAMBDA(4))
  1007. EPSI(1) = EPSI0(1)
  1008. WORK = GAMBDA(ICAS)*ISGNTC*RESIST(ICAS)
  1009. EPSI(2) = EPSI0(2) + DEGRAD*WORK
  1010. CALL VALEUR(EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  1011. &XKISO(2),XHH2,KERRE)
  1012. C
  1013. PHI = ISGNTC*(-UNDEMI*(SITPRO(1)+SITPRO(2))+SIGHP(IHARD)
  1014. & + RESIST(ICAS))
  1015. C
  1016. C
  1017. SOMLAM = GAMBDA(1) + GAMBDA(2) + GAMBDA(3) + GAMBDA(4)
  1018. AINV(1) = (DEUX*PHI+GGG*SOMLAM)/(DEUX*PHI+DEUX*GGG*SOMLAM)
  1019. AINV(2) = GGG*SOMLAM/(DEUX*PHI + DEUX*GGG*SOMLAM)
  1020. AINV(3) = PHI/(PHI + GGG*SOMLAM)
  1021. C
  1022. C Calcul des contraintes a l'aide de AINV
  1023. C
  1024. SIGT(1) = AINV(1)*SITPRO(1)+AINV(2)*SITPRO(2)
  1025. SIGT(2) = AINV(2)*SITPRO(1)+AINV(1)*SITPRO(2)
  1026. SIGT(4) = AINV(3)*SITPRO(4)
  1027. C
  1028. C Nouvelles contraintes principales
  1029. C
  1030. CALL CONPRI (SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1031. C
  1032. NMEC = 0
  1033. C
  1034. IF (ISGNTC.GT.0) THEN
  1035. FG(1) = XKISO(1)*RESIST(5)
  1036. ELSE
  1037. FG(1) = (XKISOI(2) + XHH2I*(EPSI(2)-EPSII(2)))*RESIST(6)
  1038. ENDIF
  1039. CALL CRITPE (SIGTP(IHARD),SIGHP(IHARD)
  1040. & ,RESIST(ICAS),ISGNTC,FFF(ICAS))
  1041. R(1) = + FFF(ICAS)
  1042. IF (ABS(FFF(ICAS)).GE.EPSIL(ICAS)) THEN
  1043. NMEC = NMEC + 1
  1044. ENDIF
  1045. CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC, XFFFS)
  1046. R(2) = XFFFS
  1047. IF (ICAS.GE.3) THEN
  1048. R(2) = R(2) + RESIST(6)*(XKISO(2) - XKISOI(2)
  1049. & - (EPSI(2) - EPSII(2))*XHH2I)
  1050. IF (ABS(EPSI(2) - EPSII(2)).LT.EPSILO2) THEN
  1051. XHH2I=XHH2
  1052. ELSE
  1053. XHH2I = (XKISO(2) - XKISOI(2))/
  1054. & (EPSI(2) - EPSII(2))
  1055. ENDIF
  1056. EPSII(2) = EPSI(2)
  1057. XKISOI(2) = XKISO(2)
  1058. ENDIF
  1059. IF (((ABS(XFFFS).GE.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  1060. &((ABS(XFFFS).GE.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  1061. NMEC = NMEC + 1
  1062. ENDIF
  1063. C
  1064. C Test de convergence
  1065. C
  1066. IF (NMEC.EQ.0) GOTO 3020
  1067. ENDDO
  1068. C
  1069. C Pas de convergence
  1070. C
  1071. SEGSUP,QUASIN
  1072. GOTO 6000
  1073. C
  1074. 3020 CONTINUE
  1075. C
  1076. C On verifie que la convergence ne s'est pas faite sur un point debile
  1077. C
  1078. SEGSUP,QUASIN
  1079. IF (((SPLITT-1.).LE.EPSILO).AND.(SPLITT.GE.(-1.*EPSILO))) THEN
  1080. C
  1081. IF (IHARD.EQ.2) THEN
  1082. IHARDC=1
  1083. ELSE
  1084. IHARDC=2
  1085. ENDIF
  1086. IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) ICASC = ICAS + 1
  1087. IF ((ICAS.EQ.2).OR.(ICAS.EQ.4)) ICASC = ICAS - 1
  1088. C
  1089. IF (ICAS.LE.2) THEN
  1090. FG(1) = XKISO(1)*RESIST(5)
  1091. FG(2) = XKISO(2)*RESIST(6)
  1092. ELSE
  1093. FG(1) = XKISO(2)*RESIST(6)
  1094. FG(2) = XKISO(1)*RESIST(5)
  1095. ENDIF
  1096. C
  1097. C Test des contacts lateraux
  1098. C Meme signe signtc
  1099. IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  1100. ICASB = ICAS + 1
  1101. ELSE
  1102. ICASB = ICAS -1
  1103. ENDIF
  1104. IF (KONTAC(ICASB).EQ.0) THEN
  1105. CALL CRITPE(SIGTP(IHARDC), SIGHP(IHARDC),
  1106. & RESIST(ICASB), ISGNTC,XFFFS)
  1107. IF (((XFFFS.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1108. &((XFFFS.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1109. IHARDB = IHARDC
  1110. ISGNTCB = ISGNTC
  1111. GOTO 3100
  1112. ELSE
  1113. IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1114. &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1115. KONTAC(ICASB)=1
  1116. ENDIF
  1117. ENDIF
  1118. ELSE
  1119. CALL CRITPE(SIGTP(IHARDC), XZER, FG(1), ISGNTC,XFFFS)
  1120. IF (ICASB.LE.2) THEN
  1121. ICASB = 5
  1122. ELSE
  1123. ICASB = 6
  1124. ENDIF
  1125. IF (XFFFS.GE.EPSIL(ICASB)) THEN
  1126. IHARDB = IHARDC
  1127. ISGNTCB = -ISGNTC
  1128. GOTO 3100
  1129. ENDIF
  1130. ENDIF
  1131. C Autre signe signtc
  1132. IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  1133. ICASB = ICAS + 1
  1134. ELSE
  1135. ICASB = ICAS -1
  1136. ENDIF
  1137. IF (ICASB.LE.2) THEN
  1138. ICASB = ICASB +2
  1139. ELSE
  1140. ICASB = ICASB -2
  1141. ENDIF
  1142. IF (KONTAC(ICASB).EQ.0) THEN
  1143. CALL CRITPE(SIGTP(IHARDC), SIGHP(IHARDC),RESIST(ICASB),
  1144. & -ISGNTC,XFFFS)
  1145. IF (((XFFFS.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1146. &((XFFFS.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1147. IHARDB = IHARDC
  1148. ISGNTCB = ISGNTC
  1149. GOTO 3100
  1150. ELSE
  1151. IF (((ABS(XFFFS).LT.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1152. &((ABS(XFFFS).LT.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1153. KONTAC(ICASB)=1
  1154. ENDIF
  1155. ENDIF
  1156. ELSE
  1157. CALL CRITPE(SIGTP(IHARDC), XZER, FG(2), -ISGNTC,XFFFS)
  1158. IF (ICASB.LE.2) THEN
  1159. ICASB = 5
  1160. ELSE
  1161. ICASB = 6
  1162. ENDIF
  1163. IF (XFFFS.GE.EPSIL(ICASB)) THEN
  1164. IHARDB = IHARDC
  1165. ISGNTCB = -ISGNTC
  1166. GOTO 3100
  1167. ENDIF
  1168. ENDIF
  1169. C
  1170. C Il faut faire suivre la petite surface si la grande se retrecit trop vite
  1171. C
  1172. IF (ICAS.LE.4) THEN
  1173. IF ((ICAS.EQ.1).OR.(ICAS.EQ.3)) THEN
  1174. ICASB = ICAS +1
  1175. ELSE
  1176. ICASB = ICAS -1
  1177. ENDIF
  1178. XLIMIT = FG(1) -RESIST(ICASB)
  1179. IF (ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  1180. SIGHP(IHARDC) = XLIMIT
  1181. KONTAC(ICASB) = 1
  1182. ENDIF
  1183. IF (ICASB.LE.2) THEN
  1184. ICASB = ICASB +2
  1185. ELSE
  1186. ICASB = ICASB -2
  1187. ENDIF
  1188. XLIMIT = FG(2) -RESIST(ICASB)
  1189. IF (-ISGNTC*(SIGHP(IHARDC)-XLIMIT).GT.XZER) THEN
  1190. SIGHP(IHARDC) = XLIMIT
  1191. KONTAC(ICASB) = 1
  1192. ENDIF
  1193. ENDIF
  1194. C
  1195. C
  1196. C Mis a jour des indicateurs de contact
  1197. C
  1198. KONTAC(ICAS) = 1
  1199. IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).GT.
  1200. &(RESIST(1)-RESIST(3))) THEN
  1201. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) KONTAC(ICAS+2)=0
  1202. IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) KONTAC(ICAS-2)=0
  1203. ELSE
  1204. KONTAC(1) = 1
  1205. KONTAC(2) = 1
  1206. KONTAC(3) = 1
  1207. KONTAC(4) = 1
  1208. ENDIF
  1209. SIGI(1) = SIGT(1)
  1210. SIGI(2) = SIGT(2)
  1211. SIGI(4) = SIGT(4)
  1212. DSIGI(1) = DSIGI(1)*SPLITT
  1213. DSIGI(2) = DSIGI(2)*SPLITT
  1214. DSIGI(4) = DSIGI(4)*SPLITT
  1215. SIGT(1) = SIGI(1) + DSIGI(1)
  1216. SIGT(2) = SIGI(2) + DSIGI(2)
  1217. SIGT(4) = SIGI(4) + DSIGI(4)
  1218. SIGHP0(1) = SIGHP(1)
  1219. SIGHP0(2) = SIGHP(2)
  1220. EPSI0(1) = EPSI(1)
  1221. EPSI0(2) = EPSI(2)
  1222. CPHI0 = CPHI
  1223. SPHI0 = SPHI
  1224. NPASS2 = 0
  1225. IF (ICAS.LE.2) THEN
  1226. ICAS = 5
  1227. ELSE
  1228. ICAS = 6
  1229. ENDIF
  1230. GOTO 1000
  1231. C
  1232. C
  1233. ELSE
  1234. GOTO 6000
  1235. ENDIF
  1236. C
  1237. C Cas ou on viole un critere lateral
  1238. C
  1239. 3100 CONTINUE
  1240. SIGHP(1) = SIGHP0(1)
  1241. SIGHP(2) = SIGHP0(2)
  1242. SIGT(1) = SIGI(1) + DSIGI(1)
  1243. SIGT(2) = SIGI(2) + DSIGI(2)
  1244. SIGT(4) = SIGI(4) + DSIGI(4)
  1245. CALL CONPRI ( SIGT ,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1246. CALL VALEUR(EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1),
  1247. &DUMM,KERRE)
  1248. CALL VALEUR(EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2),
  1249. &DUMM,KERRE)
  1250. CALL CRITPE(SIGTP(IHARD), SIGHP(IHARD), RESIST(ICAS),
  1251. & ISGNTC, FFF(ICAS))
  1252. IF (ICASB.LE.4) THEN
  1253. CALL CRITPE(SIGTP(IHARDB), SIGHP(IHARDB), RESIST(ICASB),
  1254. & ISGNTCB, FFF(ICASB))
  1255. ELSE
  1256. FG(1) = XKISO(ICASB-4)*RESIST(ICASB)
  1257. CALL CRITPE (SIGTP(IHARDB), XZER, FG(1),
  1258. & ISGNTCB,FFF(ICASB))
  1259. ENDIF
  1260. IF (ICAS.LE.2) THEN
  1261. FG(1) = XKISO(1)*RESIST(5)
  1262. ELSE
  1263. FG(1) = XKISO(2)*RESIST(6)
  1264. ENDIF
  1265. CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC, XFFFS)
  1266. GOTO 5000
  1267. C================================================================
  1268. C CAS 4
  1269. C Deux mecanismes
  1270. C Petite et Grande surfaces (18/07)
  1271. C================================================================
  1272. 4000 CONTINUE
  1273. C
  1274. CALL CONPRI ( SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1275. C
  1276. C
  1277. C Calcul de la normale a l'estimation elastique et du jacobien initial
  1278. C
  1279. IF (IHARD.EQ.1) THEN
  1280. IF (SIGTP(1).GE.SIGTP(2)) THEN
  1281. ISGN2 = 1
  1282. ELSE
  1283. ISGN2 = -1
  1284. ENDIF
  1285. ELSE
  1286. IF (IHARD.EQ.2) THEN
  1287. IF (SIGTP(2).GE.SIGTP(1)) THEN
  1288. ISGN2 = 1
  1289. ELSE
  1290. ISGN2 = -1
  1291. ENDIF
  1292. ENDIF
  1293. ENDIF
  1294. C
  1295. IF (IHARDB.EQ.1) THEN
  1296. IF (SIGTP(1).GE.SIGTP(2)) THEN
  1297. ISGN2B = 1
  1298. ELSE
  1299. ISGN2B = -1
  1300. ENDIF
  1301. ELSE
  1302. IF (IHARDB.EQ.2) THEN
  1303. IF (SIGTP(1).GE.SIGTP(2)) THEN
  1304. ISGN2B = -1
  1305. ELSE
  1306. ISGN2B = 1
  1307. ENDIF
  1308. ENDIF
  1309. ENDIF
  1310. C
  1311. C
  1312. C
  1313. 4010 CONTINUE
  1314. C
  1315. C Definition du segment
  1316. C
  1317. ITM = 30
  1318. MN = 2
  1319. SEGINI,QUASIN
  1320. IT = -1
  1321. CALL ZERO(D,1,MN)
  1322. C
  1323. C XKISO et XHH2
  1324. C
  1325. IF (ICAS.LE.4) THEN
  1326. IF (ICASB.LE.4) THEN
  1327. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  1328. &XKISO(1),DUMM,KERRE)
  1329. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  1330. &XKISO(2),DUMM,KERRE)
  1331. ELSE
  1332. IF (ICASB.EQ.5) THEN
  1333. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1334. & ,XHH2B,KERRE)
  1335. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1336. & ,DUMM,KERRE)
  1337. ELSE
  1338. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1339. & ,DUMM,KERRE)
  1340. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1341. & ,XHH2B,KERRE)
  1342. ENDIF
  1343. ENDIF
  1344. ELSE
  1345. IF (ICASB.LE.4) THEN
  1346. IF (ICAS.EQ.5) THEN
  1347. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1348. & ,XHH2,KERRE)
  1349. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1350. & ,DUMM,KERRE)
  1351. ELSE
  1352. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1353. & ,DUMM,KERRE)
  1354. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1355. & ,XHH2,KERRE)
  1356. ENDIF
  1357. ELSE
  1358. IF (ICASB.EQ.5) THEN
  1359. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1360. & ,XHH2B,KERRE)
  1361. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1362. & ,XHH2,KERRE)
  1363. ELSE
  1364. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1365. & ,XHH2,KERRE)
  1366. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1367. & ,XHH2B,KERRE)
  1368. ENDIF
  1369. ENDIF
  1370. ENDIF
  1371. EPSI(1) = EPSI0(1)
  1372. EPSI(2) = EPSI0(2)
  1373. XHH2I = XHH2
  1374. XHH2IB = XHH2B
  1375. XKISO0(1) = XKISO(1)
  1376. XKISO0(2) = XKISO(2)
  1377. EPSII(1) = EPSI0(1)
  1378. EPSII(2) = EPSI0(2)
  1379. XKISOI(1) = XKISO(1)
  1380. XKISOI(2) = XKISO(2)
  1381. C
  1382. DO 4011,IE1=1,6
  1383. FFF(IE1) = XZER
  1384. 4011 CONTINUE
  1385. C
  1386. IF (ICAS.LE.4) THEN
  1387. CALL CRITPE (SIGTP(IHARD),SIGHP0(IHARD),RESIST(ICAS),
  1388. & ISGNTC,FFF(ICAS))
  1389. ELSE
  1390. FG(1) = XKISOI(ICAS-4)*RESIST(ICAS)
  1391. CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,FFF(ICAS))
  1392. ENDIF
  1393. C
  1394. IF (ICASB.LE.4) THEN
  1395. CALL CRITPE (SIGTP(IHARDB),SIGHP0(IHARDB),RESIST(ICASB),
  1396. & ISGNTCB,FFF(ICASB))
  1397. ELSE
  1398. FG(1) = XKISOI(ICASB-4)*RESIST(ICASB)
  1399. CALL CRITPE (SIGTP(IHARDB),XZER,FG(1),ISGNTCB,FFF(ICASB))
  1400. ENDIF
  1401. C
  1402. C Jcobien initial
  1403. C
  1404. CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,1,XNOR1)
  1405. CALL NORMBL(SIGT,RESIST(ICASB),ISGNTCB,ISGN2B,1,XNOR2)
  1406. C
  1407. XJ00(1,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU)
  1408. XJ00(1,2) = -XXCCYY(XNOR1, XNOR2, YOUN, XNU)
  1409. XJ00(2,1) = -XXCCYY(XNOR2, XNOR1, YOUN, XNU)
  1410. XJ00(2,2) = -XXCCYY(XNOR2, XNOR2, YOUN, XNU)
  1411. C
  1412. C Premier Mecanisme
  1413. C
  1414. IF (ICAS.LE.4) THEN
  1415. XJ00(1,1) = XJ00(1,1) - XHH1
  1416. ELSE
  1417. XJ00(1,1) = XJ00(1,1)-ISGNTC*RESIST(ICAS)*XHH2
  1418. IF ((ICASB.LE.4).AND.(ICAS.EQ.6)) THEN
  1419. XJ00(1,2) = XJ00(1,2)
  1420. & +RESIST(6)*XHH2*DEGRAD*ISGNTCB*RESIST(ICASB)
  1421. ENDIF
  1422. ENDIF
  1423. C
  1424. C Second Mecanisme
  1425. C
  1426. IF (ICASB.LE.4) THEN
  1427. XJ00(2,2) = XJ00(2,2) - XHH1
  1428. ELSE
  1429. XJ00(2,2) = XJ00(2,2)-ISGNTCB*RESIST(ICASB)*XHH2B
  1430. IF ((ICAS.LE.4).AND.(ICASB.EQ.6)) THEN
  1431. XJ00(2,1) = XJ00(2,1)
  1432. & +RESIST(6)*XHH2B*DEGRAD*ISGNTC*RESIST(ICAS)
  1433. ENDIF
  1434. ENDIF
  1435. C
  1436. C
  1437. R(1) = + FFF(ICAS)
  1438. R(2) = + FFF(ICASB)
  1439. GAMBDA(1) = XZER
  1440. GAMBDA(2) = XZER
  1441. GAMBDA(3) = XZER
  1442. GAMBDA(4) = XZER
  1443. GAMBDA(5) = XZER
  1444. GAMBDA(6) = XZER
  1445. C
  1446. C Inverse du jacobien
  1447. C
  1448. XH0(1,1) = XJ00(1,1)
  1449. XH0(1,2) = XJ00(1,2)
  1450. XH0(2,1) = XJ00(2,1)
  1451. XH0(2,2) = XJ00(2,2)
  1452. CALL INVALM (XH0, MN , MN , IRD, 1.D-12)
  1453. IF (IRD.NE.0) THEN
  1454. SEGSUP,QUASIN
  1455. GOTO 6000
  1456. ENDIF
  1457. C
  1458. C On rentre dans les iterations internes
  1459. C
  1460. DO I=0,ITM
  1461. C
  1462. C Appel de Broyden
  1463. C
  1464. CALL BROYDE (QUASIN)
  1465. C
  1466. GAMBDA(ICAS) = D(1)
  1467. GAMBDA(ICASB) = D(2)
  1468. C
  1469. C Calcul du phi et des termes de la matrice A et A-1
  1470. C
  1471. SITPRO(1) = SIGI(1) + DSIGI(1)
  1472. & - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  1473. & + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  1474. SITPRO(2) = SIGI(2) + DSIGI(2)
  1475. & - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  1476. & + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  1477. SITPRO(4) = SIGI(4) + DSIGI(4)
  1478. C
  1479. C PHI et AINV
  1480. C
  1481. SIGHP(1) = SIGHP0(1) + XHH1*(GAMBDA(1)-GAMBDA(3))
  1482. SIGHP(2) = SIGHP0(2) + XHH1*(GAMBDA(2)-GAMBDA(4))
  1483. EPSI(1) = EPSI0(1) + GAMBDA(5)
  1484. IF (ICAS.LE.4) THEN
  1485. WORK = GAMBDA(ICAS)*RESIST(ICAS)*ISGNTC
  1486. ELSE
  1487. WORK = XZER
  1488. ENDIF
  1489. IF (ICASB.LE.4) THEN
  1490. WORK = WORK + GAMBDA(ICASB)*RESIST(ICASB)*ISGNTCB
  1491. ENDIF
  1492. EPSI(2) = EPSI0(2) + GAMBDA(6) + DEGRAD*WORK
  1493. C
  1494. IF (ICAS.LE.4) THEN
  1495. IF (ICASB.LE.4) THEN
  1496. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  1497. &XKISO(1),DUMM,KERRE)
  1498. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  1499. &XKISO(2),DUMM,KERRE)
  1500. ELSE
  1501. IF (ICASB.EQ.5) THEN
  1502. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1503. & ,XHH2B,KERRE)
  1504. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1505. & ,DUMM,KERRE)
  1506. ELSE
  1507. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1508. & ,DUMM,KERRE)
  1509. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1510. & ,XHH2B,KERRE)
  1511. ENDIF
  1512. ENDIF
  1513. ELSE
  1514. IF (ICASB.LE.4) THEN
  1515. IF (ICAS.EQ.5) THEN
  1516. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1517. & ,XHH2,KERRE)
  1518. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1519. & ,DUMM,KERRE)
  1520. ELSE
  1521. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1522. & ,DUMM,KERRE)
  1523. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1524. & ,XHH2,KERRE)
  1525. ENDIF
  1526. ELSE
  1527. IF (ICASB.EQ.5) THEN
  1528. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1529. & ,XHH2B,KERRE)
  1530. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1531. & ,XHH2,KERRE)
  1532. ELSE
  1533. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1534. & ,XHH2,KERRE)
  1535. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1536. & ,XHH2B,KERRE)
  1537. ENDIF
  1538. ENDIF
  1539. ENDIF
  1540. C
  1541. C Calcul des phi
  1542. C
  1543. IF (ICAS.LE.4) THEN
  1544. PHI1=ISGNTC*(-UNDEMI*(SITPRO(1) + SITPRO(2))
  1545. & + SIGHP(IHARD) + RESIST(ICAS))
  1546. ELSE
  1547. PHI1= ISGNTC*(-UNDEMI*(SITPRO(1)+SITPRO(2))
  1548. & +(XKISOI(ICAS-4)+XHH2I*(EPSI(ICAS-4)-EPSII(ICAS-4)))
  1549. & *RESIST(ICAS))
  1550. ENDIF
  1551. C
  1552. IF (ICASB.LE.4) THEN
  1553. PHI2=ISGNTCB*(-UNDEMI*(SITPRO(1) + SITPRO(2))
  1554. & + SIGHP(IHARDB)+ RESIST(ICASB))
  1555. ELSE
  1556. PHI2=ISGNTCB*(-UNDEMI*(SITPRO(1)+SITPRO(2))
  1557. & +(XKISOI(ICASB-4)+XHH2IB*(EPSI(ICASB-4)-EPSII(ICASB-4)))
  1558. & *RESIST(ICASB))
  1559. ENDIF
  1560. C
  1561. C Test sur phi1 et phi2
  1562. C
  1563. IF ((ABS(PHI2)).GT.(ABS(EPSILO*RESIST(ICASB)))) THEN
  1564. SOMLAM = GAMBDA(ICAS) + GAMBDA(ICASB)*PHI1/PHI2
  1565. PHI = PHI1
  1566. ELSE
  1567. IF ((ABS(PHI1)).GT.(ABS(EPSILO*RESIST(ICAS)))) THEN
  1568. SOMLAM = GAMBDA(ICASB) + GAMBDA(ICAS)*PHI2/PHI1
  1569. PHI = PHI2
  1570. ELSE
  1571. SOMLAM = GAMBDA(ICASB) + GAMBDA(ICAS)
  1572. PHI = PHI1
  1573. ENDIF
  1574. ENDIF
  1575. C
  1576. AINV(1) = (DEUX*PHI+GGG*SOMLAM)/(DEUX*PHI+DEUX*GGG*SOMLAM)
  1577. AINV(2) = (GGG*SOMLAM)/(DEUX*PHI + DEUX*GGG*SOMLAM)
  1578. AINV(3) = PHI/(PHI + GGG*SOMLAM)
  1579. C
  1580. C Calcul des contraintes a l'aide de AINV
  1581. C
  1582. SIGT(1) = AINV(1)*SITPRO(1)+AINV(2)*SITPRO(2)
  1583. SIGT(2) = AINV(2)*SITPRO(1)+AINV(1)*SITPRO(2)
  1584. SIGT(4) = AINV(3)*SITPRO(4)
  1585. C
  1586. C Nouvelles contraintes principales
  1587. C
  1588. CALL CONPRI ( SIGT,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1589. C
  1590. C
  1591. IF (ICAS.LE.4) THEN
  1592. CALL CRITPE (SIGTP(IHARD),SIGHP(IHARD),RESIST(ICAS),
  1593. & ISGNTC,FFF(ICAS))
  1594. R(1)= + FFF(ICAS)
  1595. ELSE
  1596. FG(1) = (XKISOI(ICAS-4)+XHH2I*(EPSI(ICAS-4)-EPSII(ICAS-4)))
  1597. & *RESIST(ICAS)
  1598. CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,FFF(ICAS))
  1599. R(1) = + FFF(ICAS)
  1600. & - ISGNTC*RESIST(ICAS)
  1601. & * (XKISO(ICAS-4) - XKISOI(ICAS-4)
  1602. & -(EPSI(ICAS-4) - EPSII(ICAS-4))*XHH2I)
  1603. IF (ABS(EPSI(ICAS-4) - EPSII(ICAS-4)).LT.EPSILO2) THEN
  1604. XHH2I=XHH2
  1605. ELSE
  1606. XHH2I = (XKISO(ICAS-4) - XKISOI(ICAS-4))/
  1607. & (EPSI(ICAS-4) - EPSII(ICAS-4))
  1608. ENDIF
  1609. EPSII(ICAS-4) = EPSI(ICAS-4)
  1610. XKISOI(ICAS-4) = XKISO(ICAS-4)
  1611. ENDIF
  1612. C
  1613. IF (ICASB.LE.4) THEN
  1614. CALL CRITPE (SIGTP(IHARDB),SIGHP(IHARDB),RESIST(ICASB),
  1615. & ISGNTCB,FFF(ICASB))
  1616. R(2)= + FFF(ICASB)
  1617. ELSE
  1618. FG(1) = (XKISOI(ICASB-4)+XHH2IB*(EPSI(ICASB-4) -
  1619. & EPSII(ICASB-4)))*RESIST(ICASB)
  1620. CALL CRITPE (SIGTP(IHARDB),XZER,FG(1),ISGNTCB,FFF(ICASB))
  1621. R(2) = + FFF(ICASB)
  1622. & - ISGNTCB*RESIST(ICASB)
  1623. & * (XKISO(ICASB-4) - XKISOI(ICASB-4)
  1624. & -(EPSI(ICASB-4) - EPSII(ICASB-4))*XHH2IB)
  1625. IF (ABS(EPSI(ICASB-4) - EPSII(ICASB-4)).LT.EPSILO2) THEN
  1626. XHH2IB=XHH2B
  1627. ELSE
  1628. XHH2IB = (XKISO(ICASB-4) - XKISOI(ICASB-4))/
  1629. & (EPSI(ICASB-4) - EPSII(ICASB-4))
  1630. ENDIF
  1631. XKISOI(ICASB-4) = XKISO(ICASB-4)
  1632. EPSII(ICASB-4) = EPSI(ICASB-4)
  1633. ENDIF
  1634. C
  1635. C Test de convergence
  1636. C
  1637. IF ((ABS(FFF(ICAS)).LE.EPSIL(ICAS))
  1638. &.AND.(ABS(FFF(ICASB)).LE.EPSIL(ICASB)).AND.
  1639. &(ABS(R(1)).LE.EPSIL(ICAS)).AND.(ABS(R(2)).LE.EPSIL(ICASB))) THEN
  1640. GOTO 4020
  1641. ENDIF
  1642. ENDDO
  1643. C
  1644. C Pas de convergence
  1645. C
  1646. SEGSUP,QUASIN
  1647. GOTO 6000
  1648. C
  1649. 4020 CONTINUE
  1650. SEGSUP,QUASIN
  1651. C
  1652. IF ((GAMBDA(ICAS).LT.EPSILO2).OR.
  1653. &(GAMBDA(ICASB).LT.EPSILO2)) THEN
  1654. 4030 CONTINUE
  1655. IF (NPASS.LE.5) THEN
  1656. NPASS = NPASS + 1
  1657. ELSE
  1658. GOTO 6000
  1659. ENDIF
  1660. IF ((GAMBDA(ICAS).LT.EPSILO2).AND.
  1661. &(GAMBDA(ICASB).LT.EPSILO2)) THEN
  1662. C Ce cas n'arrive que si les normales ont ete mal estimes:
  1663. C (increment elastique trop grand) --> on bissecte
  1664. GOTO 6000
  1665. ELSE
  1666. IF (GAMBDA(ICAS).LT.EPSILO2) THEN
  1667. ICAS = ICASB
  1668. IHARD = IHARDB
  1669. ISGNTC = ISGNTCB
  1670. ELSE
  1671. CONTINUE
  1672. ENDIF
  1673. ENDIF
  1674. SIGHP(1) = SIGHP0(1)
  1675. SIGHP(2) = SIGHP0(2)
  1676. SIGT(1) = SIGI(1) + DSIGI(1)
  1677. SIGT(2) = SIGI(2) + DSIGI(2)
  1678. SIGT(4) = SIGI(4) + DSIGI(4)
  1679. GOTO 1000
  1680. ENDIF
  1681. C On teste les mecanismes de la grande surface
  1682. C Si un des criteres est viole (et donc s'il n'y a pas contact)
  1683. C il faut trouver le point de contact
  1684. C
  1685. CALL CONPRI ( SIGT ,CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1686. C
  1687. NCON = 0
  1688. FINT1 = 0.
  1689. FINT2 = 0.
  1690. C
  1691. XFFFS = XZER
  1692. XFFFSB = XZER
  1693. IF (ICAS.LE.4) THEN
  1694. IF (ICAS.LE.2) THEN
  1695. FG(1) = XKISO(1)*RESIST(5)
  1696. ELSE
  1697. FG(1) = XKISO(2)*RESIST(6)
  1698. ENDIF
  1699. CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,XFFFS)
  1700. IF (((XFFFS.GE.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  1701. & ((XFFFS.GE.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  1702. NCON = NCON +1
  1703. ELSE
  1704. IF (((ICAS.LE.2).AND.(ABS(XFFFS).LT.EPSIL(5))).OR.
  1705. &((ICAS.GE.3).AND.(ABS(XFFFS).LT.EPSIL(6)))) THEN
  1706. KONTAC(ICAS) = 1
  1707. ELSE
  1708. CONTINUE
  1709. ENDIF
  1710. ENDIF
  1711. ENDIF
  1712. C
  1713. IF (ICASB.LE.4) THEN
  1714. IF (ICASB.LE.2) THEN
  1715. FG(1) = XKISO(1)*RESIST(5)
  1716. ELSE
  1717. FG(1) = XKISO(2)*RESIST(6)
  1718. ENDIF
  1719. CALL CRITPE (SIGTP(IHARDB),XZER,FG(1),ISGNTCB,XFFFSB)
  1720. IF (((XFFFSB.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  1721. & ((XFFFSB.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  1722. NCON = NCON +1
  1723. ELSE
  1724. IF (((ICASB.LE.2).AND.(ABS(XFFFSB).LT.EPSIL(5))).OR.
  1725. &((ICASB.GT.3).AND.(ABS(XFFFSB).LT.EPSIL(6)))) THEN
  1726. KONTAC(ICASB) = 1
  1727. ELSE
  1728. CONTINUE
  1729. ENDIF
  1730. ENDIF
  1731. ENDIF
  1732. C
  1733. IF (NCON.EQ.1) THEN
  1734. IF (((XFFFS.GE.EPSIL(5)).AND.(ICAS.LE.2)).OR.
  1735. & ((XFFFS.GE.EPSIL(6)).AND.(ICAS.GE.3))) THEN
  1736. GOTO 4040
  1737. ELSE
  1738. ICASI = ICAS
  1739. ISGNTCI = ISGNTC
  1740. IHARDI = IHARD
  1741. ICAS = ICASB
  1742. ISGNTC = ISGNTCB
  1743. IHARD = IHARDB
  1744. ICASB = ICASI
  1745. ISGNTCB = ISGNTCI
  1746. IHARDB = IHARDI
  1747. GOTO 4040
  1748. ENDIF
  1749. ELSE
  1750. IF (NCON.EQ.2) THEN
  1751. GOTO 4040
  1752. ELSE
  1753. CONTINUE
  1754. ENDIF
  1755. CONTINUE
  1756. ENDIF
  1757. C
  1758. C Regles de suivi avant de sortir normalement
  1759. C Direction de l'ecoulement
  1760. IF (ICAS.GE.5) THEN
  1761. SIGHP(IHARD) = XKISO(ICAS-4)*RESIST(ICAS)
  1762. & -RESIST(2*ICAS-10+IHARD)
  1763. ENDIF
  1764. IF (ICASB.GE.5) THEN
  1765. SIGHP(IHARDB) = XKISO(ICASB-4)*RESIST(ICASB)
  1766. & -RESIST(2*ICASB-10+IHARDB)
  1767. ENDIF
  1768. C
  1769. C Mis a jour des indicateurs de contact
  1770. C
  1771. IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).GT.
  1772. &(RESIST(1)-RESIST(3))) THEN
  1773. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) KONTAC(ICAS+2)=0
  1774. IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) KONTAC(ICAS-2)=0
  1775. IF ((ICASB.EQ.1).OR.(ICASB.EQ.2)) KONTAC(ICASB+2)=0
  1776. IF ((ICASB.EQ.3).OR.(ICASB.EQ.4)) KONTAC(ICASB-2)=0
  1777. ELSE
  1778. KONTAC(1) = 1
  1779. KONTAC(2) = 1
  1780. KONTAC(3) = 1
  1781. KONTAC(4) = 1
  1782. ENDIF
  1783. C
  1784. GOTO 9999
  1785. C
  1786. C On va au contact
  1787. 4040 CONTINUE
  1788. SIGHP(1) = SIGHP0(1)
  1789. SIGHP(2) = SIGHP0(2)
  1790. SIGT(1) = SIGI(1) + DSIGI(1)
  1791. SIGT(2) = SIGI(2) + DSIGI(2)
  1792. SIGT(4) = SIGI(4) + DSIGI(4)
  1793. EPSI(1) = EPSI0(1)
  1794. EPSI(2) = EPSI0(2)
  1795. GOTO 5000
  1796. C================================================================
  1797. C CAS 5
  1798. C Deux mecanismes avec
  1799. C Entree en contact de 2 mecanismes
  1800. C================================================================
  1801. 5000 CONTINUE
  1802. C
  1803. CALL CONPRI ( SIGT , CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  1804. C
  1805. C ISGN2 et ISGN2B
  1806. IF (IHARD.EQ.1) THEN
  1807. IF (SIGTP(1).GE.SIGTP(2)) THEN
  1808. ISGN2 = 1
  1809. ELSE
  1810. ISGN2 = -1
  1811. ENDIF
  1812. ELSE
  1813. IF (IHARD.EQ.2) THEN
  1814. IF (SIGTP(1).GE.SIGTP(2)) THEN
  1815. ISGN2 = -1
  1816. ELSE
  1817. ISGN2 = 1
  1818. ENDIF
  1819. ENDIF
  1820. ENDIF
  1821. C
  1822. IF (IHARDB.EQ.1) THEN
  1823. IF (SIGTP(1).GE.SIGTP(2)) THEN
  1824. ISGN2B = 1
  1825. ELSE
  1826. ISGN2B = -1
  1827. ENDIF
  1828. ELSE
  1829. IF (IHARDB.EQ.2) THEN
  1830. IF (SIGTP(1).GE.SIGTP(2)) THEN
  1831. ISGN2B = -1
  1832. ELSE
  1833. ISGN2B = 1
  1834. ENDIF
  1835. ENDIF
  1836. ENDIF
  1837. C
  1838. 5010 CONTINUE
  1839. C
  1840. C Definition du segment
  1841. C
  1842. ITM = 30
  1843. MN = 3
  1844. SEGINI,QUASIN
  1845. IT = -1
  1846. CALL ZERO(D,1,MN)
  1847. C XKISO et XHH2
  1848. IF (ICASB.LE.4) THEN
  1849. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  1850. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  1851. &XKISO(1),XHH2,KERRE)
  1852. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  1853. &XKISO(2),DUMM,KERRE)
  1854. ELSE
  1855. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,
  1856. &XKISO(1),DUMM,KERRE)
  1857. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,
  1858. &XKISO(2),XHH2,KERRE)
  1859. ENDIF
  1860. ELSE
  1861. IF (ICASB.EQ.5) THEN
  1862. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  1863. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1864. & ,XHH2,KERRE)
  1865. XHH2B=XHH2
  1866. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1867. & ,DUMM,KERRE)
  1868. ELSE
  1869. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1870. & ,XHH2B,KERRE)
  1871. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1872. & ,XHH2,KERRE)
  1873. ENDIF
  1874. ELSE
  1875. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  1876. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1877. & ,XHH2,KERRE)
  1878. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1879. & ,XHH2B,KERRE)
  1880. ELSE
  1881. CALL VALEUR (EPSI0(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  1882. & ,DUMM,KERRE)
  1883. CALL VALEUR (EPSI0(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  1884. & ,XHH2,KERRE)
  1885. XHH2B = XHH2
  1886. ENDIF
  1887. ENDIF
  1888. ENDIF
  1889. C
  1890. XKISO0(1)=XKISO(1)
  1891. EPSI(1)=EPSI0(1)
  1892. XKISO0(2)=XKISO(2)
  1893. EPSI(2)=EPSI0(2)
  1894. XHH2I = XHH2
  1895. XHH2IB = XHH2B
  1896. XKISOI(1)=XKISO(1)
  1897. XKISOI(2)=XKISO(2)
  1898. EPSII(1)=EPSI0(1)
  1899. EPSII(2)=EPSI0(2)
  1900. C
  1901. DO 5011,IE1=1,6
  1902. FFF(IE1) = XZER
  1903. 5011 CONTINUE
  1904. C
  1905. CALL CRITPE (SIGTP(IHARD),SIGHP0(IHARD),RESIST(ICAS),
  1906. & ISGNTC,FFF(ICAS))
  1907. C
  1908. IF (ICASB.LE.4) THEN
  1909. CALL CRITPE (SIGTP(IHARDB),SIGHP0(IHARDB),RESIST(ICASB),
  1910. & ISGNTCB,FFF(ICASB))
  1911. ELSE
  1912. FG(1) = XKISOI(ICASB -4)*RESIST(ICASB)
  1913. CALL CRITPE (SIGTP(IHARDB), XZER, FG(1),ISGNTCB,FFF(ICASB))
  1914. ENDIF
  1915. C
  1916. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  1917. FG(1) = XKISO(1)*RESIST(5)
  1918. CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC,XFFFS)
  1919. ELSE
  1920. FG(1) = XKISO(2)*RESIST(6)
  1921. CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,XFFFS)
  1922. ENDIF
  1923. C
  1924. CALL NORMBL(SIGT,RESIST(ICAS),ISGNTC,ISGN2,1,XNOR1)
  1925. CALL NORMBL(SIGT,RESIST(ICASB),ISGNTCB,ISGN2B,1,XNOR2)
  1926. C
  1927. C Jacobien
  1928. C
  1929. XJ00(1,1) = -XXCCYY(XNOR1, XNOR1, YOUN, XNU) - XHH1
  1930. XJ00(1,2) = -XXCCYY(XNOR1, XNOR2, YOUN, XNU)
  1931. XJ00(1,3) = -XXYY(XNOR1 , DSIGI)
  1932. XJ00(2,1) = -XXCCYY(XNOR2, XNOR1, YOUN, XNU)
  1933. XJ00(2,2) = -XXCCYY(XNOR2, XNOR2, YOUN, XNU)
  1934. XJ00(2,3) = -XXYY(XNOR2,DSIGI)
  1935. XJ00(3,1) = -XXCCYY(XNOR1,XNOR1,YOUN,XNU)
  1936. XJ00(3,2) = -XXCCYY(XNOR1,XNOR2,YOUN,XNU)
  1937. XJ00(3,3) = -XXYY(XNOR1,DSIGI)
  1938. C
  1939. IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) THEN
  1940. XJ00(3,1)=XJ00(3,1)+RESIST(6)*XHH2*DEGRAD*ISGNTC*RESIST(ICAS)
  1941. ENDIF
  1942. C
  1943. IF (ICASB.LE.4) THEN
  1944. XJ00(2,2) = XJ00(2,2) - XHH1
  1945. IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) THEN
  1946. XJ00(3,2)=XJ00(3,2)+RESIST(6)*XHH2
  1947. &*DEGRAD*ISGNTCB*RESIST(ICASB)
  1948. ENDIF
  1949. ELSE
  1950. XJ00(2,2) = XJ00(2,2) - RESIST(ICASB)*XHH2B*ISGNTCB
  1951. ENDIF
  1952. C
  1953. C
  1954. R(1) = FFF(ICAS)
  1955. R(2) = FFF(ICASB)
  1956. R(3) = XFFFS
  1957. C
  1958. SPLITT = XZER
  1959. GAMBDA(1) = XZER
  1960. GAMBDA(2) = XZER
  1961. GAMBDA(3) = XZER
  1962. GAMBDA(4) = XZER
  1963. GAMBDA(5) = XZER
  1964. GAMBDA(6) = XZER
  1965. C
  1966. C Inverse du jacobien
  1967. C
  1968. XH0(1,1) = XJ00(1,1)
  1969. XH0(1,2) = XJ00(1,2)
  1970. XH0(1,3) = XJ00(1,3)
  1971. XH0(2,1) = XJ00(2,1)
  1972. XH0(2,2) = XJ00(2,2)
  1973. XH0(2,3) = XJ00(2,3)
  1974. XH0(3,1) = XJ00(3,1)
  1975. XH0(3,2) = XJ00(3,2)
  1976. XH0(3,3) = XJ00(3,3)
  1977. CALL INVALM (XH0, MN , MN , IRD, 1.D-12)
  1978. IF (IRD.NE.0) THEN
  1979. SEGSUP,QUASIN
  1980. GOTO 6000
  1981. ENDIF
  1982. C
  1983. C On rentre dans les iterations internes
  1984. C
  1985. DO I=0,ITM
  1986. C
  1987. C Appel de Broyden
  1988. C
  1989. CALL BROYDE (QUASIN)
  1990. C
  1991. GAMBDA(ICAS) = D(1)
  1992. GAMBDA(ICASB) = D(2)
  1993. SPLITT = D(3)
  1994. C
  1995. C Calcul du phi et des termes de la matrice A et A-1
  1996. C
  1997. SITPRO(1) = SIGI(1) + (UN-SPLITT)*DSIGI(1)
  1998. & - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  1999. & + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  2000. SITPRO(2) = SIGI(2) + (UN-SPLITT)*DSIGI(2)
  2001. & - UNDEMI*(GAMBDA(1) + GAMBDA(2) + GAMBDA(5))*CPI
  2002. & + UNDEMI*(GAMBDA(3) + GAMBDA(4) + GAMBDA(6))*CPI
  2003. SITPRO(4) = SIGI(4) + (UN-SPLITT)*DSIGI(4)
  2004. C
  2005. C PHI et AINV
  2006. C
  2007. SIGHP(1) = SIGHP0(1) + XHH1*(GAMBDA(1)-GAMBDA(3))
  2008. SIGHP(2) = SIGHP0(2) + XHH1*(GAMBDA(2)-GAMBDA(4))
  2009. EPSI(1) = EPSI0(1) + GAMBDA(5)
  2010. WORK = ISGNTC*RESIST(ICAS)*GAMBDA(ICAS)
  2011. IF (ICASB.LE.4) THEN
  2012. WORK = WORK + ISGNTCB*RESIST(ICASB)*GAMBDA(ICASB)
  2013. ENDIF
  2014. EPSI(2)=EPSI0(2)+ GAMBDA(6) + DEGRAD*WORK
  2015. C
  2016. IF (ICASB.LE.4) THEN
  2017. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2018. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  2019. &XKISO(1),XHH2,KERRE)
  2020. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  2021. &XKISO(2),DUMM,KERRE)
  2022. ELSE
  2023. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,
  2024. &XKISO(1),DUMM,KERRE)
  2025. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,
  2026. &XKISO(2),XHH2,KERRE)
  2027. ENDIF
  2028. ELSE
  2029. IF (ICASB.EQ.5) THEN
  2030. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2031. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  2032. & ,XHH2,KERRE)
  2033. XHH2B=XHH2
  2034. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  2035. & ,DUMM,KERRE)
  2036. ELSE
  2037. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  2038. & ,XHH2B,KERRE)
  2039. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  2040. & ,XHH2,KERRE)
  2041. ENDIF
  2042. ELSE
  2043. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2044. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  2045. & ,XHH2,KERRE)
  2046. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  2047. & ,XHH2B,KERRE)
  2048. ELSE
  2049. CALL VALEUR (EPSI(1),TRAC(IPT),NCURVT,XLCAT,XKISO(1)
  2050. & ,DUMM,KERRE)
  2051. CALL VALEUR (EPSI(2),TRAC(IPC),NCURVC,XLCAC,XKISO(2)
  2052. & ,XHH2,KERRE)
  2053. XHH2B = XHH2
  2054. ENDIF
  2055. ENDIF
  2056. ENDIF
  2057. C
  2058. C Calcul des phi
  2059. C
  2060. PHI1=ISGNTC*(-UNDEMI*(SITPRO(1) + SITPRO(2))
  2061. & + SIGHP(IHARD) + RESIST(ICAS))
  2062. C
  2063. IF (ICASB.LE.4) THEN
  2064. PHI2=ISGNTCB*(-UNDEMI*(SITPRO(1) + SITPRO(2))
  2065. & + SIGHP(IHARDB)+ RESIST(ICASB))
  2066. ELSE
  2067. PHI2=ISGNTCB*(-UNDEMI*(SITPRO(1)+SITPRO(2))
  2068. & + (XKISOI(ICASB-4)+XHH2IB*(EPSI(ICASB-4)-EPSII(ICASB-4)))
  2069. & *RESIST(ICASB))
  2070. ENDIF
  2071. C
  2072. C Test sur phi1 et phi2
  2073. C
  2074. IF ((ABS(PHI2)).GT.(ABS(EPSILO*RESIST(ICASB)))) THEN
  2075. SOMLAM = GAMBDA(ICAS) + GAMBDA(ICASB)*PHI1/PHI2
  2076. PHI = PHI1
  2077. ELSE
  2078. IF ((ABS(PHI1)).GT.(ABS(EPSILO*RESIST(ICAS)))) THEN
  2079. SOMLAM = GAMBDA(ICASB) + GAMBDA(ICAS)*PHI2/PHI1
  2080. PHI = PHI2
  2081. ELSE
  2082. SOMLAM = GAMBDA(ICASB) + GAMBDA(ICAS)
  2083. PHI = PHI1
  2084. ENDIF
  2085. ENDIF
  2086. C
  2087. AINV(1) = (DEUX*PHI+GGG*SOMLAM)/(DEUX*PHI+DEUX*GGG*SOMLAM)
  2088. AINV(2) = (GGG*SOMLAM)/(DEUX*PHI + DEUX*GGG*SOMLAM)
  2089. AINV(3) = PHI/(PHI + GGG*SOMLAM)
  2090. C
  2091. C Calcul des contraintes a l'aide de AINV
  2092. C
  2093. SIGT(1) = AINV(1)*SITPRO(1)+AINV(2)*SITPRO(2)
  2094. SIGT(2) = AINV(2)*SITPRO(1)+AINV(1)*SITPRO(2)
  2095. SIGT(4) = AINV(3)*SITPRO(4)
  2096. C
  2097. C Nouvelles contraintes principales
  2098. C
  2099. CALL CONPRI ( SIGT , CPHI0,SPHI0,SIGTP,CPHI,SPHI)
  2100. C
  2101. CALL CRITPE (SIGTP(IHARD),SIGHP(IHARD),RESIST(ICAS),
  2102. & ISGNTC,FFF(ICAS))
  2103. R(1)= + FFF(ICAS)
  2104. C
  2105. IF (ICASB.LE.4) THEN
  2106. CALL CRITPE (SIGTP(IHARDB),SIGHP(IHARDB),RESIST(ICASB),
  2107. & ISGNTCB,FFF(ICASB))
  2108. R(2)= + FFF(ICASB)
  2109. ELSE
  2110. FG(1) = (XKISOI(ICASB -4)+XHH2IB*(EPSI(ICASB-4)
  2111. & - EPSII(ICASB-4)))*RESIST(ICASB)
  2112. CALL CRITPE (SIGTP(IHARDB), XZER, FG(1),ISGNTCB,FFF(ICASB))
  2113. R(2) = + FFF(ICASB)
  2114. & - ISGNTCB*RESIST(ICASB)
  2115. & * (XKISO(ICASB-4) - XKISOI(ICASB-4)
  2116. & -(EPSI(ICASB-4) - EPSII(ICASB-4))*XHH2IB)
  2117. IF (ABS(EPSI(ICASB-4) - EPSII(ICASB-4)).LT.EPSILO2) THEN
  2118. XHH2IB=XHH2B
  2119. ELSE
  2120. XHH2IB = (XKISO(ICASB-4) - XKISOI(ICASB-4))/
  2121. & (EPSI(ICASB-4) - EPSII(ICASB-4))
  2122. ENDIF
  2123. ENDIF
  2124. C
  2125. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2126. FG(1) = XKISO(1)*RESIST(5)
  2127. CALL CRITPE (SIGTP(IHARD), XZER, FG(1), ISGNTC,XFFFS)
  2128. R(3) = XFFFS
  2129. ELSE
  2130. FG(1) = (XKISOI(2)+XHH2I*(EPSI(2)-EPSII(2)))*RESIST(6)
  2131. CALL CRITPE (SIGTP(IHARD),XZER,FG(1),ISGNTC,XFFFS)
  2132. R(3) = XFFFS + RESIST(6)*
  2133. & (XKISO(2) - XKISOI(2)-(EPSI(2) - EPSII(2))*XHH2I)
  2134. IF (ABS(EPSI(2) - EPSII(2)).LT.EPSILO2) THEN
  2135. XHH2I=XHH2
  2136. ELSE
  2137. XHH2I = (XKISO(2) - XKISOI(2))/
  2138. & (EPSI(2) - EPSII(2))
  2139. ENDIF
  2140. ENDIF
  2141. XKISOI(1) = XKISO(1)
  2142. EPSII(1) = EPSI(1)
  2143. XKISOI(2) = XKISO(2)
  2144. EPSII(2) = EPSI(2)
  2145. C
  2146. C Test de convergence
  2147. C
  2148. IF (ICAS.LE.2) THEN
  2149. IF ((ABS(FFF(ICAS)).LT.EPSIL(ICAS)).AND.
  2150. & (ABS(R(1)).LT.EPSIL(ICAS)).AND.
  2151. & (ABS(FFF(ICASB)).LT.EPSIL(ICASB)).AND.
  2152. & (ABS(R(2)).LT.EPSIL(ICASB)).AND.
  2153. & (ABS(XFFFS).LT.EPSIL(5)).AND.
  2154. & (ABS(R(3)).LT.EPSIL(5))) THEN
  2155. GOTO 5020
  2156. ENDIF
  2157. ELSE
  2158. IF ((ABS(FFF(ICAS)).LT.EPSIL(ICAS)).AND.
  2159. & (ABS(R(1)).LT.EPSIL(ICAS)).AND.
  2160. & (ABS(FFF(ICASB)).LT.EPSIL(ICASB)).AND.
  2161. & (ABS(R(2)).LT.EPSIL(ICASB)).AND.
  2162. & (ABS(XFFFS).LT.EPSIL(6)).AND.
  2163. & (ABS(R(3)).LT.EPSIL(6))) THEN
  2164. GOTO 5020
  2165. ENDIF
  2166. ENDIF
  2167. ENDDO
  2168. C
  2169. C Pas de convergence
  2170. C
  2171. SEGSUP,QUASIN
  2172. GOTO 6000
  2173. C
  2174. 5020 CONTINUE
  2175. SEGSUP,QUASIN
  2176. C
  2177. C
  2178. IF (((SPLITT-1.).LE.EPSILO).AND.(SPLITT.GE.(-1.*EPSILO)).AND.
  2179. &(GAMBDA(ICAS).GE.EPSILO2)
  2180. &.AND.(GAMBDA(ICASB).GE.EPSILO2)) THEN
  2181. C
  2182. C On a les contraintes, les variables d'ecrouissage cinemat. et tgphi
  2183. C Mis a jour des variables d'ecrouissage isotrope ou regles de suivi
  2184. C
  2185.  
  2186. IF (NCON.EQ.2) THEN
  2187. IF (ICASB.LE.4) THEN
  2188. C
  2189. IF (ICASB.LE.2) THEN
  2190. FG(1) = XKISO(1)*RESIST(5)
  2191. ELSE
  2192. FG(1) = XKISO(2)*RESIST(6)
  2193. ENDIF
  2194. CALL CRITPE (SIGTP(IHARDB),XZER,FG(1),ISGNTCB,XFFFSB)
  2195. C
  2196. C On regarde si on ne s'est pas trompe d'ordre
  2197. C
  2198. IF (((XFFFSB.GE.EPSIL(5)).AND.(ICASB.LE.2)).OR.
  2199. & ((XFFFSB.GE.EPSIL(6)).AND.(ICASB.GE.3))) THEN
  2200. ICASI = ICAS
  2201. ISGNTCI = ISGNTC
  2202. IHARDI = IHARD
  2203. ICAS = ICASB
  2204. ISGNTC = ISGNTCB
  2205. IHARD = IHARDB
  2206. ICASB = ICASI
  2207. ISGNTCB = ISGNTCI
  2208. IHARDB = IHARDI
  2209. GOTO 4040
  2210. ELSE
  2211. IF ((ICASB.LE.2).AND.(ABS(XFFFSB).LT.EPSIL(5))) THEN
  2212. KONTAC(ICASB) = 1
  2213. ICASB=5
  2214. ENDIF
  2215. IF (((ICASB.EQ.3).OR.(ICASB.EQ.4))
  2216. & .AND.(ABS(XFFFSB).LT.EPSIL(6))) THEN
  2217. KONTAC(ICASB) = 1
  2218. ICASB=6
  2219. ENDIF
  2220. ENDIF
  2221. C
  2222. ELSE
  2223. CONTINUE
  2224. ENDIF
  2225. ELSE
  2226. CONTINUE
  2227. ENDIF
  2228. C
  2229. C Le contact est realise
  2230. C
  2231. KONTAC(ICAS) = 1
  2232. IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).GT.
  2233. &(RESIST(1)-RESIST(3))) THEN
  2234. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) KONTAC(ICAS+2)=0
  2235. IF ((ICAS.EQ.3).OR.(ICAS.EQ.4)) KONTAC(ICAS-2)=0
  2236. IF ((ICASB.EQ.1).OR.(ICASB.EQ.2)) KONTAC(ICASB+2)=0
  2237. IF ((ICASB.EQ.3).OR.(ICASB.EQ.4)) KONTAC(ICASB-2)=0
  2238. ELSE
  2239. KONTAC(1) = 1
  2240. KONTAC(2) = 1
  2241. KONTAC(3) = 1
  2242. KONTAC(4) = 1
  2243. ENDIF
  2244. C
  2245. SIGI(1) = SIGT(1)
  2246. SIGI(2) = SIGT(2)
  2247. SIGI(4) = SIGT(4)
  2248. DSIGI(1) = DSIGI(1)*SPLITT
  2249. DSIGI(2) = DSIGI(2)*SPLITT
  2250. DSIGI(4) = DSIGI(4)*SPLITT
  2251. SIGT(1) = SIGI(1) + DSIGI(1)
  2252. SIGT(2) = SIGI(2) + DSIGI(2)
  2253. SIGT(4) = SIGI(4) + DSIGI(4)
  2254. SIGHP0(1) = SIGHP(1)
  2255. SIGHP0(2) = SIGHP(2)
  2256. EPSI0(1) = EPSI(1)
  2257. EPSI0(2) = EPSI(2)
  2258. NPASS2 = 0
  2259. CPHI0 = CPHI
  2260. SPHI0 = SPHI
  2261. IF ((ICAS.EQ.1).OR.(ICAS.EQ.2)) THEN
  2262. ICAS = 5
  2263. ELSE
  2264. ICAS = 6
  2265. ENDIF
  2266. IF (ICASB.LE.4) THEN
  2267. GOTO 4000
  2268. ELSE
  2269. IF (ICASB.EQ.ICAS) THEN
  2270. C On passe de 2 mecanismes a 1 mecanisme
  2271. GOTO 1000
  2272. ELSE
  2273. GOTO 4000
  2274. ENDIF
  2275. ENDIF
  2276. C
  2277. C
  2278. ELSE
  2279. 5040 CONTINUE
  2280. IF ((GAMBDA(ICASB).LT.EPSILO2).AND.(NPASS.LE.3)) THEN
  2281. SIGHP(1) = SIGHP0(1)
  2282. SIGHP(2) = SIGHP0(2)
  2283. SIGT(1) = SIGI(1) + DSIGI(1)
  2284. SIGT(2) = SIGI(2) + DSIGI(2)
  2285. SIGT(4) = SIGI(4) + DSIGI(4)
  2286. IF (NPASS.LT.5) THEN
  2287. NPASS = NPASS + 1
  2288. GOTO 3000
  2289. ELSE
  2290. GOTO 6000
  2291. ENDIF
  2292. ELSE
  2293. IF ((GAMBDA(ICAS).LT.EPSILO2).AND.(NPASS2.EQ.0)) THEN
  2294. NPASS2 = 1
  2295. SIGHP(1) = SIGHP0(1)
  2296. SIGHP(2) = SIGHP0(2)
  2297. SIGT(1) = SIGI(1) + DSIGI(1)
  2298. SIGT(2) = SIGI(2) + DSIGI(2)
  2299. SIGT(4) = SIGI(4) + DSIGI(4)
  2300. ICAS = ICASB
  2301. IHARD = IHARDB
  2302. ISGNTC = ISGNTCB
  2303. IF (NPASS.LT.5) THEN
  2304. NPASS = NPASS + 1
  2305. GOTO 1000
  2306. ELSE
  2307. GOTO 6000
  2308. ENDIF
  2309. ELSE
  2310. GOTO 6000
  2311. ENDIF
  2312. ENDIF
  2313. ENDIF
  2314. C
  2315. C On bissecte
  2316. C
  2317. 6000 CONTINUE
  2318. NDICHO = NDICHO + 1
  2319. IF (NDICHO.LE.90) THEN
  2320. IDICHO = 1
  2321. DSIGIDI(1) = DSIGIDI(1) + UNDEMI*DSIGI(1)
  2322. DSIGIDI(2) = DSIGIDI(2) + UNDEMI*DSIGI(2)
  2323. DSIGIDI(4) = DSIGIDI(4) + UNDEMI*DSIGI(4)
  2324. DSIGI(1) = DSIGI(1)*UNDEMI
  2325. DSIGI(2) = DSIGI(2)*UNDEMI
  2326. DSIGI(4) = DSIGI(4)*UNDEMI
  2327. CPHI0 = CPHI00
  2328. SPHI0 = SPHI00
  2329. GOTO 10
  2330. ELSE
  2331. KERRE = 2
  2332. RETURN
  2333. ENDIF
  2334. C
  2335. C=================================================================
  2336. C MIS A JOUR DE TOUTES LES VARIABLES ET CONTRAINTES AVANT
  2337. C DE SORTIR
  2338. C=================================================================
  2339. 9999 CONTINUE
  2340. C Cas de la bissection
  2341. IF (IDICHO.EQ.1) THEN
  2342. IDICHO = 0
  2343. NDICHO = 0
  2344. NPASS = 0
  2345. DSIGI(1) = DSIGIDI(1)
  2346. DSIGI(2) = DSIGIDI(2)
  2347. DSIGI(4) = DSIGIDI(4)
  2348. SIGI(1) = SIGT(1)
  2349. SIGI(2) = SIGT(2)
  2350. SIGI(4) = SIGT(4)
  2351. DSIGIDI(1) = XZER
  2352. DSIGIDI(2) = XZER
  2353. DSIGIDI(4) = XZER
  2354. SIGHP0(1) = SIGHP(1)
  2355. SIGHP0(2) = SIGHP(2)
  2356. EPSI0(1) = EPSI(1)
  2357. EPSI0(2) = EPSI(2)
  2358. CPHI0 = CPHI
  2359. CPHI00 = CPHI0
  2360. SPHI0 = SPHI
  2361. SPHI00 = SPHI0
  2362. KONTA0(1) = KONTAC(1)
  2363. KONTA0(2) = KONTAC(2)
  2364. KONTA0(3) = KONTAC(3)
  2365. KONTA0(4) = KONTAC(4)
  2366. GOTO 10
  2367. ENDIF
  2368. C STRESS
  2369. SIGF(1) = SIGT(1)
  2370. SIGF(2) = SIGT(2)
  2371. SIGF(4) = SIGT(4)
  2372. SIGF(3) = XZER
  2373. C PLASTIC STRAIN
  2374. CINV(1) = UN/YOUN
  2375. CINV(2) = -XNU/YOUN
  2376. CINV(3) = UN/GGG
  2377. DEFP(1) =CINV(1)*(SIG0(1) + DSIGT(1) - SIGT(1))
  2378. & +CINV(2)*(SIG0(2) + DSIGT(2) - SIGT(2))
  2379. DEFP(2) =CINV(2)*(SIG0(1) + DSIGT(1) - SIGT(1))
  2380. & +CINV(1)*(SIG0(2) + DSIGT(2) - SIGT(2))
  2381. DEFP(4) =CINV(3)*(SIG0(4) + DSIGT(4) - SIGT(4))
  2382. DEFP(3) = XZER
  2383. C
  2384. C Cas ou la grande surface est a l'interieur de la petite
  2385. C
  2386. IF ((XKISO(1)*RESIST(5)-XKISO(2)*RESIST(6)).LE.
  2387. &(RESIST(1)-RESIST(3))) THEN
  2388. KONTAC(1) = 1
  2389. KONTAC(2) = 1
  2390. KONTAC(3) = 1
  2391. KONTAC(4) = 1
  2392. ENDIF
  2393. VARF(1) = SIGHP(1)
  2394. VARF(2) = SIGHP(2)
  2395. VARF(3) = EPSI(1)
  2396. VARF(4) = EPSI(2)
  2397. VARF(5) = CPHI - 1.D0
  2398. VARF(6) = SPHI
  2399. VARF(7) = KONTAC(1)
  2400. VARF(8) = KONTAC(2)
  2401. VARF(9) = KONTAC(3)
  2402. VARF(10) = KONTAC(4)
  2403. VARF(11) = ICAS
  2404. VARF(12) = ICASB
  2405. VARF(13) = GAMBDA(ICAS)
  2406. VARF(14) = GAMBDA(ICASB)
  2407. C=====================================================
  2408. C FIN DE LA ROUTINE DU MODELE
  2409. C======================================================
  2410. RETURN
  2411. END
  2412.  
  2413.  
  2414.  
  2415.  

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