Télécharger iwpr3d.eso

Retour à la liste

Numérotation des lignes :

  1. C IWPR3D SOURCE AF221230 16/11/23 21:15:07 9204
  2. C IWPR3D SOURCE AF221230 12/10/04 21:15:18 7520
  3. SUBROUTINE IWPR3D (XMAT,VAR0,VARF,SIG0,SIGF,DEPST,XCAR,
  4. & EPIN0,EPINF,EPST0,EPSTF)
  5. C
  6. C====&===1=========2=========3=========4=========5=========6=========7==
  7. C Commentaires : Loi d Iwan generalisée en 3D et basee sur les travaux
  8. C de Prevost
  9. C Traits : -
  10. C
  11. C Auteur : A. Frau (Ing.) - CEA/DEN/DANS/DM2S/SEMT/EMSI
  12. C====&===1=========2=========3=========4=========5=========6=========7==
  13. C
  14. C----DECLARATION GENERALES----------------------------------------------
  15. C
  16. IMPLICIT INTEGER(I-N)
  17. IMPLICIT REAL*8(A-H,O-Z)
  18. -INC CCOPTIO
  19. C
  20. C----DECLARATION VARIABLES----------------------------------------------
  21. C
  22. REAL*8 XID1(6),XID2(6,6)
  23. C
  24. C
  25. REAL*8 XMAT(*),SIG0(*),SIGF(*),VAR0(*),VARF(*),DEPST(*),XCAR(*)
  26. REAL*8 EPIN0(*),EPINF(*),EPST0(*),EPSTF(*)
  27. C
  28. C Parametres materiaux
  29. REAL*8 XE1,XNU1,XE0,XA0,XC1,XN1,XPREF1,XK1
  30. REAL*8 XM1, XM2, XM3, XM4, XM5, XM6, XM7, XM8, XM9, XM10
  31. REAL*8 XHC1, XHC2, XHC3, XHC4, XHC5, XHC6, XHC7, XHC8, XHC9,XHC10
  32. REAL*8 XHE1, XHE2, XHE3, XHE4, XHE5, XHE6, XHE7, XHE8, XHE9,XHE10
  33. REAL*8 XD1, XD2, XD3, XD4, XD5, XD6, XD7, XD8, XD9, XD10
  34. REAL*8 XETC1,XETC2,XETC3,XETC4,XETC5,XETC6,XETC7,XETC8,XETC9
  35. REAL*8 XETC10
  36. REAL*8 XETE1,XETE2,XETE3,XETE4,XETE5,XETE6,XETE7,XETE8,XETE9
  37. REAL*8 XETE10
  38. C
  39. C Valeurs elastiques
  40. REAL*8 GMOD_E1,KMOD_E1,GMOD,KMOD
  41. C
  42. C tollerances
  43. REAL*8 XTOL1,XTOL2
  44. C
  45. C indicateurs de controle pour le substep
  46. INTEGER IND_TEST,IND_STEP1,IND_STEP2,IND_STEP3,N_STEP,N_STEP2
  47. C
  48. C variables deplacement
  49. REAL*8 DDET_ST(6),DDET_ST2(6),DE_TOT(6),DSE_TOT(6),DVE_TOT
  50. c
  51. c Valeurs pour contraintes, deformations et autres valeurs initiaux
  52. REAL*8 SIG_INI(6),S_INI(6),P_INI
  53. REAL*8 X1_INI(6),X2_INI(6),X3_INI(6),X4_INI(6),X5_INI(6)
  54. REAL*8 X6_INI(6),X7_INI(6),X8_INI(6),X9_INI(6),X10_INI(6)
  55. REAL*8 NX1_INI,NX2_INI,NX3_INI,NX4_INI,NX5_INI
  56. REAL*8 NX6_INI,NX7_INI,NX8_INI,NX9_INI,NX10_INI
  57. REAL*8 NN1_INI(6),NN2_INI(6),NN3_INI(6),NN4_INI(6),NN5_INI(6)
  58. REAL*8 NN6_INI(6),NN7_INI(6),NN8_INI(6),NN9_INI(6),NN10_INI(6)
  59. REAL*8 NN21_INI(6),NN22_INI(6),NN23_INI(6),NN24_INI(6),NN25_INI(6)
  60. REAL*8 NN26_INI(6),NN27_INI(6),NN28_INI(6),NN29_INI(6)
  61. REAL*8 NN210_INI(6)
  62. REAL*8 COS3TH1_INI,COS3TH2_INI,COS3TH3_INI,COS3TH4_INI
  63. REAL*8 COS3TH5_INI
  64. REAL*8 COS3TH6_INI,COS3TH7_INI,COS3TH8_INI,COS3TH9_INI
  65. REAL*8 COS3TH10_INI
  66. REAL*8 RTH1_INI,RTH2_INI,RTH3_INI,RTH4_INI,RTH5_INI
  67. REAL*8 RTH6_INI,RTH7_INI,RTH8_INI,RTH9_INI,RTH10_INI
  68. REAL*8 F1_INI,F2_INI,F3_INI,F4_INI,F5_INI
  69. REAL*8 F6_INI,F7_INI,F8_INI,F9_INI,F10_INI
  70. REAL*8 BET1_INI(6),BET2_INI(6),BET3_INI(6),BET4_INI(6),BET5_INI(6)
  71. REAL*8 BET6_INI(6),BET7_INI(6),BET8_INI(6),BET9_INI(6)
  72. REAL*8 BET10_INI(6)
  73. REAL*8 EPLA_INI(6),SEPLA_INI(6),VEPLA_INI,ETOT_INI(6)
  74. INTEGER IINI_CRIT
  75. C
  76. C
  77. INTEGER N_ITER1
  78. c
  79. c Valeurs pour contraintes, deformations et autres valeurs pour
  80. c l'instant n+1 au demarrage des boucles
  81. REAL*8 X1_0(6),X2_0(6),X3_0(6),X4_0(6),X5_0(6)
  82. REAL*8 X6_0(6),X7_0(6),X8_0(6),X9_0(6),X10_0(6)
  83. REAL*8 NX1_0,NX2_0,NX3_0,NX4_0,NX5_0
  84. REAL*8 NX6_0,NX7_0,NX8_0,NX9_0,NX10_0
  85. REAL*8 NN1_0(6),NN2_0(6),NN3_0(6),NN4_0(6),NN5_0(6)
  86. REAL*8 NN6_0(6),NN7_0(6),NN8_0(6),NN9_0(6),NN10_0(6)
  87. REAL*8 NN21_0(6),NN22_0(6),NN23_0(6),NN24_0(6),NN25_0(6)
  88. REAL*8 NN26_0(6),NN27_0(6),NN28_0(6),NN29_0(6)
  89. REAL*8 NN210_0(6)
  90. REAL*8 COS3TH1_0,COS3TH2_0,COS3TH3_0,COS3TH4_0
  91. REAL*8 COS3TH5_0
  92. REAL*8 COS3TH6_0,COS3TH7_0,COS3TH8_0,COS3TH9_0
  93. REAL*8 COS3TH10_0
  94. REAL*8 RTH1_0,RTH2_0,RTH3_0,RTH4_0,RTH5_0
  95. REAL*8 RTH6_0,RTH7_0,RTH8_0,RTH9_0,RTH10_0
  96. REAL*8 F1_0,F2_0,F3_0,F4_0,F5_0
  97. REAL*8 F6_0,F7_0,F8_0,F9_0,F10_0
  98. REAL*8 BET1_0(6),BET2_0(6),BET3_0(6),BET4_0(6),BET5_0(6)
  99. REAL*8 BET6_0(6),BET7_0(6),BET8_0(6),BET9_0(6)
  100. REAL*8 BET10_0(6)
  101. REAL*8 VET_TEST(10)
  102. INTEGER VET_TEST2(10)
  103. INTEGER I1_CRIT,I2_CRIT
  104. REAL*8 XCONT_CRI
  105.  
  106. c Valeurs pour contraintes, deformations et autres valeurs pour
  107. c l'instant n+1 au cours des iterations une fois que le critere
  108. C est detecté
  109. REAL*8 XM,XHC,XHE,XETC,XETE,XD
  110. REAL*8 S_0(6),P_0
  111. REAL*8 SEPLA_0(6),VEPLA_0
  112. REAL*8 DGAM_0
  113. REAL*8 BET_0(6),X_0(6),NX_0,NN_0(6),N2N_0(6)
  114. REAL*8 COS3TH_0,RTH_0
  115. REAL*8 HMOD_0
  116. REAL*8 ETA_0,P1_0
  117. REAL*8 S1M(6),S2M(6),S2SCAL,MM_0(6),NMM_0,XVAL1(6),XVAL2
  118. REAL*8 F_0,F_00
  119. C
  120. C criteres
  121. REAL*8 CRIT1
  122. REAL*8 CRIT2
  123. C
  124. C partie elastiques
  125. REAL*8 GMOD_0,KMOD_0,DGDP_0,DKDP_0
  126. REAL*8 R_0(7),DTOT(7)
  127. REAL*8 AE_0I(7,7),AE_0(7,7)
  128. C
  129. C partie iteration
  130. REAL*8 D2F(6,6),DS_0(6),DP_0
  131. REAL*8 DFDS_0(6),DFDS_0D(6),DFDP_0,VECT1,DFDBET_0(6)
  132. REAL*8 XVAL3(6),XVAL4
  133. REAL*8 VAL1,VAL2,VAL3
  134. REAL*8 D2GAM,DSEPLA_0(6),DVEPLA_0,DBET_0(6)
  135. REAL*8 SEPLA_1(6),VEPLA_1,BET_1(6),DGAM_1
  136.  
  137. c Valeurs suite à la convergence à chaque sous-step n+1
  138. REAL*8 EPLA_FIN(6),ETOT_FIN(6)
  139. REAL*8 SIG_FIN(6)
  140. REAL*8 A_CORR,XXI_CORR(6)
  141. REAL*8 BET1_FIN(6),BET2_FIN(6),BET3_FIN(6),BET4_FIN(6),BET5_FIN(6)
  142. REAL*8 BET6_FIN(6),BET7_FIN(6),BET8_FIN(6),BET9_FIN(6)
  143. REAL*8 BET10_FIN(6)
  144. C
  145. C Definition des matrices identité
  146. XID1(1:6) = 0.0D0
  147. XID2(1:6,1:6) = 0.0D0
  148. XID1(1:3) = 1.0D0
  149. XID2(1,1) = 1.0D0
  150. XID2(2,2) = 1.0D0
  151. XID2(3,3) = 1.0D0
  152. XID2(4,4) = 1.0D0
  153. XID2(5,5) = 1.0D0
  154. XID2(6,6) = 1.0D0
  155. C
  156. C----CONTROLE OPTION CALCUL --------------------------------------------
  157. IF ((IFOUR.NE.2).AND.(IFOUR.NE.0).AND.(IFOUR.NE.-1)) THEN
  158. WRITE(IOIMP,10000)
  159. 10000 FORMAT('Option calcul mauvais');
  160. WRITE(IOIMP,10001)
  161. 10001 FORMAT(' seulement:');
  162. WRITE(IOIMP,10002)
  163. 10002 FORMAT(' 3D MODE PLAN');
  164. WRITE(IOIMP,10003)
  165. 10003 FORMAT(' 2D MODE AXIS');
  166. WRITE(IOIMP,10004)
  167. 10004 FORMAT(' 2D MODE PLAN');
  168. STOP
  169. ENDIF
  170. C
  171. C----MISE EN DONNEES----------------------------------------------------
  172. C
  173. C Module de young et coef poisson
  174. XE1 = XMAT(1);
  175. XNU1 = XMAT(2);
  176. C
  177. C parametres alph0, c1, n1, pref, k1
  178. XE0 = XMAT(5);
  179. XA0 = XMAT(6);
  180. XC1 = XMAT(7);
  181. XN1 = XMAT(8);
  182. XPREF1 = XMAT(9);
  183. XK1 = XMAT(10);
  184. C
  185. C seuils m1, m2, m3, m4, m5, m6, m7, m8, m9, m10
  186. XM1 = XMAT(11)
  187. XM2 = XMAT(12)
  188. XM3 = XMAT(13)
  189. XM4 = XMAT(14)
  190. XM5 = XMAT(15)
  191. XM6 = XMAT(16)
  192. XM7 = XMAT(17)
  193. XM8 = XMAT(18)
  194. XM9 = XMAT(19)
  195. XM10 = XMAT(20)
  196. C
  197. C Modules plastiques Hc1, Hc2, Hc3, Hc4, Hc5, Hc6, Hc7, Hc8, Hc9
  198. C Hc10
  199. XHC1 = XMAT(21);
  200. XHC2 = XMAT(22);
  201. XHC3 = XMAT(23);
  202. XHC4 = XMAT(24);
  203. XHC5 = XMAT(25);
  204. XHC6 = XMAT(26);
  205. XHC7 = XMAT(27);
  206. XHC8 = XMAT(28);
  207. XHC9 = XMAT(29);
  208. XHC10 = 0.01D0*XMAT(29);
  209. C
  210. C Modules plastiques He1, He2, He3, He4, He5, He6, He7, He8, He9
  211. C He10
  212. XHE1 = XMAT(30);
  213. XHE2 = XMAT(31);
  214. XHE3 = XMAT(32);
  215. XHE4 = XMAT(33);
  216. XHE5 = XMAT(34);
  217. XHE6 = XMAT(35);
  218. XHE7 = XMAT(36);
  219. XHE8 = XMAT(37);
  220. XHE9 = XMAT(38);
  221. XHE10 = 0.01D0*XMAT(38);
  222. C
  223. C Modules plastiques D1, D2, D3, D4, D5, D6, D7, D8, D9, D10
  224. XD1 = XMAT(39);
  225. XD2 = XMAT(40);
  226. XD3 = XMAT(41);
  227. XD4 = XMAT(42);
  228. XD5 = XMAT(43);
  229. XD6 = XMAT(44);
  230. XD7 = XMAT(45);
  231. XD8 = XMAT(46);
  232. XD9 = XMAT(47);
  233. XD10 = 0.0D0;
  234. C
  235. C Modules plastiques eta_c1, eta_c2, eta_c3, eta_c4, eta_c5
  236. C eta_c6, eta_c7, eta_c8, eta_c9, eta_c10
  237. XETC1 = XMAT(48);
  238. XETC2 = XMAT(49);
  239. XETC3 = XMAT(50);
  240. XETC4 = XMAT(51);
  241. XETC5 = XMAT(52);
  242. XETC6 = XMAT(53);
  243. XETC7 = XMAT(54);
  244. XETC8 = XMAT(55);
  245. XETC9 = XMAT(56);
  246. XETC10 = 1000.0D0;
  247. C
  248. C Modules plastiques eta_e1, eta_e2, eta_e3, eta_e4, eta_e5
  249. C eta_e6, eta_e7, eta_e8, eta_e9, eta_e10
  250. XETE1 = XMAT(57);
  251. XETE2 = XMAT(58);
  252. XETE3 = XMAT(59);
  253. XETE4 = XMAT(60);
  254. XETE5 = XMAT(61);
  255. XETE6 = XMAT(62);
  256. XETE7 = XMAT(63);
  257. XETE8 = XMAT(64);
  258. XETE9 = XMAT(65);
  259. XETE10 = 1000.0D0;
  260. C
  261. C----CONTROLE SUR LE VALEURS MATERIAUX ---------------------------------
  262. IF ((XM2.LT.XM1).OR.(XM3.LT.XM2).OR.(XM4.LT.XM3).OR.(XM5.LT.XM4)
  263. &.OR.(XM6.LT.XM5).OR.(XM7.LT.XM6).OR.(XM8.LT.XM7).OR.(XM9.LT.XM8)
  264. &.OR.(XM10.LT.XM9)) THEN
  265. WRITE(IOIMP,10005)
  266. 10005 FORMAT(' ERREUR DANS LES VALEURS M');
  267. STOP
  268. ENDIF
  269. C
  270. C---- DEFINITION DES VALEURS ELASTIQUES --------------------------------
  271. GMOD_E1 = ((XE1)/((2.0D0)*(1.0D0 + XNU1)))
  272. KMOD_E1 = ((XE1)/((3.0D0)*(1.0D0 - ((2.0D0)*(XNU1)))))
  273. GMOD = (((XA0)*(XE0))/((2.0D0)*(1.0D0 + XNU1)))
  274. KMOD = (((XA0)*(XE0))/((3.0D0)*(1.0D0 - ((2.0D0)*(XNU1)))))
  275. C
  276. C---- DEFINITION DES TOLLERANCES ---------------------------------------
  277. XTOL1 = 1.0E-8
  278. XTOL2 = 1.0E-8
  279. C
  280. C---- INITIALLIZZATION INTEGER POUR SUBSTEPPING ------------------------
  281. IND_TEST = INT(0)
  282. IND_STEP1 = INT(0)
  283. IND_STEP2 = INT(0)
  284. IND_STEP3 = INT(0)
  285. C
  286. C---- VALEUR INITIALES -------------------------------------------------
  287. SIG_INI(1:6) = SIG0(1:6)
  288. BET1_INI(1:6) = VAR0(2:7)
  289. BET2_INI(1:6) = VAR0(8:13)
  290. BET3_INI(1:6) = VAR0(14:19)
  291. BET4_INI(1:6) = VAR0(20:25)
  292. BET5_INI(1:6) = VAR0(26:31)
  293. BET6_INI(1:6) = VAR0(32:37)
  294. BET7_INI(1:6) = VAR0(38:43)
  295. BET8_INI(1:6) = VAR0(44:49)
  296. BET9_INI(1:6) = VAR0(50:55)
  297. BET10_INI(1:6) = VAR0(56:61)
  298. C
  299. ETOT_INI(1:6) = EPST0(1:6)
  300. EPLA_INI(1:6) = EPIN0(1:6)
  301. EPLA_INI(4:6) = (EPLA_INI(4:6))/(2.0D0)
  302. ETOT_INI(4:6) = (ETOT_INI(4:6))/(2.0D0)
  303. C
  304. N_STEP = INT(1)
  305. C
  306. C Modification 17-05-2016
  307. IF (IFOUR.NE.2) THEN
  308. DEPST(5:6) = 0.0D0
  309. DDET_ST(5:6) = 0.0D0
  310. ETOT_INI(5:6) = 0.0D0
  311. EPLA_INI(5:6) = 0.0D0
  312. ENDIF
  313. DDET_ST(1:6) = (DEPST(1:6))/(N_STEP)
  314. DDET_ST(4:6) = (DDET_ST(4:6))/(2.0D0)
  315.  
  316. C
  317. C---- ON COMMENCE LA BOUCLE WHILE POUR L ALGORITHME ADAPTATIVE ---------
  318.  
  319. 10 CONTINUE
  320. DO 100 I_STEP=1,N_STEP
  321. C
  322. C---- AU DEBUT DE CHAQUE STEP ON INTIALISE LES VARIABLE ET ------------
  323. C---- ON RECHERCHE LE CRITERE ACTIVE AU DERNIER PAS DE TEMP -----------
  324. c
  325. C Delta deformation totale pour le substep
  326. DE_TOT(1:6) = DDET_ST(1:6)
  327. c
  328. C Partie deviatore et volumique du Deps_tot
  329. DSE_TOT = DE_TOT(1:6) -
  330. & (((DE_TOT(1) + DE_TOT(2) + DE_TOT(3))/(3.0D0))*(XID1(1:6)))
  331. DVE_TOT = (DE_TOT(1) + DE_TOT(2) + DE_TOT(3))
  332. c
  333. C Partie deviatoire et hydrostatique de sig_ini
  334. S_INI = SIG_INI(1:6) -
  335. & (((SIG_INI(1) + SIG_INI(2) + SIG_INI(3))/(3.0D0))*(XID1(1:6)))
  336. P_INI = ((SIG_INI(1) + SIG_INI(2) + SIG_INI(3))/(3.0D0))
  337. c
  338. C PArtie deviatoire et volumiquue de e_pla_ini
  339. SEPLA_INI(1:6) = EPLA_INI -
  340. & (((EPLA_INI(1) + EPLA_INI(2) + EPLA_INI(3))/(3.0D0))*(XID1(1:6)))
  341. VEPLA_INI = (EPLA_INI(1) + EPLA_INI(2) + EPLA_INI(3))
  342. C
  343. C On determine le termes xi_i_ini = s_ini +(p_ini - c)*beta_i_ini
  344. C pour chaque critere
  345. X1_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET1_INI(1:6)))
  346. X2_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET2_INI(1:6)))
  347. X3_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET3_INI(1:6)))
  348. X4_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET4_INI(1:6)))
  349. X5_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET5_INI(1:6)))
  350. X6_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET6_INI(1:6)))
  351. X7_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET7_INI(1:6)))
  352. X8_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET8_INI(1:6)))
  353. X9_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET9_INI(1:6)))
  354. X10_INI(1:6) = S_INI(1:6) + ((P_INI - XC1)*(BET10_INI(1:6)))
  355. C
  356. C On determine lA norme |x_i_ini|
  357. NX1_INI = (((X1_INI(1)*X1_INI(1)) +
  358. & (X1_INI(2)*X1_INI(2)) +
  359. & (X1_INI(3)*X1_INI(3)) +
  360. & ((2.0D0)*(X1_INI(4)*X1_INI(4))) +
  361. & ((2.0D0)*(X1_INI(5)*X1_INI(5))) +
  362. & ((2.0D0)*(X1_INI(6)*X1_INI(6))))**(0.5D0))
  363. NX2_INI = (((X2_INI(1)*X2_INI(1)) +
  364. & (X2_INI(2)*X2_INI(2)) +
  365. & (X2_INI(3)*X2_INI(3)) +
  366. & ((2.0D0)*(X2_INI(4)*X2_INI(4))) +
  367. & ((2.0D0)*(X2_INI(5)*X2_INI(5))) +
  368. & ((2.0D0)*(X2_INI(6)*X2_INI(6))))**(0.5D0))
  369. NX3_INI = (((X3_INI(1)*X3_INI(1)) +
  370. & (X3_INI(2)*X3_INI(2)) +
  371. & (X3_INI(3)*X3_INI(3)) +
  372. & ((2.0D0)*(X3_INI(4)*X3_INI(4))) +
  373. & ((2.0D0)*(X3_INI(5)*X3_INI(5))) +
  374. & ((2.0D0)*(X3_INI(6)*X3_INI(6))))**(0.5D0))
  375. NX4_INI = (((X4_INI(1)*X4_INI(1)) +
  376. & (X4_INI(2)*X4_INI(2)) +
  377. & (X4_INI(3)*X4_INI(3)) +
  378. & ((2.0D0)*(X4_INI(4)*X4_INI(4))) +
  379. & ((2.0D0)*(X4_INI(5)*X4_INI(5))) +
  380. & ((2.0D0)*(X4_INI(6)*X4_INI(6))))**(0.5D0))
  381. NX5_INI = (((X5_INI(1)*X5_INI(1)) +
  382. & (X5_INI(2)*X5_INI(2)) +
  383. & (X5_INI(3)*X5_INI(3)) +
  384. & ((2.0D0)*(X5_INI(4)*X5_INI(4))) +
  385. & ((2.0D0)*(X5_INI(5)*X5_INI(5))) +
  386. & ((2.0D0)*(X5_INI(6)*X5_INI(6))))**(0.5D0))
  387. NX6_INI = (((X6_INI(1)*X6_INI(1)) +
  388. & (X6_INI(2)*X6_INI(2)) +
  389. & (X6_INI(3)*X6_INI(3)) +
  390. & ((2.0D0)*(X6_INI(4)*X6_INI(4))) +
  391. & ((2.0D0)*(X6_INI(5)*X6_INI(5))) +
  392. & ((2.0D0)*(X6_INI(6)*X6_INI(6))))**(0.5D0))
  393. NX7_INI = (((X7_INI(1)*X7_INI(1)) +
  394. & (X7_INI(2)*X7_INI(2)) +
  395. & (X7_INI(3)*X7_INI(3)) +
  396. & ((2.0D0)*(X7_INI(4)*X7_INI(4))) +
  397. & ((2.0D0)*(X7_INI(5)*X7_INI(5))) +
  398. & ((2.0D0)*(X7_INI(6)*X7_INI(6))))**(0.5D0))
  399. NX8_INI = (((X8_INI(1)*X8_INI(1)) +
  400. & (X8_INI(2)*X8_INI(2)) +
  401. & (X8_INI(3)*X8_INI(3)) +
  402. & ((2.0D0)*(X8_INI(4)*X8_INI(4))) +
  403. & ((2.0D0)*(X8_INI(5)*X8_INI(5))) +
  404. & ((2.0D0)*(X8_INI(6)*X8_INI(6))))**(0.5D0))
  405. NX9_INI = (((X9_INI(1)*X9_INI(1)) +
  406. & (X9_INI(2)*X9_INI(2)) +
  407. & (X9_INI(3)*X9_INI(3)) +
  408. & ((2.0D0)*(X9_INI(4)*X9_INI(4))) +
  409. & ((2.0D0)*(X9_INI(5)*X9_INI(5))) +
  410. & ((2.0D0)*(X9_INI(6)*X9_INI(6))))**(0.5D0))
  411. NX10_INI = (((X10_INI(1)*X10_INI(1)) +
  412. & (X10_INI(2)*X10_INI(2)) +
  413. & (X10_INI(3)*X10_INI(3)) +
  414. & ((2.0D0)*(X10_INI(4)*X10_INI(4))) +
  415. & ((2.0D0)*(X10_INI(5)*X10_INI(5))) +
  416. & ((2.0D0)*(X10_INI(6)*X10_INI(6))))**(0.5D0))
  417. C
  418. C Check sur |x_i_ini| et determination de nn_i_ini de nn2_i_ini
  419. IF ((NX1_INI).LT.(1.0E-12) ) THEN
  420. NN1_INI(1:6) = 0.0D0
  421. NN21_INI(1:6) = 0.0D0
  422. ELSE
  423. NN1_INI(1:6) = X1_INI(1:6)/(NX1_INI)
  424. NN21_INI(1) = ((NN1_INI(1))*(NN1_INI(1))) +
  425. & ((NN1_INI(4))*(NN1_INI(4))) +
  426. & ((NN1_INI(5))*(NN1_INI(5)))
  427. NN21_INI(2) = ((NN1_INI(4))*(NN1_INI(4))) +
  428. & ((NN1_INI(2))*(NN1_INI(2))) +
  429. & ((NN1_INI(6))*(NN1_INI(6)))
  430. NN21_INI(3) = ((NN1_INI(5))*(NN1_INI(5))) +
  431. & ((NN1_INI(6))*(NN1_INI(6))) +
  432. & ((NN1_INI(3))*(NN1_INI(3)))
  433. NN21_INI(4) = ((NN1_INI(1))*(NN1_INI(4))) +
  434. & ((NN1_INI(4))*(NN1_INI(2))) +
  435. & ((NN1_INI(5))*(NN1_INI(6)))
  436. NN21_INI(5) = ((NN1_INI(1))*(NN1_INI(5))) +
  437. & ((NN1_INI(4))*(NN1_INI(6))) +
  438. & ((NN1_INI(5))*(NN1_INI(3)))
  439. NN21_INI(6) = ((NN1_INI(4))*(NN1_INI(5))) +
  440. & ((NN1_INI(2))*(NN1_INI(6))) +
  441. & ((NN1_INI(6))*(NN1_INI(3)))
  442. ENDIF
  443. IF ((NX2_INI).LT.(1.0E-12) ) THEN
  444. NN2_INI(1:6) = 0.0D0
  445. NN22_INI(1:6) = 0.0D0
  446. ELSE
  447. NN2_INI(1:6) = X2_INI(1:6)/(NX2_INI)
  448. NN22_INI(1) = ((NN2_INI(1))*(NN2_INI(1))) +
  449. & ((NN2_INI(4))*(NN2_INI(4))) +
  450. & ((NN2_INI(5))*(NN2_INI(5)))
  451. NN22_INI(2) = ((NN2_INI(4))*(NN2_INI(4))) +
  452. & ((NN2_INI(2))*(NN2_INI(2))) +
  453. & ((NN2_INI(6))*(NN2_INI(6)))
  454. NN22_INI(3) = ((NN2_INI(5))*(NN2_INI(5))) +
  455. & ((NN2_INI(6))*(NN2_INI(6))) +
  456. & ((NN2_INI(3))*(NN2_INI(3)))
  457. NN22_INI(4) = ((NN2_INI(1))*(NN2_INI(4))) +
  458. & ((NN2_INI(4))*(NN2_INI(2))) +
  459. & ((NN2_INI(5))*(NN2_INI(6)))
  460. NN22_INI(5) = ((NN2_INI(1))*(NN2_INI(5))) +
  461. & ((NN2_INI(4))*(NN2_INI(6))) +
  462. & ((NN2_INI(5))*(NN2_INI(3)))
  463. NN22_INI(6) = ((NN2_INI(4))*(NN2_INI(5))) +
  464. & ((NN2_INI(2))*(NN2_INI(6))) +
  465. & ((NN2_INI(6))*(NN2_INI(3)))
  466. ENDIF
  467. IF ((NX3_INI).LT.(1.0E-12) ) THEN
  468. NN3_INI(1:6) = 0.0D0
  469. NN23_INI(1:6) = 0.0D0
  470. ELSE
  471. NN3_INI(1:6) = X3_INI(1:6)/(NX3_INI)
  472. NN23_INI(1) = ((NN3_INI(1))*(NN3_INI(1))) +
  473. & ((NN3_INI(4))*(NN3_INI(4))) +
  474. & ((NN3_INI(5))*(NN3_INI(5)))
  475. NN23_INI(2) = ((NN3_INI(4))*(NN3_INI(4))) +
  476. & ((NN3_INI(2))*(NN3_INI(2))) +
  477. & ((NN3_INI(6))*(NN3_INI(6)))
  478. NN23_INI(3) = ((NN3_INI(5))*(NN3_INI(5))) +
  479. & ((NN3_INI(6))*(NN3_INI(6))) +
  480. & ((NN3_INI(3))*(NN3_INI(3)))
  481. NN23_INI(4) = ((NN3_INI(1))*(NN3_INI(4))) +
  482. & ((NN3_INI(4))*(NN3_INI(2))) +
  483. & ((NN3_INI(5))*(NN3_INI(6)))
  484. NN23_INI(5) = ((NN3_INI(1))*(NN3_INI(5))) +
  485. & ((NN3_INI(4))*(NN3_INI(6))) +
  486. & ((NN3_INI(5))*(NN3_INI(3)))
  487. NN23_INI(6) = ((NN3_INI(4))*(NN3_INI(5))) +
  488. & ((NN3_INI(2))*(NN3_INI(6))) +
  489. & ((NN3_INI(6))*(NN3_INI(3)))
  490. ENDIF
  491. IF ((NX4_INI).LT.(1.0E-12) ) THEN
  492. NN4_INI(1:6) = 0.0D0
  493. NN24_INI(1:6) = 0.0D0
  494. ELSE
  495. NN4_INI(1:6) = X4_INI(1:6)/(NX4_INI)
  496. NN24_INI(1) = ((NN4_INI(1))*(NN4_INI(1))) +
  497. & ((NN4_INI(4))*(NN4_INI(4))) +
  498. & ((NN4_INI(5))*(NN4_INI(5)))
  499. NN24_INI(2) = ((NN4_INI(4))*(NN4_INI(4))) +
  500. & ((NN4_INI(2))*(NN4_INI(2))) +
  501. & ((NN4_INI(6))*(NN4_INI(6)))
  502. NN24_INI(3) = ((NN4_INI(5))*(NN4_INI(5))) +
  503. & ((NN4_INI(6))*(NN4_INI(6))) +
  504. & ((NN4_INI(3))*(NN4_INI(3)))
  505. NN24_INI(4) = ((NN4_INI(1))*(NN4_INI(4))) +
  506. & ((NN4_INI(4))*(NN4_INI(2))) +
  507. & ((NN4_INI(5))*(NN4_INI(6)))
  508. NN24_INI(5) = ((NN4_INI(1))*(NN4_INI(5))) +
  509. & ((NN4_INI(4))*(NN4_INI(6))) +
  510. & ((NN4_INI(5))*(NN4_INI(3)))
  511. NN24_INI(6) = ((NN4_INI(4))*(NN4_INI(5))) +
  512. & ((NN4_INI(2))*(NN4_INI(6))) +
  513. & ((NN4_INI(6))*(NN4_INI(3)))
  514. ENDIF
  515. IF ((NX5_INI).LT.(1.0E-12) ) THEN
  516. NN5_INI(1:6) = 0.0D0
  517. NN25_INI(1:6) = 0.0D0
  518. ELSE
  519. NN5_INI(1:6) = X5_INI(1:6)/(NX5_INI)
  520. NN25_INI(1) = ((NN5_INI(1))*(NN5_INI(1))) +
  521. & ((NN5_INI(4))*(NN5_INI(4))) +
  522. & ((NN5_INI(5))*(NN5_INI(5)))
  523. NN25_INI(2) = ((NN5_INI(4))*(NN5_INI(4))) +
  524. & ((NN5_INI(2))*(NN5_INI(2))) +
  525. & ((NN5_INI(6))*(NN5_INI(6)))
  526. NN25_INI(3) = ((NN5_INI(5))*(NN5_INI(5))) +
  527. & ((NN5_INI(6))*(NN5_INI(6))) +
  528. & ((NN5_INI(3))*(NN5_INI(3)))
  529. NN25_INI(4) = ((NN5_INI(1))*(NN5_INI(4))) +
  530. & ((NN5_INI(4))*(NN5_INI(2))) +
  531. & ((NN5_INI(5))*(NN5_INI(6)))
  532. NN25_INI(5) = ((NN5_INI(1))*(NN5_INI(5))) +
  533. & ((NN5_INI(4))*(NN5_INI(6))) +
  534. & ((NN5_INI(5))*(NN5_INI(3)))
  535. NN25_INI(6) = ((NN5_INI(4))*(NN5_INI(5))) +
  536. & ((NN5_INI(2))*(NN5_INI(6))) +
  537. & ((NN5_INI(6))*(NN5_INI(3)))
  538. ENDIF
  539. IF ((NX6_INI).LT.(1.0E-12) ) THEN
  540. NN6_INI(1:6) = 0.0D0
  541. NN26_INI(1:6) = 0.0D0
  542. ELSE
  543. NN6_INI(1:6) = X6_INI(1:6)/(NX6_INI)
  544. NN26_INI(1) = ((NN6_INI(1))*(NN6_INI(1))) +
  545. & ((NN6_INI(4))*(NN6_INI(4))) +
  546. & ((NN6_INI(5))*(NN6_INI(5)))
  547. NN26_INI(2) = ((NN6_INI(4))*(NN6_INI(4))) +
  548. & ((NN6_INI(2))*(NN6_INI(2))) +
  549. & ((NN6_INI(6))*(NN6_INI(6)))
  550. NN26_INI(3) = ((NN6_INI(5))*(NN6_INI(5))) +
  551. & ((NN6_INI(6))*(NN6_INI(6))) +
  552. & ((NN6_INI(3))*(NN6_INI(3)))
  553. NN26_INI(4) = ((NN6_INI(1))*(NN6_INI(4))) +
  554. & ((NN6_INI(4))*(NN6_INI(2))) +
  555. & ((NN6_INI(5))*(NN6_INI(6)))
  556. NN26_INI(5) = ((NN6_INI(1))*(NN6_INI(5))) +
  557. & ((NN6_INI(4))*(NN6_INI(6))) +
  558. & ((NN6_INI(5))*(NN6_INI(3)))
  559. NN26_INI(6) = ((NN6_INI(4))*(NN6_INI(5))) +
  560. & ((NN6_INI(2))*(NN6_INI(6))) +
  561. & ((NN6_INI(6))*(NN6_INI(3)))
  562. ENDIF
  563. IF ((NX7_INI).LT.(1.0E-12) ) THEN
  564. NN7_INI(1:6) = 0.0D0
  565. NN27_INI(1:6) = 0.0D0
  566. ELSE
  567. NN7_INI(1:6) = X7_INI(1:6)/(NX7_INI)
  568. NN27_INI(1) = ((NN7_INI(1))*(NN7_INI(1))) +
  569. & ((NN7_INI(4))*(NN7_INI(4))) +
  570. & ((NN7_INI(5))*(NN7_INI(5)))
  571. NN27_INI(2) = ((NN7_INI(4))*(NN7_INI(4))) +
  572. & ((NN7_INI(2))*(NN7_INI(2))) +
  573. & ((NN7_INI(6))*(NN7_INI(6)))
  574. NN27_INI(3) = ((NN7_INI(5))*(NN7_INI(5))) +
  575. & ((NN7_INI(6))*(NN7_INI(6))) +
  576. & ((NN7_INI(3))*(NN7_INI(3)))
  577. NN27_INI(4) = ((NN7_INI(1))*(NN7_INI(4))) +
  578. & ((NN7_INI(4))*(NN7_INI(2))) +
  579. & ((NN7_INI(5))*(NN7_INI(6)))
  580. NN27_INI(5) = ((NN7_INI(1))*(NN7_INI(5))) +
  581. & ((NN7_INI(4))*(NN7_INI(6))) +
  582. & ((NN7_INI(5))*(NN7_INI(3)))
  583. NN27_INI(6) = ((NN7_INI(4))*(NN7_INI(5))) +
  584. & ((NN7_INI(2))*(NN7_INI(6))) +
  585. & ((NN7_INI(6))*(NN7_INI(3)))
  586. ENDIF
  587. IF ((NX8_INI).LT.(1.0E-12) ) THEN
  588. NN8_INI(1:6) = 0.0D0
  589. NN28_INI(1:6) = 0.0D0
  590. ELSE
  591. NN8_INI(1:6) = X8_INI(1:6)/(NX8_INI)
  592. NN28_INI(1) = ((NN8_INI(1))*(NN8_INI(1))) +
  593. & ((NN8_INI(4))*(NN8_INI(4))) +
  594. & ((NN8_INI(5))*(NN8_INI(5)))
  595. NN28_INI(2) = ((NN8_INI(4))*(NN8_INI(4))) +
  596. & ((NN8_INI(2))*(NN8_INI(2))) +
  597. & ((NN8_INI(6))*(NN8_INI(6)))
  598. NN28_INI(3) = ((NN8_INI(5))*(NN8_INI(5))) +
  599. & ((NN8_INI(6))*(NN8_INI(6))) +
  600. & ((NN8_INI(3))*(NN8_INI(3)))
  601. NN28_INI(4) = ((NN8_INI(1))*(NN8_INI(4))) +
  602. & ((NN8_INI(4))*(NN8_INI(2))) +
  603. & ((NN8_INI(5))*(NN8_INI(6)))
  604. NN28_INI(5) = ((NN8_INI(1))*(NN8_INI(5))) +
  605. & ((NN8_INI(4))*(NN8_INI(6))) +
  606. & ((NN8_INI(5))*(NN8_INI(3)))
  607. NN28_INI(6) = ((NN8_INI(4))*(NN8_INI(5))) +
  608. & ((NN8_INI(2))*(NN8_INI(6))) +
  609. & ((NN8_INI(6))*(NN8_INI(3)))
  610. ENDIF
  611. IF ((NX9_INI).LT.(1.0E-12) ) THEN
  612. NN9_INI(1:6) = 0.0D0
  613. NN29_INI(1:6) = 0.0D0
  614. ELSE
  615. NN9_INI(1:6) = X9_INI(1:6)/(NX9_INI)
  616. NN29_INI(1) = ((NN9_INI(1))*(NN9_INI(1))) +
  617. & ((NN9_INI(4))*(NN9_INI(4))) +
  618. & ((NN9_INI(5))*(NN9_INI(5)))
  619. NN29_INI(2) = ((NN9_INI(4))*(NN9_INI(4))) +
  620. & ((NN9_INI(2))*(NN9_INI(2))) +
  621. & ((NN9_INI(6))*(NN9_INI(6)))
  622. NN29_INI(3) = ((NN9_INI(5))*(NN9_INI(5))) +
  623. & ((NN9_INI(6))*(NN9_INI(6))) +
  624. & ((NN9_INI(3))*(NN9_INI(3)))
  625. NN29_INI(4) = ((NN9_INI(1))*(NN9_INI(4))) +
  626. & ((NN9_INI(4))*(NN9_INI(2))) +
  627. & ((NN9_INI(5))*(NN9_INI(6)))
  628. NN29_INI(5) = ((NN9_INI(1))*(NN9_INI(5))) +
  629. & ((NN9_INI(4))*(NN9_INI(6))) +
  630. & ((NN9_INI(5))*(NN9_INI(3)))
  631. NN29_INI(6) = ((NN9_INI(4))*(NN9_INI(5))) +
  632. & ((NN9_INI(2))*(NN9_INI(6))) +
  633. & ((NN9_INI(6))*(NN9_INI(3)))
  634. ENDIF
  635. IF ((NX10_INI).LT.(1.0E-12) ) THEN
  636. NN10_INI(1:6) = 0.0D0
  637. NN210_INI(1:6) = 0.0D0
  638. ELSE
  639. NN10_INI(1:6) = X10_INI(1:6)/(NX10_INI)
  640. NN210_INI(1) = ((NN10_INI(1))*(NN10_INI(1))) +
  641. & ((NN10_INI(4))*(NN10_INI(4))) +
  642. & ((NN10_INI(5))*(NN10_INI(5)))
  643. NN210_INI(2) = ((NN10_INI(4))*(NN10_INI(4))) +
  644. & ((NN10_INI(2))*(NN10_INI(2))) +
  645. & ((NN10_INI(6))*(NN10_INI(6)))
  646. NN210_INI(3) = ((NN10_INI(5))*(NN10_INI(5))) +
  647. & ((NN10_INI(6))*(NN10_INI(6))) +
  648. & ((NN10_INI(3))*(NN10_INI(3)))
  649. NN210_INI(4) = ((NN10_INI(1))*(NN10_INI(4))) +
  650. & ((NN10_INI(4))*(NN10_INI(2))) +
  651. & ((NN10_INI(5))*(NN10_INI(6)))
  652. NN210_INI(5) = ((NN10_INI(1))*(NN10_INI(5))) +
  653. & ((NN10_INI(4))*(NN10_INI(6))) +
  654. & ((NN10_INI(5))*(NN10_INI(3)))
  655. NN210_INI(6) = ((NN10_INI(4))*(NN10_INI(5))) +
  656. & ((NN10_INI(2))*(NN10_INI(6))) +
  657. & ((NN10_INI(6))*(NN10_INI(3)))
  658. ENDIF
  659. C
  660. C Calcul de cos(3th_i_ini)=-sqrt(6)*trace(nn_i_ini^3) -
  661. c depandance angle de Lode
  662. COS3TH1_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  663. & (((NN1_INI(1))*(NN21_INI(1))) +
  664. & ((NN1_INI(2))*(NN21_INI(2))) +
  665. & ((NN1_INI(3))*(NN21_INI(3))) +
  666. & (2.0D0)*((NN1_INI(4))*(NN21_INI(4))) +
  667. & (2.0D0)*((NN1_INI(5))*(NN21_INI(5))) +
  668. & (2.0D0)*((NN1_INI(6))*(NN21_INI(6))))
  669. COS3TH2_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  670. & (((NN2_INI(1))*(NN22_INI(1))) +
  671. & ((NN2_INI(2))*(NN22_INI(2))) +
  672. & ((NN2_INI(3))*(NN22_INI(3))) +
  673. & (2.0D0)*((NN2_INI(4))*(NN22_INI(4))) +
  674. & (2.0D0)*((NN2_INI(5))*(NN22_INI(5))) +
  675. & (2.0D0)*((NN2_INI(6))*(NN22_INI(6))))
  676. COS3TH3_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  677. & (((NN3_INI(1))*(NN23_INI(1))) +
  678. & ((NN3_INI(2))*(NN23_INI(2))) +
  679. & ((NN3_INI(3))*(NN23_INI(3))) +
  680. & (2.0D0)*((NN3_INI(4))*(NN23_INI(4))) +
  681. & (2.0D0)*((NN3_INI(5))*(NN23_INI(5))) +
  682. & (2.0D0)*((NN3_INI(6))*(NN23_INI(6))))
  683. COS3TH4_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  684. & (((NN4_INI(1))*(NN24_INI(1))) +
  685. & ((NN4_INI(2))*(NN24_INI(2))) +
  686. & ((NN4_INI(3))*(NN24_INI(3))) +
  687. & (2.0D0)*((NN4_INI(4))*(NN24_INI(4))) +
  688. & (2.0D0)*((NN4_INI(5))*(NN24_INI(5))) +
  689. & (2.0D0)*((NN4_INI(6))*(NN24_INI(6))))
  690. COS3TH5_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  691. & (((NN5_INI(1))*(NN25_INI(1))) +
  692. & ((NN5_INI(2))*(NN25_INI(2))) +
  693. & ((NN5_INI(3))*(NN25_INI(3))) +
  694. & (2.0D0)*((NN5_INI(4))*(NN25_INI(4))) +
  695. & (2.0D0)*((NN5_INI(5))*(NN25_INI(5))) +
  696. & (2.0D0)*((NN5_INI(6))*(NN25_INI(6))))
  697. COS3TH6_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  698. & (((NN6_INI(1))*(NN26_INI(1))) +
  699. & ((NN6_INI(2))*(NN26_INI(2))) +
  700. & ((NN6_INI(3))*(NN26_INI(3))) +
  701. & (2.0D0)*((NN6_INI(4))*(NN26_INI(4))) +
  702. & (2.0D0)*((NN6_INI(5))*(NN26_INI(5))) +
  703. & (2.0D0)*((NN6_INI(6))*(NN26_INI(6))))
  704. COS3TH7_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  705. & (((NN7_INI(1))*(NN27_INI(1))) +
  706. & ((NN7_INI(2))*(NN27_INI(2))) +
  707. & ((NN7_INI(3))*(NN27_INI(3))) +
  708. & (2.0D0)*((NN7_INI(4))*(NN27_INI(4))) +
  709. & (2.0D0)*((NN7_INI(5))*(NN27_INI(5))) +
  710. & (2.0D0)*((NN7_INI(6))*(NN27_INI(6))))
  711. COS3TH8_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  712. & (((NN8_INI(1))*(NN28_INI(1))) +
  713. & ((NN8_INI(2))*(NN28_INI(2))) +
  714. & ((NN8_INI(3))*(NN28_INI(3))) +
  715. & (2.0D0)*((NN8_INI(4))*(NN28_INI(4))) +
  716. & (2.0D0)*((NN8_INI(5))*(NN28_INI(5))) +
  717. & (2.0D0)*((NN8_INI(6))*(NN28_INI(6))))
  718. COS3TH9_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  719. & (((NN9_INI(1))*(NN29_INI(1))) +
  720. & ((NN9_INI(2))*(NN29_INI(2))) +
  721. & ((NN9_INI(3))*(NN29_INI(3))) +
  722. & (2.0D0)*((NN9_INI(4))*(NN29_INI(4))) +
  723. & (2.0D0)*((NN9_INI(5))*(NN29_INI(5))) +
  724. & (2.0D0)*((NN9_INI(6))*(NN29_INI(6))))
  725. COS3TH10_INI = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  726. & (((NN10_INI(1))*(NN210_INI(1))) +
  727. & ((NN10_INI(2))*(NN210_INI(2))) +
  728. & ((NN10_INI(3))*(NN210_INI(3))) +
  729. & (2.0D0)*((NN10_INI(4))*(NN210_INI(4))) +
  730. & (2.0D0)*((NN10_INI(5))*(NN210_INI(5))) +
  731. & (2.0D0)*((NN10_INI(6))*(NN210_INI(6))))
  732. C
  733. C Calcul du terme R(th)_i_ini=2k1/((1+k1)-(1-k1)*(cos(3th_i_ini)))
  734. RTH1_INI = ((2.0D0)*(XK1))/
  735. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH1_INI)))
  736. RTH2_INI = ((2.0D0)*(XK1))/
  737. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH2_INI)))
  738. RTH3_INI = ((2.0D0)*(XK1))/
  739. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH3_INI)))
  740. RTH4_INI = ((2.0D0)*(XK1))/
  741. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH4_INI)))
  742. RTH5_INI = ((2.0D0)*(XK1))/
  743. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH5_INI)))
  744. RTH6_INI = ((2.0D0)*(XK1))/
  745. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH6_INI)))
  746. RTH7_INI = ((2.0D0)*(XK1))/
  747. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH7_INI)))
  748. RTH8_INI = ((2.0D0)*(XK1))/
  749. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH8_INI)))
  750. RTH9_INI = ((2.0D0)*(XK1))/
  751. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH9_INI)))
  752. RTH10_INI = ((2.0D0)*(XK1))/
  753. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH10_INI)))
  754. C
  755. C Calcul des criteres f_i_ini
  756. F1_INI = NX1_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM1))*
  757. & ((P_INI - XC1)*(RTH1_INI))
  758. F2_INI = NX2_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM2))*
  759. & ((P_INI - XC1)*(RTH2_INI))
  760. F3_INI = NX3_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM3))*
  761. & ((P_INI - XC1)*(RTH3_INI))
  762. F4_INI = NX4_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM4))*
  763. & ((P_INI - XC1)*(RTH4_INI))
  764. F5_INI = NX5_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM5))*
  765. & ((P_INI - XC1)*(RTH5_INI))
  766. F6_INI = NX6_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM6))*
  767. & ((P_INI - XC1)*(RTH6_INI))
  768. F7_INI = NX7_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM7))*
  769. & ((P_INI - XC1)*(RTH7_INI))
  770. F8_INI = NX8_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM8))*
  771. & ((P_INI - XC1)*(RTH8_INI))
  772. F9_INI = NX9_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM9))*
  773. & ((P_INI - XC1)*(RTH9_INI))
  774. F10_INI = NX10_INI + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM10))*
  775. & ((P_INI - XC1)*(RTH10_INI))
  776. C
  777. C Determination du critere activé, si
  778. IINI_CRIT = INT(0)
  779. IF ((F1_INI.GT.(0.0D0)).OR.((ABS(F1_INI)).LT.(1.0E-8))) THEN
  780. IINI_CRIT = INT(1);
  781. ENDIF
  782. IF ((F2_INI.GT.(0.0D0)).OR.((ABS(F2_INI)).LT.(1.0E-8))) THEN
  783. IINI_CRIT = INT(2);
  784. ENDIF
  785. IF ((F3_INI.GT.(0.0D0)).OR.((ABS(F3_INI)).LT.(1.0E-8))) THEN
  786. IINI_CRIT = INT(3);
  787. ENDIF
  788. IF ((F4_INI.GT.(0.0D0)).OR.((ABS(F4_INI)).LT.(1.0E-8))) THEN
  789. IINI_CRIT = INT(4);
  790. ENDIF
  791. IF ((F5_INI.GT.(0.0D0)).OR.((ABS(F5_INI)).LT.(1.0E-8))) THEN
  792. IINI_CRIT = INT(5);
  793. ENDIF
  794. IF ((F6_INI.GT.(0.0D0)).OR.((ABS(F6_INI)).LT.(1.0E-8))) THEN
  795. IINI_CRIT = INT(6);
  796. ENDIF
  797. IF ((F7_INI.GT.(0.0D0)).OR.((ABS(F7_INI)).LT.(1.0E-8))) THEN
  798. IINI_CRIT = INT(7);
  799. ENDIF
  800. IF ((F8_INI.GT.(0.0D0)).OR.((ABS(F8_INI)).LT.(1.0E-8))) THEN
  801. IINI_CRIT = INT(8);
  802. ENDIF
  803. IF ((F9_INI.GT.(0.0D0)).OR.((ABS(F9_INI)).LT.(1.0E-8))) THEN
  804. IINI_CRIT = INT(9);
  805. ENDIF
  806. IF ((F10_INI.GT.(0.0D0)).OR.((ABS(F10_INI)).LT.(1.0E-8))) THEN
  807. IINI_CRIT = INT(10);
  808. ENDIF
  809. C---- FIN DE L INTIALISE LES VARIABLE ET DE LA -------------------------
  810. C---- ON RECHERCHE LE CRITERE ACTIVE AU DERNIER PAS DE TEMP -----------
  811. C
  812. C------ ITERATION POUR LA RECHERCHE DE L'ETAT PLASTIQUE ADMISSIBLE -----
  813. C
  814. C Initialissation deformations plastique
  815. SEPLA_0(1:6) = SEPLA_INI(1:6)
  816. VEPLA_0 = VEPLA_INI
  817. C
  818. C Initialissation back stress
  819. BET1_0(1:6) = BET1_INI(1:6)
  820. BET2_0(1:6) = BET2_INI(1:6)
  821. BET3_0(1:6) = BET3_INI(1:6)
  822. BET4_0(1:6) = BET4_INI(1:6)
  823. BET5_0(1:6) = BET5_INI(1:6)
  824. BET6_0(1:6) = BET6_INI(1:6)
  825. BET7_0(1:6) = BET7_INI(1:6)
  826. BET8_0(1:6) = BET8_INI(1:6)
  827. BET9_0(1:6) = BET9_INI(1:6)
  828. BET10_0(1:6) = BET10_INI(1:6)
  829. C
  830. C Initialissation Gamma_plastique
  831. DGAM_0 = 0.0D0
  832. C
  833. C Prevision Elastique des contraintes (etat trial)
  834. S_0(1:6) = S_INI(1:6) +
  835. & ((2.0D0)*(GMOD_E1))*(DSE_TOT(1:6) -
  836. & (SEPLA_0(1:6) - SEPLA_INI(1:6)))
  837. IF (IFOUR.NE.2) THEN
  838. S_0(5:6) = 0.0D0
  839. ENDIF
  840. P_0 = P_INI +
  841. & ((1.0D0)*(KMOD_E1))*(DVE_TOT -
  842. & (VEPLA_0 - VEPLA_INI))
  843. C
  844. N_ITER1 = 1000;
  845. C
  846. C------ CICLE DE CONVERGENCE POUR CHAQUE SUBSTEP -----------------------
  847. DO 101 I_ITER1=1,N_ITER1
  848. C Cas apres la premiere iteration k=1
  849. IF (I_ITER1.GT.(INT(1))) THEN
  850. C Mise a jour pour l'iteration k+1
  851. SEPLA_0(1:6) = SEPLA_1(1:6)
  852. VEPLA_0 = VEPLA_1
  853. BET_0(1:6) = BET_1(1:6)
  854. DGAM_0 = DGAM_1
  855. ELSE
  856. C Choix du critere activé des variables pour l'evolution
  857. C de l'etat plastique
  858. C
  859. C On determine le termes xi_i_0 = s_0 +(p_0 - c)*beta_i_0
  860. C pour chaque critere pour k=1
  861. X1_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET1_0(1:6)))
  862. X2_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET2_0(1:6)))
  863. X3_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET3_0(1:6)))
  864. X4_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET4_0(1:6)))
  865. X5_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET5_0(1:6)))
  866. X6_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET6_0(1:6)))
  867. X7_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET7_0(1:6)))
  868. X8_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET8_0(1:6)))
  869. X9_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET9_0(1:6)))
  870. X10_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET10_0(1:6)))
  871. C
  872. C On determine lA norme |x_i_ini|
  873. NX1_0 = (((X1_0(1)*X1_0(1)) +
  874. & (X1_0(2)*X1_0(2)) +
  875. & (X1_0(3)*X1_0(3)) +
  876. & ((2.0D0)*(X1_0(4)*X1_0(4))) +
  877. & ((2.0D0)*(X1_0(5)*X1_0(5))) +
  878. & ((2.0D0)*(X1_0(6)*X1_0(6))))**(0.5D0))
  879. NX2_0 = (((X2_0(1)*X2_0(1)) +
  880. & (X2_0(2)*X2_0(2)) +
  881. & (X2_0(3)*X2_0(3)) +
  882. & ((2.0D0)*(X2_0(4)*X2_0(4))) +
  883. & ((2.0D0)*(X2_0(5)*X2_0(5))) +
  884. & ((2.0D0)*(X2_0(6)*X2_0(6))))**(0.5D0))
  885. NX3_0 = (((X3_0(1)*X3_0(1)) +
  886. & (X3_0(2)*X3_0(2)) +
  887. & (X3_0(3)*X3_0(3)) +
  888. & ((2.0D0)*(X3_0(4)*X3_0(4))) +
  889. & ((2.0D0)*(X3_0(5)*X3_0(5))) +
  890. & ((2.0D0)*(X3_0(6)*X3_0(6))))**(0.5D0))
  891. NX4_0 = (((X4_0(1)*X4_0(1)) +
  892. & (X4_0(2)*X4_0(2)) +
  893. & (X4_0(3)*X4_0(3)) +
  894. & ((2.0D0)*(X4_0(4)*X4_0(4))) +
  895. & ((2.0D0)*(X4_0(5)*X4_0(5))) +
  896. & ((2.0D0)*(X4_0(6)*X4_0(6))))**(0.5D0))
  897. NX5_0 = (((X5_0(1)*X5_0(1)) +
  898. & (X5_0(2)*X5_0(2)) +
  899. & (X5_0(3)*X5_0(3)) +
  900. & ((2.0D0)*(X5_0(4)*X5_0(4))) +
  901. & ((2.0D0)*(X5_0(5)*X5_0(5))) +
  902. & ((2.0D0)*(X5_0(6)*X5_0(6))))**(0.5D0))
  903. NX6_0 = (((X6_0(1)*X6_0(1)) +
  904. & (X6_0(2)*X6_0(2)) +
  905. & (X6_0(3)*X6_0(3)) +
  906. & ((2.0D0)*(X6_0(4)*X6_0(4))) +
  907. & ((2.0D0)*(X6_0(5)*X6_0(5))) +
  908. & ((2.0D0)*(X6_0(6)*X6_0(6))))**(0.5D0))
  909. NX7_0 = (((X7_0(1)*X7_0(1)) +
  910. & (X7_0(2)*X7_0(2)) +
  911. & (X7_0(3)*X7_0(3)) +
  912. & ((2.0D0)*(X7_0(4)*X7_0(4))) +
  913. & ((2.0D0)*(X7_0(5)*X7_0(5))) +
  914. & ((2.0D0)*(X7_0(6)*X7_0(6))))**(0.5D0))
  915. NX8_0 = (((X8_0(1)*X8_0(1)) +
  916. & (X8_0(2)*X8_0(2)) +
  917. & (X8_0(3)*X8_0(3)) +
  918. & ((2.0D0)*(X8_0(4)*X8_0(4))) +
  919. & ((2.0D0)*(X8_0(5)*X8_0(5))) +
  920. & ((2.0D0)*(X8_0(6)*X8_0(6))))**(0.5D0))
  921. NX9_0 = (((X9_0(1)*X9_0(1)) +
  922. & (X9_0(2)*X9_0(2)) +
  923. & (X9_0(3)*X9_0(3)) +
  924. & ((2.0D0)*(X9_0(4)*X9_0(4))) +
  925. & ((2.0D0)*(X9_0(5)*X9_0(5))) +
  926. & ((2.0D0)*(X9_0(6)*X9_0(6))))**(0.5D0))
  927. NX10_0 = (((X10_0(1)*X10_0(1)) +
  928. & (X10_0(2)*X10_0(2)) +
  929. & (X10_0(3)*X10_0(3)) +
  930. & ((2.0D0)*(X10_0(4)*X10_0(4))) +
  931. & ((2.0D0)*(X10_0(5)*X10_0(5))) +
  932. & ((2.0D0)*(X10_0(6)*X10_0(6))))**(0.5D0))
  933. C
  934. C Check sur |x_i_0| et determination de nn_i_0 de nn2_i_0
  935. c pour l'iteration k=0
  936. IF ((NX1_0).LT.(1.0E-12) ) THEN
  937. NN1_0(1:6) = 0.0D0
  938. NN21_0(1:6) = 0.0D0
  939. ELSE
  940. NN1_0(1:6) = X1_0(1:6)/(NX1_0)
  941. NN21_0(1) = ((NN1_0(1))*(NN1_0(1))) +
  942. & ((NN1_0(4))*(NN1_0(4))) +
  943. & ((NN1_0(5))*(NN1_0(5)))
  944. NN21_0(2) = ((NN1_0(4))*(NN1_0(4))) +
  945. & ((NN1_0(2))*(NN1_0(2))) +
  946. & ((NN1_0(6))*(NN1_0(6)))
  947. NN21_0(3) = ((NN1_0(5))*(NN1_0(5))) +
  948. & ((NN1_0(6))*(NN1_0(6))) +
  949. & ((NN1_0(3))*(NN1_0(3)))
  950. NN21_0(4) = ((NN1_0(1))*(NN1_0(4))) +
  951. & ((NN1_0(4))*(NN1_0(2))) +
  952. & ((NN1_0(5))*(NN1_0(6)))
  953. NN21_0(5) = ((NN1_0(1))*(NN1_0(5))) +
  954. & ((NN1_0(4))*(NN1_0(6))) +
  955. & ((NN1_0(5))*(NN1_0(3)))
  956. NN21_0(6) = ((NN1_0(4))*(NN1_0(5))) +
  957. & ((NN1_0(2))*(NN1_0(6))) +
  958. & ((NN1_0(6))*(NN1_0(3)))
  959. ENDIF
  960. IF ((NX2_0).LT.(1.0E-12) ) THEN
  961. NN2_0(1:6) = 0.0D0
  962. NN22_0(1:6) = 0.0D0
  963. ELSE
  964. NN2_0(1:6) = X2_0(1:6)/(NX2_0)
  965. NN22_0(1) = ((NN2_0(1))*(NN2_0(1))) +
  966. & ((NN2_0(4))*(NN2_0(4))) +
  967. & ((NN2_0(5))*(NN2_0(5)))
  968. NN22_0(2) = ((NN2_0(4))*(NN2_0(4))) +
  969. & ((NN2_0(2))*(NN2_0(2))) +
  970. & ((NN2_0(6))*(NN2_0(6)))
  971. NN22_0(3) = ((NN2_0(5))*(NN2_0(5))) +
  972. & ((NN2_0(6))*(NN2_0(6))) +
  973. & ((NN2_0(3))*(NN2_0(3)))
  974. NN22_0(4) = ((NN2_0(1))*(NN2_0(4))) +
  975. & ((NN2_0(4))*(NN2_0(2))) +
  976. & ((NN2_0(5))*(NN2_0(6)))
  977. NN22_0(5) = ((NN2_0(1))*(NN2_0(5))) +
  978. & ((NN2_0(4))*(NN2_0(6))) +
  979. & ((NN2_0(5))*(NN2_0(3)))
  980. NN22_0(6) = ((NN2_0(4))*(NN2_0(5))) +
  981. & ((NN2_0(2))*(NN2_0(6))) +
  982. & ((NN2_0(6))*(NN2_0(3)))
  983. ENDIF
  984. IF ((NX3_0).LT.(1.0E-12) ) THEN
  985. NN3_0(1:6) = 0.0D0
  986. NN23_0(1:6) = 0.0D0
  987. ELSE
  988. NN3_0(1:6) = X3_0(1:6)/(NX3_0)
  989. NN23_0(1) = ((NN3_0(1))*(NN3_0(1))) +
  990. & ((NN3_0(4))*(NN3_0(4))) +
  991. & ((NN3_0(5))*(NN3_0(5)))
  992. NN23_0(2) = ((NN3_0(4))*(NN3_0(4))) +
  993. & ((NN3_0(2))*(NN3_0(2))) +
  994. & ((NN3_0(6))*(NN3_0(6)))
  995. NN23_0(3) = ((NN3_0(5))*(NN3_0(5))) +
  996. & ((NN3_0(6))*(NN3_0(6))) +
  997. & ((NN3_0(3))*(NN3_0(3)))
  998. NN23_0(4) = ((NN3_0(1))*(NN3_0(4))) +
  999. & ((NN3_0(4))*(NN3_0(2))) +
  1000. & ((NN3_0(5))*(NN3_0(6)))
  1001. NN23_0(5) = ((NN3_0(1))*(NN3_0(5))) +
  1002. & ((NN3_0(4))*(NN3_0(6))) +
  1003. & ((NN3_0(5))*(NN3_0(3)))
  1004. NN23_0(6) = ((NN3_0(4))*(NN3_0(5))) +
  1005. & ((NN3_0(2))*(NN3_0(6))) +
  1006. & ((NN3_0(6))*(NN3_0(3)))
  1007. ENDIF
  1008. IF ((NX4_0).LT.(1.0E-12) ) THEN
  1009. NN4_0(1:6) = 0.0D0
  1010. NN24_0(1:6) = 0.0D0
  1011. ELSE
  1012. NN4_0(1:6) = X4_0(1:6)/(NX4_0)
  1013. NN24_0(1) = ((NN4_0(1))*(NN4_0(1))) +
  1014. & ((NN4_0(4))*(NN4_0(4))) +
  1015. & ((NN4_0(5))*(NN4_0(5)))
  1016. NN24_0(2) = ((NN4_0(4))*(NN4_0(4))) +
  1017. & ((NN4_0(2))*(NN4_0(2))) +
  1018. & ((NN4_0(6))*(NN4_0(6)))
  1019. NN24_0(3) = ((NN4_0(5))*(NN4_0(5))) +
  1020. & ((NN4_0(6))*(NN4_0(6))) +
  1021. & ((NN4_0(3))*(NN4_0(3)))
  1022. NN24_0(4) = ((NN4_0(1))*(NN4_0(4))) +
  1023. & ((NN4_0(4))*(NN4_0(2))) +
  1024. & ((NN4_0(5))*(NN4_0(6)))
  1025. NN24_0(5) = ((NN4_0(1))*(NN4_0(5))) +
  1026. & ((NN4_0(4))*(NN4_0(6))) +
  1027. & ((NN4_0(5))*(NN4_0(3)))
  1028. NN24_0(6) = ((NN4_0(4))*(NN4_0(5))) +
  1029. & ((NN4_0(2))*(NN4_0(6))) +
  1030. & ((NN4_0(6))*(NN4_0(3)))
  1031. ENDIF
  1032. IF ((NX5_0).LT.(1.0E-12) ) THEN
  1033. NN5_0(1:6) = 0.0D0
  1034. NN25_0(1:6) = 0.0D0
  1035. ELSE
  1036. NN5_0(1:6) = X5_0(1:6)/(NX5_0)
  1037. NN25_0(1) = ((NN5_0(1))*(NN5_0(1))) +
  1038. & ((NN5_0(4))*(NN5_0(4))) +
  1039. & ((NN5_0(5))*(NN5_0(5)))
  1040. NN25_0(2) = ((NN5_0(4))*(NN5_0(4))) +
  1041. & ((NN5_0(2))*(NN5_0(2))) +
  1042. & ((NN5_0(6))*(NN5_0(6)))
  1043. NN25_0(3) = ((NN5_0(5))*(NN5_0(5))) +
  1044. & ((NN5_0(6))*(NN5_0(6))) +
  1045. & ((NN5_0(3))*(NN5_0(3)))
  1046. NN25_0(4) = ((NN5_0(1))*(NN5_0(4))) +
  1047. & ((NN5_0(4))*(NN5_0(2))) +
  1048. & ((NN5_0(5))*(NN5_0(6)))
  1049. NN25_0(5) = ((NN5_0(1))*(NN5_0(5))) +
  1050. & ((NN5_0(4))*(NN5_0(6))) +
  1051. & ((NN5_0(5))*(NN5_0(3)))
  1052. NN25_0(6) = ((NN5_0(4))*(NN5_0(5))) +
  1053. & ((NN5_0(2))*(NN5_0(6))) +
  1054. & ((NN5_0(6))*(NN5_0(3)))
  1055. ENDIF
  1056. IF ((NX6_0).LT.(1.0E-12) ) THEN
  1057. NN6_0(1:6) = 0.0D0
  1058. NN26_0(1:6) = 0.0D0
  1059. ELSE
  1060. NN6_0(1:6) = X6_0(1:6)/(NX6_0)
  1061. NN26_0(1) = ((NN6_0(1))*(NN6_0(1))) +
  1062. & ((NN6_0(4))*(NN6_0(4))) +
  1063. & ((NN6_0(5))*(NN6_0(5)))
  1064. NN26_0(2) = ((NN6_0(4))*(NN6_0(4))) +
  1065. & ((NN6_0(2))*(NN6_0(2))) +
  1066. & ((NN6_0(6))*(NN6_0(6)))
  1067. NN26_0(3) = ((NN6_0(5))*(NN6_0(5))) +
  1068. & ((NN6_0(6))*(NN6_0(6))) +
  1069. & ((NN6_0(3))*(NN6_0(3)))
  1070. NN26_0(4) = ((NN6_0(1))*(NN6_0(4))) +
  1071. & ((NN6_0(4))*(NN6_0(2))) +
  1072. & ((NN6_0(5))*(NN6_0(6)))
  1073. NN26_0(5) = ((NN6_0(1))*(NN6_0(5))) +
  1074. & ((NN6_0(4))*(NN6_0(6))) +
  1075. & ((NN6_0(5))*(NN6_0(3)))
  1076. NN26_0(6) = ((NN6_0(4))*(NN6_0(5))) +
  1077. & ((NN6_0(2))*(NN6_0(6))) +
  1078. & ((NN6_0(6))*(NN6_0(3)))
  1079. ENDIF
  1080. IF ((NX7_0).LT.(1.0E-12) ) THEN
  1081. NN7_0(1:6) = 0.0D0
  1082. NN27_0(1:6) = 0.0D0
  1083. ELSE
  1084. NN7_0(1:6) = X7_0(1:6)/(NX7_0)
  1085. NN27_0(1) = ((NN7_0(1))*(NN7_0(1))) +
  1086. & ((NN7_0(4))*(NN7_0(4))) +
  1087. & ((NN7_0(5))*(NN7_0(5)))
  1088. NN27_0(2) = ((NN7_0(4))*(NN7_0(4))) +
  1089. & ((NN7_0(2))*(NN7_0(2))) +
  1090. & ((NN7_0(6))*(NN7_0(6)))
  1091. NN27_0(3) = ((NN7_0(5))*(NN7_0(5))) +
  1092. & ((NN7_0(6))*(NN7_0(6))) +
  1093. & ((NN7_0(3))*(NN7_0(3)))
  1094. NN27_0(4) = ((NN7_0(1))*(NN7_0(4))) +
  1095. & ((NN7_0(4))*(NN7_0(2))) +
  1096. & ((NN7_0(5))*(NN7_0(6)))
  1097. NN27_0(5) = ((NN7_0(1))*(NN7_0(5))) +
  1098. & ((NN7_0(4))*(NN7_0(6))) +
  1099. & ((NN7_0(5))*(NN7_0(3)))
  1100. NN27_0(6) = ((NN7_0(4))*(NN7_0(5))) +
  1101. & ((NN7_0(2))*(NN7_0(6))) +
  1102. & ((NN7_0(6))*(NN7_0(3)))
  1103. ENDIF
  1104. IF ((NX8_0).LT.(1.0E-12) ) THEN
  1105. NN8_0(1:6) = 0.0D0
  1106. NN28_0(1:6) = 0.0D0
  1107. ELSE
  1108. NN8_0(1:6) = X8_0(1:6)/(NX8_0)
  1109. NN28_0(1) = ((NN8_0(1))*(NN8_0(1))) +
  1110. & ((NN8_0(4))*(NN8_0(4))) +
  1111. & ((NN8_0(5))*(NN8_0(5)))
  1112. NN28_0(2) = ((NN8_0(4))*(NN8_0(4))) +
  1113. & ((NN8_0(2))*(NN8_0(2))) +
  1114. & ((NN8_0(6))*(NN8_0(6)))
  1115. NN28_0(3) = ((NN8_0(5))*(NN8_0(5))) +
  1116. & ((NN8_0(6))*(NN8_0(6))) +
  1117. & ((NN8_0(3))*(NN8_0(3)))
  1118. NN28_0(4) = ((NN8_0(1))*(NN8_0(4))) +
  1119. & ((NN8_0(4))*(NN8_0(2))) +
  1120. & ((NN8_0(5))*(NN8_0(6)))
  1121. NN28_0(5) = ((NN8_0(1))*(NN8_0(5))) +
  1122. & ((NN8_0(4))*(NN8_0(6))) +
  1123. & ((NN8_0(5))*(NN8_0(3)))
  1124. NN28_0(6) = ((NN8_0(4))*(NN8_0(5))) +
  1125. & ((NN8_0(2))*(NN8_0(6))) +
  1126. & ((NN8_0(6))*(NN8_0(3)))
  1127. ENDIF
  1128. IF ((NX9_0).LT.(1.0E-12) ) THEN
  1129. NN9_0(1:6) = 0.0D0
  1130. NN29_0(1:6) = 0.0D0
  1131. ELSE
  1132. NN9_0(1:6) = X9_0(1:6)/(NX9_0)
  1133. NN29_0(1) = ((NN9_0(1))*(NN9_0(1))) +
  1134. & ((NN9_0(4))*(NN9_0(4))) +
  1135. & ((NN9_0(5))*(NN9_0(5)))
  1136. NN29_0(2) = ((NN9_0(4))*(NN9_0(4))) +
  1137. & ((NN9_0(2))*(NN9_0(2))) +
  1138. & ((NN9_0(6))*(NN9_0(6)))
  1139. NN29_0(3) = ((NN9_0(5))*(NN9_0(5))) +
  1140. & ((NN9_0(6))*(NN9_0(6))) +
  1141. & ((NN9_0(3))*(NN9_0(3)))
  1142. NN29_0(4) = ((NN9_0(1))*(NN9_0(4))) +
  1143. & ((NN9_0(4))*(NN9_0(2))) +
  1144. & ((NN9_0(5))*(NN9_0(6)))
  1145. NN29_0(5) = ((NN9_0(1))*(NN9_0(5))) +
  1146. & ((NN9_0(4))*(NN9_0(6))) +
  1147. & ((NN9_0(5))*(NN9_0(3)))
  1148. NN29_0(6) = ((NN9_0(4))*(NN9_0(5))) +
  1149. & ((NN9_0(2))*(NN9_0(6))) +
  1150. & ((NN9_0(6))*(NN9_0(3)))
  1151. ENDIF
  1152. IF ((NX10_0).LT.(1.0E-12) ) THEN
  1153. NN10_0(1:6) = 0.0D0
  1154. NN210_0(1:6) = 0.0D0
  1155. ELSE
  1156. NN10_0(1:6) = X10_0(1:6)/(NX10_0)
  1157. NN210_0(1) = ((NN10_0(1))*(NN10_0(1))) +
  1158. & ((NN10_0(4))*(NN10_0(4))) +
  1159. & ((NN10_0(5))*(NN10_0(5)))
  1160. NN210_0(2) = ((NN10_0(4))*(NN10_0(4))) +
  1161. & ((NN10_0(2))*(NN10_0(2))) +
  1162. & ((NN10_0(6))*(NN10_0(6)))
  1163. NN210_0(3) = ((NN10_0(5))*(NN10_0(5))) +
  1164. & ((NN10_0(6))*(NN10_0(6))) +
  1165. & ((NN10_0(3))*(NN10_0(3)))
  1166. NN210_0(4) = ((NN10_0(1))*(NN10_0(4))) +
  1167. & ((NN10_0(4))*(NN10_0(2))) +
  1168. & ((NN10_0(5))*(NN10_0(6)))
  1169. NN210_0(5) = ((NN10_0(1))*(NN10_0(5))) +
  1170. & ((NN10_0(4))*(NN10_0(6))) +
  1171. & ((NN10_0(5))*(NN10_0(3)))
  1172. NN210_0(6) = ((NN10_0(4))*(NN10_0(5))) +
  1173. & ((NN10_0(2))*(NN10_0(6))) +
  1174. & ((NN10_0(6))*(NN10_0(3)))
  1175. ENDIF
  1176. C
  1177. C Calcul de cos(3th_i_0)=-sqrt(6)*trace(nn_i_0^3) -
  1178. c depandance angle de Lode pour k=0
  1179. COS3TH1_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1180. & (((NN1_0(1))*(NN21_0(1))) +
  1181. & ((NN1_0(2))*(NN21_0(2))) +
  1182. & ((NN1_0(3))*(NN21_0(3))) +
  1183. & (2.0D0)*((NN1_0(4))*(NN21_0(4))) +
  1184. & (2.0D0)*((NN1_0(5))*(NN21_0(5))) +
  1185. & (2.0D0)*((NN1_0(6))*(NN21_0(6))))
  1186. COS3TH2_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1187. & (((NN2_0(1))*(NN22_0(1))) +
  1188. & ((NN2_0(2))*(NN22_0(2))) +
  1189. & ((NN2_0(3))*(NN22_0(3))) +
  1190. & (2.0D0)*((NN2_0(4))*(NN22_0(4))) +
  1191. & (2.0D0)*((NN2_0(5))*(NN22_0(5))) +
  1192. & (2.0D0)*((NN2_0(6))*(NN22_0(6))))
  1193. COS3TH3_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1194. & (((NN3_0(1))*(NN23_0(1))) +
  1195. & ((NN3_0(2))*(NN23_0(2))) +
  1196. & ((NN3_0(3))*(NN23_0(3))) +
  1197. & (2.0D0)*((NN3_0(4))*(NN23_0(4))) +
  1198. & (2.0D0)*((NN3_0(5))*(NN23_0(5))) +
  1199. & (2.0D0)*((NN3_0(6))*(NN23_0(6))))
  1200. COS3TH4_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1201. & (((NN4_0(1))*(NN24_0(1))) +
  1202. & ((NN4_0(2))*(NN24_0(2))) +
  1203. & ((NN4_0(3))*(NN24_0(3))) +
  1204. & (2.0D0)*((NN4_0(4))*(NN24_0(4))) +
  1205. & (2.0D0)*((NN4_0(5))*(NN24_0(5))) +
  1206. & (2.0D0)*((NN4_0(6))*(NN24_0(6))))
  1207. COS3TH5_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1208. & (((NN5_0(1))*(NN25_0(1))) +
  1209. & ((NN5_0(2))*(NN25_0(2))) +
  1210. & ((NN5_0(3))*(NN25_0(3))) +
  1211. & (2.0D0)*((NN5_0(4))*(NN25_0(4))) +
  1212. & (2.0D0)*((NN5_0(5))*(NN25_0(5))) +
  1213. & (2.0D0)*((NN5_0(6))*(NN25_0(6))))
  1214. COS3TH6_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1215. & (((NN6_0(1))*(NN26_0(1))) +
  1216. & ((NN6_0(2))*(NN26_0(2))) +
  1217. & ((NN6_0(3))*(NN26_0(3))) +
  1218. & (2.0D0)*((NN6_0(4))*(NN26_0(4))) +
  1219. & (2.0D0)*((NN6_0(5))*(NN26_0(5))) +
  1220. & (2.0D0)*((NN6_0(6))*(NN26_0(6))))
  1221. COS3TH7_0 = ((-1.0D0)*((6.0D0)**(0.5)))*
  1222. & (((NN7_0(1))*(NN27_0(1))) +
  1223. & ((NN7_0(2))*(NN27_0(2))) +
  1224. & ((NN7_0(3))*(NN27_0(3))) +
  1225. & (2.0D0)*((NN7_0(4))*(NN27_0(4))) +
  1226. & (2.0D0)*((NN7_0(5))*(NN27_0(5))) +
  1227. & (2.0D0)*((NN7_0(6))*(NN27_0(6))))
  1228. COS3TH8_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1229. & (((NN8_0(1))*(NN28_0(1))) +
  1230. & ((NN8_0(2))*(NN28_0(2))) +
  1231. & ((NN8_0(3))*(NN28_0(3))) +
  1232. & (2.0D0)*((NN8_0(4))*(NN28_0(4))) +
  1233. & (2.0D0)*((NN8_0(5))*(NN28_0(5))) +
  1234. & (2.0D0)*((NN8_0(6))*(NN28_0(6))))
  1235. COS3TH9_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1236. & (((NN9_0(1))*(NN29_0(1))) +
  1237. & ((NN9_0(2))*(NN29_0(2))) +
  1238. & ((NN9_0(3))*(NN29_0(3))) +
  1239. & (2.0D0)*((NN9_0(4))*(NN29_0(4))) +
  1240. & (2.0D0)*((NN9_0(5))*(NN29_0(5))) +
  1241. & (2.0D0)*((NN9_0(6))*(NN29_0(6))))
  1242. COS3TH10_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1243. & (((NN10_0(1))*(NN210_0(1))) +
  1244. & ((NN10_0(2))*(NN210_0(2))) +
  1245. & ((NN10_0(3))*(NN210_0(3))) +
  1246. & (2.0D0)*((NN10_0(4))*(NN210_0(4))) +
  1247. & (2.0D0)*((NN10_0(5))*(NN210_0(5))) +
  1248. & (2.0D0)*((NN10_0(6))*(NN210_0(6))))
  1249. C
  1250. C Calcul du terme
  1251. C R(th)_i_0=2k1/((1+k1)-(1-k1)*(cos(3th_i_0))) pour
  1252. C l'iteration k=0
  1253. RTH1_0 = ((2.0D0)*(XK1))/
  1254. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH1_0)))
  1255. RTH2_0 = ((2.0D0)*(XK1))/
  1256. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH2_0)))
  1257. RTH3_0 = ((2.0D0)*(XK1))/
  1258. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH3_0)))
  1259. RTH4_0 = ((2.0D0)*(XK1))/
  1260. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH4_0)))
  1261. RTH5_0 = ((2.0D0)*(XK1))/
  1262. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH5_0)))
  1263. RTH6_0 = ((2.0D0)*(XK1))/
  1264. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH6_0)))
  1265. RTH7_0 = ((2.0D0)*(XK1))/
  1266. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH7_0)))
  1267. RTH8_0 = ((2.0D0)*(XK1))/
  1268. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH8_0)))
  1269. RTH9_0 = ((2.0D0)*(XK1))/
  1270. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH9_0)))
  1271. RTH10_0 = ((2.0D0)*(XK1))/
  1272. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH10_0)))
  1273. C
  1274. C Calcul des criteres f_i_0 pour l'iteration k=0
  1275. F1_0 = NX1_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM1))*
  1276. & ((P_0 - XC1)*(RTH1_0))
  1277. F2_0 = NX2_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM2))*
  1278. & ((P_0 - XC1)*(RTH2_0))
  1279. F3_0 = NX3_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM3))*
  1280. & ((P_0 - XC1)*(RTH3_0))
  1281. F4_0 = NX4_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM4))*
  1282. & ((P_0 - XC1)*(RTH4_0))
  1283. F5_0 = NX5_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM5))*
  1284. & ((P_0 - XC1)*(RTH5_0))
  1285. F6_0 = NX6_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM6))*
  1286. & ((P_0 - XC1)*(RTH6_0))
  1287. F7_0 = NX7_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM7))*
  1288. & ((P_0 - XC1)*(RTH7_0))
  1289. F8_0 = NX8_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM8))*
  1290. & ((P_0 - XC1)*(RTH8_0))
  1291. F9_0 = NX9_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM9))*
  1292. & ((P_0 - XC1)*(RTH9_0))
  1293. F10_0 = NX10_0 + ((((2.0D0)/(3.0D0))**(0.5D0))*(XM10))*
  1294. & ((P_0 - XC1)*(RTH10_0))
  1295.  
  1296. C
  1297. C Detectetion du critere activé pour l'iteration k=0
  1298. VET_TEST(1) = F1_0
  1299. VET_TEST(2) = F2_0
  1300. VET_TEST(3) = F3_0
  1301. VET_TEST(4) = F4_0
  1302. VET_TEST(5) = F5_0
  1303. VET_TEST(6) = F6_0
  1304. VET_TEST(7) = F7_0
  1305. VET_TEST(8) = F8_0
  1306. VET_TEST(9) = F9_0
  1307. VET_TEST(10) = F10_0
  1308. DO 102 II1 = 1,10
  1309. IF ((VET_TEST(II1).GT.(0.0D0)).OR.
  1310. & ((ABS(VET_TEST(II1))).LT.(1.0E-8))) THEN
  1311. VET_TEST2(II1) = INT(1)
  1312. ELSE
  1313. VET_TEST2(II1) = INT(0)
  1314. ENDIF
  1315. 102 CONTINUE
  1316. I1_CRIT = INT(0)
  1317. DO 103 II1 = 1,10
  1318. I1_CRIT = I1_CRIT + (INT(VET_TEST2(II1)))
  1319. 103 CONTINUE
  1320. I2_CRIT = I1_CRIT + INT(1);
  1321. C
  1322. C Choix des parametres selon le critere activé pour la
  1323. C recherche de l'etat plastiquement admissible
  1324. C XM m (inclinaison de la pente du critere activé)
  1325. C XHC Hc (module plastique en compression pour crit activ)
  1326. C XHE He (module plastique en exstension pour crit activ)
  1327. C XETC eta_c (pente surf caracteristique en compress)
  1328. C XETE eta_e (pente surf caracteristique en exstension)
  1329. C XD D (parametre de controle de la dilatance)
  1330. IF (I1_CRIT.EQ.(INT(1))) THEN
  1331. BET_0(1:6) = BET1_0(1:6)
  1332. XM = XM1
  1333. XHC = XHC1
  1334. XHE = XHE1
  1335. XETC = XETC1
  1336. XETE = XETE1
  1337. XD = XD1
  1338. ELSEIF (I1_CRIT.EQ.(INT(2))) THEN
  1339. BET_0(1:6) = BET2_0(1:6)
  1340. XM = XM2
  1341. XHC = XHC2
  1342. XHE = XHE2
  1343. XETC = XETC2
  1344. XETE = XETE2
  1345. XD = XD2
  1346. ELSEIF (I1_CRIT.EQ.(INT(3))) THEN
  1347. BET_0(1:6) = BET3_0(1:6)
  1348. XM = XM3
  1349. XHC = XHC3
  1350. XHE = XHE3
  1351. XETC = XETC3
  1352. XETE = XETE3
  1353. XD = XD3
  1354. ELSEIF (I1_CRIT.EQ.(INT(4))) THEN
  1355. BET_0(1:6) = BET4_0(1:6)
  1356. XM = XM4
  1357. XHC = XHC4
  1358. XHE = XHE4
  1359. XETC = XETC4
  1360. XETE = XETE4
  1361. XD = XD4
  1362. ELSEIF (I1_CRIT.EQ.(INT(5))) THEN
  1363. BET_0(1:6) = BET5_0(1:6)
  1364. XM = XM5
  1365. XHC = XHC5
  1366. XHE = XHE5
  1367. XETC = XETC5
  1368. XETE = XETE5
  1369. XD = XD5
  1370. ELSEIF (I1_CRIT.EQ.(INT(6))) THEN
  1371. BET_0(1:6) = BET6_0(1:6)
  1372. XM = XM6
  1373. XHC = XHC6
  1374. XHE = XHE6
  1375. XETC = XETC6
  1376. XETE = XETE6
  1377. XD = XD6
  1378. ELSEIF (I1_CRIT.EQ.(INT(7))) THEN
  1379. BET_0(1:6) = BET7_0(1:6)
  1380. XM = XM7
  1381. XHC = XHC7
  1382. XHE = XHE7
  1383. XETC = XETC7
  1384. XETE = XETE7
  1385. XD = XD7
  1386. ELSEIF (I1_CRIT.EQ.(INT(8))) THEN
  1387. BET_0(1:6) = BET8_0(1:6)
  1388. XM = XM8
  1389. XHC = XHC8
  1390. XHE = XHE8
  1391. XETC = XETC8
  1392. XETE = XETE8
  1393. XD = XD8
  1394. ELSEIF (I1_CRIT.EQ.(INT(9))) THEN
  1395. BET_0(1:6) = BET9_0(1:6)
  1396. XM = XM9
  1397. XHC = XHC9
  1398. XHE = XHE9
  1399. XETC = XETC9
  1400. XETE = XETE9
  1401. XD = XD9
  1402. ELSEIF (I1_CRIT.EQ.(INT(10))) THEN
  1403. BET_0(1:6) = BET10_0(1:6)
  1404. XM = XM10
  1405. XHC = XHC10
  1406. XHE = XHE10
  1407. XETC = XETC10
  1408. XETE = XETE10
  1409. XD = XD10
  1410. ELSEIF (I1_CRIT.EQ.(INT(0))) THEN
  1411. BET_0(1:6) = BET1_0(1:6)
  1412. XM = XM1
  1413. XHC = XHC1
  1414. XHE = XHE1
  1415. XETC = XETC1
  1416. XETE = XETE1
  1417. XD = XD1
  1418. ENDIF
  1419. ENDIF
  1420. C
  1421. C Calcul des residues des contraintes et des derivé pour
  1422. C l'instant n+1 à l'iteration k
  1423. C Moidification V4 22-03-2016
  1424. IF ((ABS(XN1 - 0.0D0)).LT.(1.E-12)) THEN
  1425. DGDP_0 = 0.0D0;
  1426. DKDP_0 = 0.0D0;
  1427. ELSE
  1428. DGDP_0 = ((GMOD)*(((P_0)/(XPREF1))**(XN1 - 1.0D0)))*
  1429. & ((XN1)/(XPREF1))
  1430. DKDP_0 = ((KMOD)*(((P_0)/(XPREF1))**(XN1 - 1.0D0)))*
  1431. & ((XN1)/(XPREF1))
  1432. ENDIF
  1433. GMOD_0 = ((GMOD)*(((P_0)/(XPREF1))**(XN1)))
  1434. KMOD_0 = ((KMOD)*(((P_0)/(XPREF1))**(XN1)))
  1435. c Determination residu pour l iteration k
  1436. R_0(1:6) = ((-1.0D0)*(S_0)) + S_INI +
  1437. & ((2.0D0)*(GMOD_0))*(DSE_TOT(1:6) -
  1438. & (SEPLA_0(1:6) - SEPLA_INI(1:6)))
  1439. IF (IFOUR.NE.2) THEN
  1440. R_0(5:6) = 0.0D0
  1441. ENDIF
  1442.  
  1443. R_0(7) = -P_0 + P_INI +
  1444. & ((1.0D0)*(KMOD_0))*(DVE_TOT -
  1445. & (VEPLA_0 - VEPLA_INI))
  1446. C Determination du terme xi pour l'instant n+1 à l'iteration k
  1447. X_0(1:6) = S_0(1:6) + ((P_0 - XC1)*(BET_0(1:6)))
  1448. C Determination du terme |xi| pour l'instant n+1 à l'iteration k
  1449. NX_0 = (((X_0(1)*X_0(1)) +
  1450. & (X_0(2)*X_0(2)) +
  1451. & (X_0(3)*X_0(3)) +
  1452. & ((2.0D0)*(X_0(4)*X_0(4))) +
  1453. & ((2.0D0)*(X_0(5)*X_0(5))) +
  1454. & ((2.0D0)*(X_0(6)*X_0(6))))**(0.5D0))
  1455. C
  1456. C Determination du terme nn=xi/|xi| et nn2 = nn*nn pour
  1457. C l'instant n+1 à l'iteration k
  1458. IF ((NX_0).LT.(1.0E-12) ) THEN
  1459. NN_0(1:6) = 0.0D0
  1460. N2N_0(1:6) = 0.0D0
  1461. ELSE
  1462. NN_0(1:6) = X_0(1:6)/(NX_0)
  1463. N2N_0(1) = ((NN_0(1))*(NN_0(1))) +
  1464. & ((NN_0(4))*(NN_0(4))) +
  1465. & ((NN_0(5))*(NN_0(5)))
  1466. N2N_0(2) = ((NN_0(4))*(NN_0(4))) +
  1467. & ((NN_0(2))*(NN_0(2))) +
  1468. & ((NN_0(6))*(NN_0(6)))
  1469. N2N_0(3) = ((NN_0(5))*(NN_0(5))) +
  1470. & ((NN_0(6))*(NN_0(6))) +
  1471. & ((NN_0(3))*(NN_0(3)))
  1472. N2N_0(4) = ((NN_0(1))*(NN_0(4))) +
  1473. & ((NN_0(4))*(NN_0(2))) +
  1474. & ((NN_0(5))*(NN_0(6)))
  1475. N2N_0(5) = ((NN_0(1))*(NN_0(5))) +
  1476. & ((NN_0(4))*(NN_0(6))) +
  1477. & ((NN_0(5))*(NN_0(3)))
  1478. N2N_0(6) = ((NN_0(4))*(NN_0(5))) +
  1479. & ((NN_0(2))*(NN_0(6))) +
  1480. & ((NN_0(6))*(NN_0(3)))
  1481. ENDIF
  1482. C
  1483. C Determination du terme cos(3th) à l'instnt n+1 pour
  1484. C l'iteration k
  1485. COS3TH_0 = ((-1.0D0)*((6.0D0)**(0.5D0)))*
  1486. & (((NN_0(1))*(N2N_0(1))) +
  1487. & ((NN_0(2))*(N2N_0(2))) +
  1488. & ((NN_0(3))*(N2N_0(3))) +
  1489. & (2.0D0)*((NN_0(4))*(N2N_0(4))) +
  1490. & (2.0D0)*((NN_0(5))*(N2N_0(5))) +
  1491. & (2.0D0)*((NN_0(6))*(N2N_0(6))))
  1492.  
  1493. C
  1494. C Determination du terme R(th) à l'instant n+1 pour
  1495. C l'iteration k
  1496. RTH_0 = ((2.0D0)*(XK1))/
  1497. & ((1.0D0 + XK1) - ((1.0D0 - XK1)*(COS3TH_0)))
  1498.  
  1499. C
  1500. C Determination du terme eta_0= (((3/2)*(s_0*s_0))^0.5)/p_0
  1501. C à l'instant n+1 pour l'iteration k
  1502. ETA_0 = (((((3.0D0)/(2.0D0))*(((S_0(1))*(S_0(1))) +
  1503. & ((S_0(2))*(S_0(2))) +
  1504. & ((S_0(3))*(S_0(3))) +
  1505. & (2.0D0)*((S_0(4))*(S_0(4))) +
  1506. & (2.0D0)*((S_0(5))*(S_0(5))) +
  1507. & (2.0D0)*((S_0(6))*(S_0(6)))))**(0.5D0)))
  1508. & /(P_0)
  1509. C
  1510. C Determination du terme de la loi d'ecoulement
  1511. C P1_0 = D*(((eta/eta*)² - 1)/((eta/eta*)² + 1)) où
  1512. C eta* = eta_c si cos3th=tr(n3)>0
  1513. C eta* = eta_e si cos3th=tr(n3)<0
  1514. C pour l'instant n+ 1 à l'iteration k
  1515. C
  1516. C Modification 22-03-2016 valeur d'exclusions
  1517. IF ((COS3TH_0.GE.(0.0D0)).OR.
  1518. & ((ABS(XK1 - 1.0D0)).LT.(1.0E-12))) THEN
  1519. P1_0 = (XD)*(((((ETA_0)/(XETC))**(2.0D0)) - 1.0D0)
  1520. & /((((ETA_0)/(XETC))**(2.0D0)) + 1.0D0))
  1521. IF (ETA_0.GT.(1.0E8)) THEN
  1522. P1_0 = (XD)*(1.0D0)
  1523. ENDIF
  1524. ELSE
  1525. P1_0 = (XD)*(((((ETA_0)/(XETE))**(2.0D0)) - 1.0D0)
  1526. & /((((ETA_0)/(XETE))**(2.0D0)) + 1.0D0))
  1527. IF (ETA_0.GT.(1.0E8)) THEN
  1528. P1_0 = (XD)*(1.0D0)
  1529. ENDIF
  1530. ENDIF
  1531. IF ((ABS(XD - 0.0D0)).LT.(1.0E-12)) THEN
  1532. P1_0 = 0.0D0
  1533. ENDIF
  1534. c
  1535. C
  1536. C Determination du module plastique en fonction
  1537. C de l'angle de lode de la pression hydrostatique:
  1538. C H = H'*(p_0/p_ref)^n1
  1539. C H' = (Hc - He)/2*cos(3th) + (Hc + He)/2
  1540. C pour l'instant n+ 1 à l'iteration k
  1541. IF ((ABS(XK1 - 1.0D0)).LT.(1.0E-12)) THEN
  1542. HMOD_0 = XHC
  1543. ELSE
  1544. HMOD_0 = ((((XHC - XHE)/(2.0D0))*(COS3TH_0))
  1545. & + ((XHC + XHE)/(2.0D0)))*(((P_0)/(XPREF1))**(XN1))
  1546. ENDIF
  1547. C
  1548. C Determination du terme D2f pour l'instant n+1 à l'iteration k
  1549. C D2f = 1/|xi|*(Id - nn_0^nn_0)
  1550. D2F(1,1) = ((NN_0(1))*(NN_0(1)))
  1551. D2F(1,2) = ((NN_0(1))*(NN_0(2)))
  1552. D2F(1,3) = ((NN_0(1))*(NN_0(3)))
  1553. D2F(1,4) = (2.0D0)*((NN_0(1))*(NN_0(4)))
  1554. D2F(1,5) = (2.0D0)*((NN_0(1))*(NN_0(5)))
  1555. D2F(1,6) = (2.0D0)*((NN_0(1))*(NN_0(6)))
  1556. C
  1557. D2F(2,1) = ((NN_0(2))*(NN_0(1)))
  1558. D2F(2,2) = ((NN_0(2))*(NN_0(2)))
  1559. D2F(2,3) = ((NN_0(2))*(NN_0(3)))
  1560. D2F(2,4) = (2.0D0)*((NN_0(2))*(NN_0(4)))
  1561. D2F(2,5) = (2.0D0)*((NN_0(2))*(NN_0(5)))
  1562. D2F(2,6) = (2.0D0)*((NN_0(2))*(NN_0(6)))
  1563. C
  1564. D2F(3,1) = ((NN_0(3))*(NN_0(1)))
  1565. D2F(3,2) = ((NN_0(3))*(NN_0(2)))
  1566. D2F(3,3) = ((NN_0(3))*(NN_0(3)))
  1567. D2F(3,4) = (2.0D0)*((NN_0(3))*(NN_0(4)))
  1568. D2F(3,5) = (2.0D0)*((NN_0(3))*(NN_0(5)))
  1569. D2F(3,6) = (2.0D0)*((NN_0(3))*(NN_0(6)))
  1570. C
  1571. D2F(4,1) = ((NN_0(4))*(NN_0(1)))
  1572. D2F(4,2) = ((NN_0(4))*(NN_0(2)))
  1573. D2F(4,3) = ((NN_0(4))*(NN_0(3)))
  1574. D2F(4,4) = (2.0D0)*((NN_0(4))*(NN_0(4)))
  1575. D2F(4,5) = (2.0D0)*((NN_0(4))*(NN_0(5)))
  1576. D2F(4,6) = (2.0D0)*((NN_0(4))*(NN_0(6)))
  1577. C
  1578. D2F(5,1) = ((NN_0(5))*(NN_0(1)))
  1579. D2F(5,2) = ((NN_0(5))*(NN_0(2)))
  1580. D2F(5,3) = ((NN_0(5))*(NN_0(3)))
  1581. D2F(5,4) = (2.0D0)*((NN_0(5))*(NN_0(4)))
  1582. D2F(5,5) = (2.0D0)*((NN_0(5))*(NN_0(5)))
  1583. D2F(5,6) = (2.0D0)*((NN_0(5))*(NN_0(6)))
  1584. C
  1585. D2F(6,1) = ((NN_0(6))*(NN_0(1)))
  1586. D2F(6,2) = ((NN_0(6))*(NN_0(2)))
  1587. D2F(6,3) = ((NN_0(6))*(NN_0(3)))
  1588. D2F(6,4) = (2.0D0)*((NN_0(6))*(NN_0(4)))
  1589. D2F(6,5) = (2.0D0)*((NN_0(6))*(NN_0(5)))
  1590. D2F(6,6) = (2.0D0)*((NN_0(6))*(NN_0(6)))
  1591. C
  1592. D2F = ((-1.0D0)*(D2F))
  1593. D2F(1:6,1:6) = (XID2(1:6,1:6) + D2F(1:6,1:6))/(NX_0)
  1594. C
  1595. C Determination du verseur (direction dans la quelle le critere
  1596. C evolue) mm_0 selon la relge d'ecrouissage de Morz
  1597. S1M(1:6) = S_0(1:6)
  1598. IF (I2_CRIT.EQ.(INT(2))) THEN
  1599. S2SCAL = (-1.0D0)*
  1600. & ((((((2.0D0)/(3.0D0))**(0.5D0))*(XM2))*
  1601. & (P_0 - XC1))*(RTH_0))
  1602. S2M(1:6) = ((S2SCAL)*(NN_0)) - ((P_0 - XC1)*(BET2_0(1:6)))
  1603. MM_0(1:6) = (S2M(1:6)) - (S1M(1:6))
  1604. NMM_0 = (((MM_0(1))*(MM_0(1))) +
  1605. & ((MM_0(2))*(MM_0(2))) +
  1606. & ((MM_0(3))*(MM_0(3))) +
  1607. & (2.0D0)*((MM_0(4))*(MM_0(4))) +
  1608. & (2.0D0)*((MM_0(5))*(MM_0(5))) +
  1609. & (2.0D0)*((MM_0(6))*(MM_0(6))))**(0.5D0)
  1610. MM_0(1:6) = ((-1.0D0)*((MM_0(1:6))/(NMM_0)))
  1611. ELSEIF (I2_CRIT.EQ.(INT(3))) THEN
  1612. S2SCAL = (-1.0D0)*
  1613. & ((((((2.0D0)/(3.0D0))**(0.5D0))*(XM3))*
  1614. & (P_0 - XC1))*(RTH_0))
  1615. S2M(1:6) = ((S2SCAL)*(NN_0)) - ((P_0 - XC1)*(BET3_0(1:6)))
  1616. MM_0(1:6) = (S2M(