Télécharger bbcal2.eso

Retour à la liste

Numérotation des lignes :

  1. C BBCAL2 SOURCE KICH 18/01/11 21:15:02 9690
  2. SUBROUTINE BBCAL2(I2B,I2GAU,IDIM,NBPGAU,IVACAR,
  3. 1 POIGAU,QSIGAU,ETAGAU,DZEGAU,MELE,MFR,NBNN,LRE,IFOUR,NSTRS,NHRM,
  4. 2 A,BB,XE,SHPTOT,SHP,BGENE,XDPGE,YDPGE,P)
  5. * SUBROUTINE BBCALC(XE,MELE,NBNO,IDIM,NBPGAU,POIGAU,
  6. * 1 QSIGAU,ETAGAU,DZEGAU,
  7. * 2 NST,LRE,IFOU,A,BB,
  8. * 3 NN,SHPTOT,SHP,XDPGE,YDPGE)
  9. *******************************************************************
  10. * projection de dN/dx dans l espace de dimension inferieure
  11. c=======================================================================
  12. c calcul des composantes de la matrice b-barre-dlatation
  13. c - en 2d: elements ict3,icq4,ict6
  14. c - en 3d: elements ....
  15. c calcul des coefficients de modification de la matrice b-barre-dilatation
  16. c - en 2d: elements icq8
  17. c - en 3d: elements ....
  18. c coefficients {a}={a1,a2,a3}
  19. c=======================================================================
  20. c input
  21. c xe = coordonnees de l'element
  22. c mele = numero de l'element dans nomtp
  23. c nbnn = nombre de noeuds
  24. c idim = dimension espace
  25. c nbpgau = nombre de points de gauss pour la rigidite
  26. c poigau = poids des points de gauss
  27. c qsigau,etagau,dzegau = coordonnees des points de gauss
  28. c nstrs = nombre de composantes de contraintes
  29. c lre = nombre de colonnes de la matrice b
  30. c ifour = ifour dans ccoptio
  31. c nhrm = numero du mode de fourier
  32. c xdpge,ydpge = coordonnees du point autour duquel se fait
  33. c le mouvement de la section
  34. c
  35. c output
  36. c tableau des coefficients de modification: a(4,*)
  37. c tableau des composantes bbarre-dilatation: bb(3,*)
  38. c=======================================================================
  39. IMPLICIT INTEGER(I-N)
  40. IMPLICIT REAL*8 (A-H,O-Z)
  41. -INC CCREEL
  42. -INC SMCHAML
  43.  
  44. PARAMETER (XZer=0., XUn=1., X1s2=0.5, X1s4 = 0.25,
  45. & X2s3 = 0.666666666666666666666666666666666666666667,
  46. & X1s3 = 0.333333333333333333333333333333333333333333,
  47. & X1s6 = 0.166666666666666666666666666666666666666667)
  48.  
  49. DIMENSION XE(3,*), POIGAU(*),QSIGAU(*),ETAGAU(*),DZEGAU(*),
  50. & SHP(6,*),SHPTOT(6,NBNN,*),
  51. & A(4,*),BB(3,*),BGENE(NSTRS,*)
  52.  
  53. DIMENSION P(4,4),XNUM1(8),XNUM2(8),XNUM3(8)
  54.  
  55. SEGMENT MPTVAL
  56. INTEGER IPOS(NS) ,NSOF(NS)
  57. INTEGER IVAL(NCOSOU)
  58. CHARACTER*16 TYVAL(NCOSOU)
  59. ENDSEGMENT
  60.  
  61. if (MELE.EQ.72.or.MELE.GE.76) then
  62. CALL ZERO(A,4,60)
  63. CALL ZERO(P,4,4)
  64. endif
  65.  
  66. IFR = IFOUR+4
  67. if (ifr.le.4.or.ifr.ge.7) then
  68. KINC = 2
  69. else
  70. KINC = 3
  71. endif
  72. * if(i2b.eq.34)
  73. * &write(6,*)'BBCAL2',xzprec,I2B,I2GAU,nbnn,NBPGAU,IFOUR,mele,kinc
  74.  
  75.  
  76. IF ( MELE.LT.69 .OR.
  77. &(MELE.GT.78.AND.MELE.NE.273.AND.MELE.NE.274)) GOTO 999
  78. C=======================================================================
  79. C= Elements incompressibles implementes :
  80. C=======================================================================
  81. C= NOM : ICT3, ICQ4, ICT6, ICQ8, ICC8, ICT4, ICP6, IC20, IC10, IC15,
  82. C= MELE : 69 , 70 , 71 , 72 , 73 , 74 , 75 , 76 , 77 , 78 ,
  83. C= NOM : ICY5, IC13
  84. C= MELE : 273, 274
  85. CC=======================================================================
  86. CALL ZERO(BB,3,60)
  87. IF(MELE.EQ.273.OR.MELE.EQ.274) GOTO 100
  88. GOTO ( 10, 20, 30, 40, 100,100,100,100,100,100) , (MELE-68)
  89. GOTO 999
  90.  
  91. 3 FORMAT (4(A,I1,A,I1,A,F10.4,2X))
  92. 4 FORMAT (3(A,I1,A,I1,A,F10.4,2X))
  93.  
  94. C=======================================================================
  95. C========== Elements MASSIFS INCOMPRESSIBLES BIDIMENSIONNELS ===========
  96. C=======================================================================
  97.  
  98. *-----------------------------------------------------------------------
  99. *----------------------------- Element ICT3 ----------------------------
  100. *-----------------------------------------------------------------------
  101. 10 CONTINUE
  102. * write(*,*) 'Element ICT3',ifour
  103. * write(*,*) 'Ecriture de la matrice xe[ ]'
  104. * write(*,3) (('xe(',i,',',j,')=',xe(i,j),i=1,2),j=1,3)
  105. IFR = ifour+4
  106. GOTO (666,101,101,103,666,999),IFR
  107. GOTO 999
  108.  
  109. *--------contraintes planes ou déformations planes (ifour=-2,-1)
  110. 101 CONTINUE
  111. *--------donnee des composantes de bb-dil
  112. r_z = - XE(1,1)*XE(2,3) - XE(1,2)*XE(2,1) + XE(1,2)*XE(2,3)
  113. & + XE(1,1)*XE(2,2) + XE(1,3)*XE(2,1) - XE(1,3)*XE(2,2)
  114. * AIRE = X1s2 * r_z
  115. r_z = XUn / r_z
  116. BB(1,1) = (XE(2,2)-XE(2,3)) * r_z
  117. BB(2,2) = (XE(1,3)-XE(1,2)) * r_z
  118. BB(1,3) = (XE(2,3)-XE(2,1)) * r_z
  119. BB(2,4) = (XE(1,1)-XE(1,3)) * r_z
  120. BB(1,5) = (XE(2,1)-XE(2,2)) * r_z
  121. BB(2,6) = (XE(1,2)-XE(1,1)) * r_z
  122. * write (*,*) 'calcul de aire pour bbdil'
  123. * write (*,*) 'aire = ',AIRE
  124. * write (*,*) 'ecriture de la matrice bb'
  125. * write (*,3) (('bb(',i,',',j,')=',bb(i,j),i=1,3),j=1,2*nbnn)
  126. GOTO 666
  127. *--------axisymétrique (ifour=0)
  128. 103 CONTINUE
  129. *--------donnee des composantes de bb-dil
  130. r_z = - XE(1,1)*XE(2,3) - XE(1,2)*XE(2,1) + XE(1,2)*XE(2,3)
  131. & + XE(1,1)*XE(2,2) + XE(1,3)*XE(2,1) - XE(1,3)*XE(2,2)
  132. * AIRE = X1s6 * r_z * (XE(1,1)+XE(1,2)+XE(1,3))
  133. r_z = XUn / r_z
  134. BB(1,1) = r_z * (XE(2,2)-XE(2,3))
  135. BB(2,2) = r_z * (XE(1,3)-XE(1,2))
  136. BB(1,3) = r_z * (XE(2,3)-XE(2,1))
  137. BB(2,4) = r_z * (XE(1,1)-XE(1,3))
  138. BB(1,5) = r_z * (XE(2,1)-XE(2,2))
  139. BB(2,6) = r_z * (XE(1,2)-XE(1,1))
  140. *-- Integration analytique des termes BB(3,1),BB(3,3),BB(3,5)
  141. BB(3,1) = XUn / (XE(1,1)+XE(1,2)+XE(1,3))
  142. BB(3,3) = BB(3,1)
  143. BB(3,5) = BB(3,1)
  144. * write(*,*) 'calcul de aire pour bbdil'
  145. * write(*,*) 'aire = ',AIRE
  146. * write(*,*) 'ecriture de la matrice bb'
  147. * write(*,4) (('bb(',i,',',j,')=',bb(i,j),i=1,3),j=1,2*nbnn)
  148. GOTO 666
  149.  
  150. *-----------------------------------------------------------------------
  151. *----------------------------- Element ICQ4 ----------------------------
  152. *-----------------------------------------------------------------------
  153. 20 CONTINUE
  154. * write(*,*) 'Element ICQ4',ifour
  155. * write(*,*) 'Ecriture de la matrice xe[ ]'
  156. * write(*,3) (('xe(',i,',',j,')=',xe(i,j),i=1,2),j=1,4)
  157. IFR = ifour+4
  158. GOTO (666,201,201,203,666,999),IFR
  159. GOTO 999
  160. *--------contraintes planes ou déformations planes (ifour=-2,-1)
  161. 201 CONTINUE
  162. *--------donnée des composantes de bb-dil
  163. r_z = - XE(1,4)*XE(2,3) - XE(1,2)*XE(2,1) + XE(1,1)*XE(2,2)
  164. & - XE(1,1)*XE(2,4) - XE(1,3)*XE(2,2) + XE(1,2)*XE(2,3)
  165. & + XE(1,3)*XE(2,4) + XE(1,4)*XE(2,1)
  166. * AIRE = X1s2 * r_z
  167. r_z = XUn / r_z
  168. BB(1,1) = r_z * (XE(2,2)-XE(2,4))
  169. BB(2,2) = r_z * (XE(1,4)-XE(1,2))
  170. BB(1,3) = r_z * (XE(2,3)-XE(2,1))
  171. BB(2,4) = r_z * (XE(1,1)-XE(1,3))
  172. BB(1,5) = r_z * (XE(2,4)-XE(2,2))
  173. BB(2,6) = r_z * (XE(1,2)-XE(1,4))
  174. BB(1,7) = r_z * (XE(2,1)-XE(2,3))
  175. BB(2,8) = r_z * (XE(1,3)-XE(1,1))
  176. * write (*,*) 'calcul de aire pour bbdil'
  177. * write (*,*) 'aire = ',aire
  178. * write (*,*) 'ecriture de la matrice bb'
  179. * write (*,4) (('bb(',i,',',j,')=',bb(i,j),i=1,3),j=1,2*nbnn)
  180. GOTO 666
  181. *
  182. *--------axisymétriques(ifour=0)
  183. 203 CONTINUE
  184. r_z = XE(1,1)*( XE(1,1)*XE(2,2)-XE(1,1)*XE(2,4)-XE(1,4)*XE(2,4)
  185. & +XE(1,4)*XE(2,1)-XE(1,2)*XE(2,1)+XE(1,2)*XE(2,2))
  186. & + XE(1,2)*(-XE(1,2)*XE(2,1)+XE(1,2)*XE(2,3)+XE(1,3)*XE(2,3)
  187. & -XE(1,3)*XE(2,2))
  188. & + XE(1,3)*( XE(1,3)*XE(2,4)+XE(1,4)*XE(2,4)-XE(1,4)*XE(2,3)
  189. & -XE(1,3)*XE(2,2))
  190. & + XE(1,4)*XE(1,4)*(XE(2,1)-XE(2,3))
  191. * AIRE = x1s6 * r_z
  192. r_z = xUn / r_z
  193. BB(1,1) = ( XE(2,2)*(XE(1,2)+XE(1,1))-XE(2,4)*(XE(1,1)+XE(1,4))
  194. & +X1s2*( XE(1,3)*(XE(2,2)-XE(2,4))
  195. & +XE(1,4)*(XE(2,2)+XE(2,3))
  196. & -XE(1,2)*(XE(2,3)+XE(2,4))) ) * r_z
  197. BB(2,2) = ( XE(1,4)*(XE(1,1)+XE(1,4))
  198. & -XE(1,2)*(XE(1,1)+XE(1,2)) ) * r_z
  199. BB(1,3) = ( XE(2,3)*(XE(1,2)+XE(1,3))-XE(2,1)*(XE(1,2)+XE(1,1))
  200. & +X1s2*( XE(1,4)*(XE(2,3)-XE(2,1))
  201. & +XE(1,1)*(XE(2,3)+XE(2,4))
  202. & -XE(1,3)*(XE(2,4)+XE(2,1))) ) * r_z
  203. BB(2,4) = ( XE(1,1)*(XE(1,2)+XE(1,1))
  204. & -XE(1,3)*(XE(1,3)+XE(1,2)) ) * r_z
  205. BB(1,5) = ( XE(2,4)*(XE(1,3)+XE(1,4))-XE(2,2)*(XE(1,3)+XE(1,2))
  206. & +X1s2*( XE(1,1)*(XE(2,4)-XE(2,2))
  207. & +XE(1,2)*(XE(2,4)+XE(2,1))
  208. & -XE(1,4)*(XE(2,1)+XE(2,2))) ) * r_z
  209. BB(2,6) = ( XE(1,2)*(XE(1,2)+XE(1,3))
  210. & -XE(1,4)*(XE(1,4)+XE(1,3)) ) * r_z
  211. BB(1,7) = ( XE(2,1)*(XE(1,4)+XE(1,1))-XE(2,3)*(XE(1,4)+XE(1,3))
  212. & +X1s2*( XE(1,2)*(XE(2,1)-XE(2,3))
  213. & +XE(1,3)*(XE(2,1)+XE(2,2))
  214. & -XE(1,1)*(XE(2,2)+XE(2,3))) ) * r_z
  215. BB(2,8) = ( XE(1,3)*(XE(1,3)+XE(1,4))
  216. & -XE(1,1)*(XE(1,1)+XE(1,4)) ) * r_z
  217. *- Integration analytique des termes BB(3,1),BB(3,3),BB(3,5),BB(3,7)
  218. BB(3,1) = ( XE(1,1)*(XE(2,2)-XE(2,4))+XE(2,1)*(XE(1,4)-XE(1,2))
  219. & +X1s2*( XE(1,2)*(XE(2,3)+XE(2,4))
  220. & +XE(1,3)*(XE(2,4)-XE(2,2))
  221. & -XE(1,4)*(XE(2,3)+XE(2,2))) ) * r_z
  222. BB(3,3) = ( XE(1,2)*(XE(2,3)-XE(2,1))+XE(2,2)*(XE(1,1)-XE(1,3))
  223. & +X1s2*( XE(1,3)*(XE(2,1)+XE(2,4))
  224. & +XE(1,4)*(XE(2,1)-XE(2,3))
  225. & -XE(1,1)*(XE(2,3)+XE(2,4))) ) * r_z
  226. BB(3,5) = ( XE(1,3)*(XE(2,4)-XE(2,2))+XE(2,3)*(XE(1,2)-XE(1,4))
  227. & +X1s2*( XE(1,4)*(XE(2,1)+XE(2,2))
  228. & +XE(1,1)*(XE(2,2)-XE(2,4))
  229. & -XE(1,2)*(XE(2,4)+XE(2,1))) ) * r_z
  230. BB(3,7) = ( XE(1,4)*(XE(2,1)-XE(2,3))+XE(2,4)*(XE(1,3)-XE(1,1))
  231. & +X1s2*( XE(1,1)*(XE(2,2)+XE(2,3))
  232. & +XE(1,2)*(XE(2,3)-XE(2,1))
  233. & -XE(1,3)*(XE(2,1)+XE(2,2))) ) * r_z
  234. * write (*,*) 'calcul de aire pour bbdil'
  235. * write (*,*) 'aire = ',aire
  236. * write (*,*) 'ecriture de la matrice bb'
  237. * write (*,4) (('bb(',i,',',j,')=',bb(i,j),i=1,3),j=1,2*nbnn)
  238. GOTO 666
  239.  
  240. *-----------------------------------------------------------------------
  241. *----------------------------- Element ICT6 ----------------------------
  242. *-----------------------------------------------------------------------
  243. 30 CONTINUE
  244. * write(*,*) 'Element ICT6',ifour
  245. * write (*,*) 'ecriture de la matrice xe[ ]'
  246. * write (*,3) (('xe(',i,',',j,')=',xe(i,j),i=1,2),j=1,6)
  247. IFR = ifour+4
  248. GOTO (666,301,301,303,666,999),IFR
  249. GOTO 999
  250. *
  251. *--------contraintes planes ou déformations planes(ifour=-2,-1)
  252. 301 CONTINUE
  253. *--------donnée des composantes de bb-dil
  254. r_z = ( XE(1,1)*(XE(2,2)-XE(2,6)) + XE(1,2)*(XE(2,3)-XE(2,1))
  255. & + XE(1,3)*(XE(2,4)-XE(2,2)) + XE(1,4)*(XE(2,5)-XE(2,3))
  256. & + XE(1,5)*(XE(2,6)-XE(2,4)) + XE(1,6)*(XE(2,1)-XE(2,5)) )
  257. & + X1s4 * ( XE(1,1)*(XE(2,5)-XE(2,3))
  258. & + XE(1,3)*(XE(2,1)-XE(2,5))
  259. & + XE(1,5)*(XE(2,3)-XE(2,1)) )
  260. AIRE = X2s3 * r_z
  261. r_z = xUn / r_z
  262. BB(1, 1) = ( (XE(2,2)-XE(2,6)) + X1s4 * (XE(2,5)-XE(2,3)) ) * r_z
  263. BB(2, 2) = ( (XE(1,6)-XE(1,2)) + X1s4 * (XE(1,3)-XE(1,5)) ) * r_z
  264. BB(1, 3) = (XE(2,3)-XE(2,1)) * r_z
  265. BB(2, 4) = (XE(1,1)-XE(1,3)) * r_z
  266. BB(1, 5) = ( (XE(2,4)-XE(2,2)) + X1s4 * (XE(2,1)-XE(2,5)) ) * r_z
  267. BB(2, 6) = ( (XE(1,2)-XE(1,4)) + X1s4 * (XE(1,5)-XE(1,1)) ) * r_z
  268. BB(1, 7) = (XE(2,5)-XE(2,3)) * r_z
  269. BB(2, 8) = (XE(1,3)-XE(1,5)) * r_z
  270. BB(1, 9) = ( (XE(2,6)-XE(2,4)) + X1s4 * (XE(2,3)-XE(2,1)) ) * r_z
  271. BB(2,10) = ( (XE(1,4)-XE(1,6)) + X1s4 * (XE(1,1)-XE(1,3)) ) * r_z
  272. BB(1,11) = (XE(2,1)-XE(2,5)) * r_z
  273. BB(2,12) = (XE(1,5)-XE(1,1)) * r_z
  274.  
  275. * write (*,*) 'calcul de aire pour bbdil'
  276. * write (*,*) 'aire = ',aire
  277. * write (*,*) 'ecriture de la matrice bb'
  278. * write (*,4) (('bb(',i,',',j,')=',bb(i,j),i=1,3),j=1,2*nbnn)
  279. GOTO 666
  280. *
  281. *--------axisymétriques (ifour=0)
  282. 303 CONTINUE
  283. *--------donnée des composantes de bb-dil
  284. *-----intégration numérique des composantes bb-dil------
  285. AIRE = XZer
  286. DO i = 1,nbnn
  287. XNUM1(i) = XZer
  288. XNUM2(i) = XZer
  289. XNUM3(i) = XZer
  290. ENDDO
  291. DO 3031 IGAU = 1,NBPGAU
  292. DO i = 1,nbnn
  293. SHP(1,i) = SHPTOT(1,i,IGAU)
  294. SHP(2,i) = SHPTOT(2,i,IGAU)
  295. SHP(3,i) = SHPTOT(3,i,IGAU)
  296. ENDDO
  297. CALL DISTRR(XE,SHP,nbnn,RR)
  298. CALL DISTR2(XE,SHP,nbnn,RR2)
  299. CALL DISTR3(XE,SHP,nbnn,RR3)
  300. CALL DISTZ2(XE,SHP,nbnn,ZZ2)
  301. CALL DISTZ3(XE,SHP,nbnn,ZZ3)
  302. CALL JACOBI(XE,SHP,2,nbnn,DJAC)
  303. r_z = POIGAU(IGAU) * RR
  304. r_zz2 = r_z * ZZ2
  305. r_zz3 = r_z * ZZ3
  306. r_rr2 = r_z * RR2
  307. r_rr3 = r_z * RR3
  308. r_z = POIGAU(IGAU) * DJAC
  309. AIRE = AIRE + (r_z * RR)
  310. DO i = 1,nbnn
  311. XNUM1(i) = XNUM1(i) + ( SHPTOT(2,i,IGAU)*r_zz3
  312. & - SHPTOT(3,i,IGAU)*r_zz2 )
  313. XNUM2(i) = XNUM2(i) + ( SHPTOT(3,i,IGAU)*r_rr2
  314. & - SHPTOT(2,i,IGAU)*r_rr3 )
  315. XNUM3(i) = XNUM3(i) + ( SHPTOT(1,i,IGAU)*r_z )
  316. ENDDO
  317. 3031 CONTINUE
  318. r_z = XUn / AIRE
  319. DO i = 1,nbnn
  320. k = 2*i
  321. BB(1,k-1) = XNUM1(i) * r_z
  322. BB(2,k) = XNUM2(i) * r_z
  323. BB(3,k-1) = XNUM3(i) * r_z
  324. ENDDO
  325. * write (*,*) 'calcul de aire pour bbdil'
  326. * write (*,*) 'aire = ',aire
  327. * write (*,*) 'ecriture de la matrice bb'
  328. * write (*,4) (('bb(',i,',',j,')=',bb(i,j),i=1,3),j=1,2*nbnn)
  329. GOTO 666
  330.  
  331. *-----------------------------------------------------------------------
  332. *----------------------------- Element ICQ8 ----------------------------
  333. *-----------------------------------------------------------------------
  334. 40 CONTINUE
  335. * write(*,*) 'Element ICQ8',ifour
  336. * write (*,*) 'Ecriture de la matrice xe[ ]'
  337. * write (*,3) (('xe(',i,',',j,')=',xe(i,j),i=1,2),j=1,6)
  338.  
  339. IFR = ifour+4
  340. GOTO (666,401,401,403,666,999),IFR
  341. GOTO 666
  342.  
  343. *--------contraintes planes ou déformations planes (ifour=-2,-1)
  344. 401 CONTINUE
  345. *--------donnée de la matrice symétrique p(3,3) (1 fois par élément)
  346. P(1,1) = X2s3 * ( XE(1,1)*(XE(2,2)-XE(2,8))
  347. 1 + XE(1,2)*(XE(2,3)-XE(2,1))
  348. 1 +XE(1,3)*(XE(2,4)-XE(2,2))+XE(1,4)*(XE(2,5)-XE(2,3))
  349. 2 +XE(1,5)*(XE(2,6)-XE(2,4))+XE(1,6)*(XE(2,7)-XE(2,5))
  350. 3 +XE(1,7)*(XE(2,8)-XE(2,6))+XE(1,8)*(XE(2,1)-XE(2,7)) )
  351. 4 + X1s6 * ( (XE(1,1)-XE(1,5))*(XE(2,7)-XE(2,3))
  352. 5 +(XE(1,3)-XE(1,7))*(XE(2,1)-XE(2,5)) )
  353. P(1,2) = 5.D0/9.D0*(-XE(2,2)*(XE(1,1)+XE(1,3))
  354. 1 +XE(1,2)*(XE(2,1)+XE(2,3))+XE(1,6)*(-XE(2,5)-XE(2,7))
  355. 1 +XE(2,6)*(XE(1,5)+XE(1,7)))
  356. 2 +4.D0/9.D0*(XE(1,8)*(XE(2,7)+XE(2,2)-XE(2,6)-XE(2,1))
  357. 2 +XE(2,8)*(-XE(1,7)-XE(1,2)+XE(1,6)+XE(1,1))
  358. 2 +XE(1,4)*(XE(2,2)-XE(2,6)+XE(2,5)-XE(2,3))
  359. 2 +XE(2,4)*(-XE(1,2)+XE(1,6)-XE(1,5)+XE(1,3)))
  360. 3 +7.D0/45.D0*(XE(1,2)*(XE(2,5)+XE(2,7))
  361. 3 -XE(2,2)*(XE(1,5)+XE(1,7))-XE(1,6)*(XE(2,1)+XE(2,3))
  362. 3 +XE(2,6)*(XE(1,1)+XE(1,3)))
  363. 4 +8.D0/15.D0*(XE(1,6)*XE(2,2)-XE(1,2)*XE(2,6))
  364. 5 +7.D0/90.D0*(XE(1,7)*XE(2,1)-XE(1,3)*XE(2,5)
  365. 5 +XE(1,5)*XE(2,3)-XE(1,1)*XE(2,7))
  366. 6 +1.D0/30.D0*(-XE(1,7)*XE(2,3)+XE(1,3)*XE(2,7)
  367. 6 -XE(1,5)*XE(2,1)+XE(1,1)*XE(2,5))
  368. P(1,3) = 5.D0/9.D0*(-XE(1,8)*(XE(2,7)+XE(2,1))
  369. 1 +XE(2,8)*(XE(1,7)+XE(1,1))
  370. 1 +XE(1,4)*(XE(2,3)+XE(2,5))
  371. 1 -XE(2,4)*(XE(1,3)+XE(1,5)))
  372. 2 +4.D0/9.D0*((XE(1,8)-XE(1,4))*(XE(2,2)+XE(2,6))
  373. 2 +(-XE(2,8)+XE(2,4))*(XE(1,2)+XE(1,6))
  374. 2 +XE(1,2)*(XE(2,1)-XE(2,3))+XE(2,2)*(-XE(1,1)+XE(1,3))
  375. 2 +XE(1,6)*(XE(2,7)-XE(2,5))+XE(2,6)*(-XE(1,7)+XE(1,5)))
  376. 3 +7.D0/45.D0*(-XE(1,8)*(XE(2,3)+XE(2,5))
  377. 3 +XE(2,8)*(XE(1,3)+XE(1,5))
  378. 3 +XE(1,4)*(XE(2,1)+XE(2,7))
  379. 3 -XE(2,4)*(XE(1,1)+XE(1,7)))
  380. 4 +1.D0/30.D0*(-XE(1,7)*XE(2,3)+XE(1,3)*XE(2,7)
  381. 4 +XE(1,5)*XE(2,1)-XE(1,1)*XE(2,5))
  382. 5 +7.D0/90.D0*(XE(1,7)*XE(2,5)-XE(1,3)*XE(2,1)
  383. 5 -XE(1,5)*XE(2,7)+XE(1,1)*XE(2,3))
  384. 6 +8.D0/15.D0*(XE(1,8)*XE(2,4)-XE(1,4)*XE(2,8))
  385. P(1,4) = XZer
  386.  
  387. P(2,1) = P(1,2)
  388. P(2,2) = 16.D0/45.D0*(XE(2,6)*(XE(1,5)-XE(1,7))
  389. 1 +XE(1,6)*(XE(2,7)-XE(2,5))+XE(1,2)*(XE(2,3)-XE(2,1))
  390. 1 +XE(2,2)*(XE(1,1)-XE(1,3)))
  391. 2 +14.D0/45.D0*(XE(2,4)*(XE(1,3)-XE(1,5))
  392. 2 +XE(1,4)*(XE(2,5)-XE(2,3))+XE(1,8)*(XE(2,1)-XE(2,7))
  393. 2 +XE(2,8)*(XE(1,7)-XE(1,1)))
  394. 3 +8.D0/45.D0*(XE(1,4)*(XE(2,2)-XE(2,6))
  395. 3 +XE(1,8)*(XE(2,6)-XE(2,2))+XE(1,2)*(XE(2,8)-XE(2,4))
  396. 3 +XE(1,6)*(XE(2,4)-XE(2,8)))
  397. 4 +2.D0/45.D0*(XE(2,2)*(XE(1,5)-XE(1,7))
  398. 4 +XE(1,2)*(XE(2,7)-XE(2,5)) +XE(1,6)*(XE(2,3)-XE(2,1))
  399. 4 +XE(2,6)*(XE(1,1)-XE(1,3)))
  400. 5 +4.D0/45.D0*(XE(2,4)*(XE(1,1)-XE(1,7))
  401. 5 +XE(1,4)*(XE(2,7)-XE(2,1)) +XE(1,8)*(XE(2,3)-XE(2,5))
  402. 5 +XE(2,8)*(XE(1,5)-XE(1,3)))
  403. 6 +17.D0/90.D0*(XE(1,7)*XE(2,5)+XE(1,3)*XE(2,1)
  404. 6 -XE(1,5)*XE(2,7)-XE(1,1)*XE(2,3))
  405. 7 +1.D0/90.D0*(-XE(1,7)*XE(2,1)-XE(1,3)*XE(2,5)
  406. 7 +XE(1,5)*XE(2,3)+XE(1,1)*XE(2,7))
  407. P(2,3) = 1.D0/3.D0*(XE(1,5)*(XE(2,6)-XE(2,4))
  408. 1 +XE(1,8)*(XE(2,7)+XE(2,1))+XE(1,7)*(XE(2,6)-XE(2,8))
  409. 1 -XE(1,2)*(XE(2,1)+XE(2,3))+XE(1,1)*(XE(2,2)-XE(2,8))
  410. 1 +XE(1,4)*(XE(2,5)+XE(2,3))-XE(1,6)*(XE(2,7)+XE(2,5))
  411. 1 +XE(1,3)*(XE(2,2)-XE(2,4)))
  412. 2 +4.D0/9.D0*(-XE(1,4)*(XE(2,2)+XE(2,6))
  413. 2 -XE(1,8)*(XE(2,6)+XE(2,2))+XE(1,2)*(XE(2,8)+XE(2,4))
  414. 2 +XE(1,6)*(XE(2,4)+XE(2,8)))
  415. 3 +1.D0/9.D0*(XE(1,7)*(XE(2,2)-XE(2,4))
  416. 3 +XE(1,8)*(XE(2,3)+XE(2,5))
  417. 3 +XE(1,5)*(XE(2,2)-XE(2,8))+XE(1,3)*(XE(2,6)-XE(2,8))
  418. 3 -XE(1,2)*(XE(2,5)+XE(2,7))+XE(1,1)*(XE(2,6)-XE(2,4))
  419. 3 +XE(1,4)*(XE(2,1)+XE(2,7))-XE(1,6)*(XE(2,1)+XE(2,3)))
  420. P(2,4) = XZer
  421.  
  422. P(3,1) = P(1,3)
  423. P(3,2) = P(2,3)
  424. P(3,3) = 14.D0/45.D0*(XE(1,6)*(XE(2,7)-XE(2,5))
  425. 1 +XE(2,6)*(XE(1,5)-XE(1,7))+XE(1,2)*(XE(2,3)-XE(2,1))
  426. 1 +XE(2,2)*(XE(1,1)-XE(1,3)))
  427. 2 +16.D0/45.D0*(XE(1,8)*(XE(2,1)-XE(2,7))
  428. 2 +XE(2,8)*(XE(1,7)-XE(1,1))+XE(1,4)*(XE(2,5)-XE(2,3))
  429. 2 +XE(2,4)*(XE(1,3)-XE(1,5)))
  430. 3 +8.D0/45.D0*(XE(1,4)*(XE(2,2)-XE(2,6))
  431. 3 +XE(1,8)*(XE(2,6)-XE(2,2))+XE(1,2)*(XE(2,8)-XE(2,4))
  432. 3 +XE(1,6)*(XE(2,4)-XE(2,8)))
  433. 4 +2.D0/45.D0*(XE(2,4)*(XE(1,7)-XE(1,1))
  434. 4 +XE(1,8)*(XE(2,5)-XE(2,3))+XE(2,8)*(-XE(1,5)+XE(1,3))
  435. 4 +XE(1,4)*(XE(2,1)-XE(2,7)))
  436. 5 +4.D0/45.D0*(XE(2,2)*(XE(1,7)-XE(1,5))
  437. 5 +XE(1,2)*(XE(2,5)-XE(2,7))+XE(1,6)*(XE(2,1)-XE(2,3))
  438. 5 +XE(2,6)*(XE(1,3)-XE(1,1)))
  439. 6 +1.D0/90.D0*(-XE(1,5)*XE(2,7)+XE(1,3)*XE(2,1)
  440. 6 +XE(1,7)*XE(2,5)-XE(1,1)*XE(2,3))
  441. 7 +17.D0/90.D0*(-XE(1,7)*XE(2,1)-XE(1,3)*XE(2,5)
  442. 7 +XE(1,5)*XE(2,3)+XE(1,1)*XE(2,7))
  443. P(3,4) = XZer
  444.  
  445. P(4,1) = XZer
  446. P(4,2) = XZer
  447. P(4,3) = XZer
  448. P(4,4) = XZer
  449.  
  450. * write (*,*) 'ecriture de la matrice p[ ]'
  451. * write (*,4) (('p(',i,',',j,')=',p(i,j),i=1,3),j=1,3)
  452.  
  453. *-----------------donnée des vecteurs {t}={t1,t2,t3}--------------------
  454. A(1,1) = 2.D0/3.D0*(XE(2,2)-XE(2,8))
  455. 1 +1.D0/6.D0*(XE(2,7)-XE(2,3))
  456. A(2,1) = 1.D0/9.D0*(4.D0*XE(2,8)-5.D0*XE(2,2))
  457. 1 +7.D0/90.D0*(-XE(2,7)+2.D0*XE(2,6))
  458. 1 +1.D0/30.D0*XE(2,5)
  459. A(3,1) = 1.D0/9.D0*(-4.D0*XE(2,2)+5.D0*XE(2,8))
  460. 1 +7.D0/90.D0*(XE(2,3)-2.D0*XE(2,4))
  461. 1 -1.D0/30.D0*XE(2,5)
  462.  
  463. A(1,2) = 2.D0/3.D0*(XE(1,8)-XE(1,2))
  464. 1 +1.D0/6.D0*(XE(1,3)-XE(1,7))
  465. A(2,2) = 1.D0/9.D0*(-4.D0*XE(1,8)+5.D0*XE(1,2))
  466. 1 +7.D0/90.D0*(XE(1,7)-2.D0*XE(1,6))
  467. 1 -1.D0/30.D0*XE(1,5)
  468. A(3,2) = 1.D0/9.D0*(4.D0*XE(1,2)-5.D0*XE(1,8))
  469. 1 +7.D0/90.D0*(-XE(1,3)+2.D0*XE(1,4))
  470. 1 +1.D0/30.D0*XE(1,5)
  471.  
  472. A(1,3) = 2.D0/3.D0*(-XE(2,1)+XE(2,3))
  473. A(2,3) = 5.D0/9.D0*(XE(2,1)+XE(2,3))
  474. 1 +7.D0/45.D0*(XE(2,5)+XE(2,7))
  475. 1 -4.D0/9.D0*(XE(2,4)+XE(2,8))
  476. 1 -8.D0/15.D0*XE(2,6)
  477. A(3,3) = 4.D0/9.D0*(XE(2,1)-XE(2,3)+XE(2,4)-XE(2,8))
  478.  
  479. A(1,4) = 2.D0/3.D0*(-XE(1,3)+XE(1,1))
  480. A(2,4) = -5.D0/9.D0*(XE(1,1)+XE(1,3))
  481. 1 -7.D0/45.D0*(XE(1,5)+XE(1,7))
  482. 1 +4.D0/9.D0*(XE(1,4)+XE(1,8))
  483. 1 +8.D0/15.D0*XE(1,6)
  484. A(3,4) = 4.D0/9.D0*(-XE(1,1)+XE(1,3)-XE(1,4)+XE(1,8))
  485.  
  486. A(1,5) = 2.D0/3.D0*(XE(2,4)-XE(2,2))
  487. 1 +1.D0/6.D0*(-XE(2,5)+XE(2,1))
  488. A(2,5) = 1.D0/9.D0*(4.D0*XE(2,4)-5.D0*XE(2,2))
  489. 1 +7.D0/90.D0*(-XE(2,5)+2.D0*XE(2,6))
  490. 1 +1.D0/30.D0*XE(2,7)
  491. A(3,5) = 1.D0/9.D0*(4.D0*XE(2,2)-5.D0*XE(2,4))
  492. 1 +7.D0/90.D0*(-XE(2,1)+2.D0*XE(2,8))
  493. 1 +1.D0/30.D0*XE(2,7)
  494.  
  495. A(1,6) = 2.D0/3.D0*(-XE(1,4)+XE(1,2))
  496. 1 +1.D0/6.D0*(XE(1,5)-XE(1,1))
  497. A(2,6) = 1.D0/9.D0*(-4.D0*XE(1,4)+5.D0*XE(1,2))
  498. 1 +7.D0/90.D0*(XE(1,5)-2.D0*XE(1,6))
  499. 1 -1.D0/30.D0*XE(1,7)
  500. A(3,6) = 1.D0/9.D0*(-4.D0*XE(1,2)+5.D0*XE(1,4))
  501. 1 +7.D0/90.D0*(XE(1,1)-2.D0*XE(1,8))
  502. 1 -1.D0/30.D0*XE(1,7)
  503.  
  504. A(1,7) = 2.D0/3.D0*(-XE(2,3)+XE(2,5))
  505. A(2,7) = 4.D0/9.D0*(XE(2,5)-XE(2,3)+XE(2,2)-XE(2,6))
  506. A(3,7) = 5.D0/9.D0*(XE(2,3)+XE(2,5))
  507. 1 +7.D0/45.D0*(XE(2,1)+XE(2,7))
  508. 1 -4.D0/9.D0*(XE(2,2)+XE(2,6))
  509. 1 -8.D0/15.D0*XE(2,8)
  510.  
  511. A(1,8) = 2.D0/3.D0*(XE(1,3)-XE(1,5))
  512. A(2,8) = 4.D0/9.D0*(-XE(1,5)+XE(1,3)-XE(1,2)+XE(1,6))
  513. A(3,8) = -5.D0/9.D0*(XE(1,3)+XE(1,5))
  514. 1 -7.D0/45.D0*(XE(1,1)+XE(1,7))
  515. 1 +4.D0/9.D0*(XE(1,2)+XE(1,6))
  516. 1 +8.D0/15.D0*XE(1,8)
  517.  
  518. A(1,9) = 2.D0/3.D0*(-XE(2,4)+XE(2,6))
  519. 1 +1.D0/6.D0*(XE(2,3)-XE(2,7))
  520. A(2,9) = 1.D0/9.D0*(-4.D0*XE(2,4)+5.D0*XE(2,6))
  521. 1 +7.D0/90.D0*(XE(2,3)-2.D0*XE(2,2))
  522. 1 -1.D0/30.D0*XE(2,1)
  523. A(3,9) = 1.D0/9.D0*(4.D0*XE(2,6)-5.D0*XE(2,4))
  524. 1 +7.D0/90.D0*(-XE(2,7)+2.D0*XE(2,8))
  525. 1 +1.D0/30.D0*XE(2,1)
  526.  
  527. A(1,10) = 2.D0/3.D0*(XE(1,4)-XE(1,6))
  528. 1 +1.D0/6.D0*(-XE(1,3)+XE(1,7))
  529. A(2,10) = 1.D0/9.D0*(4.D0*XE(1,4)-5.D0*XE(1,6))
  530. 1 +7.D0/90.D0*(-XE(1,3)+2.D0*XE(1,2))
  531. 1 +1.D0/30.D0*XE(1,1)
  532. A(3,10) = 1.D0/9.D0*(-4.D0*XE(1,6)+5.D0*XE(1,4))
  533. 1 +7.D0/90.D0*(XE(1,7)-2.D0*XE(1,8))
  534. 1 -1.D0/30.D0*XE(1,1)
  535.  
  536. A(1,11) = 2.D0/3.D0*(-XE(2,5)+XE(2,7))
  537. A(2,11) = -5.D0/9.D0*(XE(2,7)+XE(2,5))
  538. 1 -7.D0/45.D0*(XE(2,1)+XE(2,3))
  539. 1 +4.D0/9.D0*(XE(2,4)+XE(2,8))
  540. 1 +8.D0/15.D0*XE(2,2)
  541. A(3,11) = 4.D0/9.D0*(XE(2,4)-XE(2,5)-XE(2,8)+XE(2,7))
  542.  
  543. A(1,12) = 2.D0/3.D0*(XE(1,5)-XE(1,7))
  544. A(2,12) = 5.D0/9.D0*(XE(1,7)+XE(1,5))
  545. 1 +7.D0/45.D0*(XE(1,1)+XE(1,3))
  546. 1 -4.D0/9.D0*(XE(1,4)+XE(1,8))
  547. 1 -8.D0/15.D0*XE(1,2)
  548. A(3,12) = 4.D0/9.D0*(-XE(1,4)+XE(1,5)+XE(1,8)-XE(1,7))
  549.  
  550. A(1,13) = 2.D0/3.D0*(-XE(2,6)+XE(2,8))
  551. 1 +1.D0/6.D0*(XE(2,5)-XE(2,1))
  552. A(2,13) = 1.D0/9.D0*(-4.D0*XE(2,8)+5.D0*XE(2,6))
  553. 1 +7.D0/90.D0*(XE(2,1)-2.D0*XE(2,2))
  554. 1 -1.D0/30.D0*XE(2,3)
  555. A(3,13) = 1.D0/9.D0*(-4.D0*XE(2,6)+5.D0*XE(2,8))
  556. 1 +7.D0/90.D0*(XE(2,5)-2.D0*XE(2,4))
  557. 1 -1.D0/30.D0*XE(2,3)
  558.  
  559. A(1,14) = 2.D0/3.D0*(XE(1,6)-XE(1,8))
  560. 1 +1.D0/6.D0*(-XE(1,5)+XE(1,1))
  561. A(2,14) = 1.D0/9.D0*(4.D0*XE(1,8)-5.D0*XE(1,6))
  562. 1 +7.D0/90.D0*(-XE(1,1)+2.D0*XE(1,2))
  563. 1 +1.D0/30.D0*XE(1,3)
  564. A(3,14) = 1.D0/9.D0*(4.D0*XE(1,6)-5.D0*XE(1,8))
  565. 1 +7.D0/90.D0*(-XE(1,5)+2.D0*XE(1,4))
  566. 1 +1.D0/30.D0*XE(1,3)
  567.  
  568. A(1,15) = 2.D0/3.D0*(-XE(2,7)+XE(2,1))
  569. A(2,15) = 4.D0/9.D0*(-XE(2,6)+XE(2,7)+XE(2,2)-XE(2,1))
  570. A(3,15) = -5.D0/9.D0*(XE(2,1)+XE(2,7))
  571. 1 -7.D0/45.D0*(XE(2,3)+XE(2,5))
  572. 1 +4.D0/9.D0*(XE(2,2)+XE(2,6))
  573. 1 +8.D0/15.D0*XE(2,4)
  574.  
  575. A(1,16) = 2.D0/3.D0*(XE(1,7)-XE(1,1))
  576. A(2,16) = 4.D0/9.D0*(XE(1,6)-XE(1,7)-XE(1,2)+XE(1,1))
  577. A(3,16) = 5.D0/9.D0*(XE(1,1)+XE(1,7))
  578. 1 +7.D0/45.D0*(XE(1,3)+XE(1,5))
  579. 1 -4.D0/9.D0*(XE(1,2)+XE(1,6))
  580. 1 -8.D0/15.D0*XE(1,4)
  581.  
  582. * write (*,*) 'ecriture de la matrice t[ ]'
  583. * write (*,4) (('t(',i,',',j,')=',a(i,j),i=1,3),j=1,16)
  584. * if (IB.EQ.34.and.mele.eq.72) then
  585. * write(*,*)'bbcalc'
  586. * do j = 1,16
  587. * write(*,*) (a(i,j),i=1,3)
  588. * enddo
  589. * do j = 1,4
  590. * write(6,*) (p(i,j),i=1,4)
  591. * enddo
  592. * endif
  593. *--------------résolution matricielle de [p]{a}={t}---------------------
  594. * write (*,*) 'Appel a GAUSSJ'
  595. N = 3
  596. NP = 4
  597. M = 16
  598. MP = 60
  599. CALL GAUSSJ(P,N,NP,A,M,MP)
  600. GOTO 666
  601.  
  602. *--------axisymétriques (ifour=0)
  603. 403 CONTINUE
  604. *--------donnée de la matrice symétrique p(3,3) (1 fois par élément)
  605. CALL ZERO(P,4,4)
  606. DO IGAU = 1,NBPGAU
  607. DO I = 1,nbnn
  608. SHP(1,I) = SHPTOT(1,I,IGAU)
  609. SHP(2,I) = SHPTOT(2,I,IGAU)
  610. SHP(3,I) = SHPTOT(3,I,IGAU)
  611. ENDDO
  612. CALL DISTRR(XE,SHP,nbnn,RR)
  613. CALL JACOBI(XE,SHP,2,nbnn,DJAC)
  614. r_z = POIGAU(IGAU)*RR*DJAC
  615. P(1,1) = P(1,1) + r_z
  616. P(2,1) = P(2,1) + (r_z * QSIGAU(IGAU))
  617. P(3,1) = P(3,1) + (r_z * ETAGAU(IGAU))
  618. P(2,2) = P(2,2) + (r_z * QSIGAU(IGAU) * QSIGAU(IGAU))
  619. P(3,2) = P(3,2) + (r_z * QSIGAU(IGAU) * ETAGAU(IGAU))
  620. P(3,3) = P(3,3) + (r_z * ETAGAU(IGAU) * ETAGAU(IGAU))
  621. ENDDO
  622. P(1,2) = P(2,1)
  623. P(1,3) = P(3,1)
  624. P(2,3) = P(3,2)
  625. * write (*,*) 'ecriture de la matrice p[ ]'
  626. * write (*,4) (('p(',i,',',j,')=',p(i,j),i=1,3),j=1,3)
  627. *-----------calcul des vecteurs {t}={t1,t2,t3}--intégration numérique---
  628. CALL ZERO(A,4,24)
  629. DO 4031 IGAU = 1,NBPGAU
  630. DO I = 1,nbnn
  631. SHP(1,I) = SHPTOT(1,I,IGAU)
  632. SHP(2,I) = SHPTOT(2,I,IGAU)
  633. SHP(3,I) = SHPTOT(3,I,IGAU)
  634. ENDDO
  635. CALL DISTRR(XE,SHP,nbnn,RR)
  636. CALL DISTR2(XE,SHP,nbnn,RR2)
  637. CALL DISTR3(XE,SHP,nbnn,RR3)
  638. CALL DISTZ2(XE,SHP,nbnn,ZZ2)
  639. CALL DISTZ3(XE,SHP,nbnn,ZZ3)
  640. CALL JACOBI(XE,SHP,2,nbnn,DJAC)
  641. K = 0
  642. DO I = 1,nbnn
  643. K = K+1
  644. r_z = POIGAU(IGAU)*RR
  645. & *(SHPTOT(2,I,IGAU)*ZZ3-SHPTOT(3,I,IGAU)*ZZ2)
  646. A(1,K) = A(1,K) + r_z
  647. A(2,K) = A(2,K) + r_z * QSIGAU(IGAU)
  648. A(3,K) = A(3,K) + r_z * ETAGAU(IGAU)
  649. K = K+1
  650. r_z = POIGAU(IGAU)*RR
  651. & *(-SHPTOT(2,I,IGAU)*RR3+SHPTOT(3,I,IGAU)*RR2)
  652. A(1,K) = A(1,K) + r_z
  653. A(2,K) = A(2,K) + r_z * QSIGAU(IGAU)
  654. A(3,K) = A(3,K) + r_z * ETAGAU(IGAU)
  655. K = K+1
  656. r_z = POIGAU(IGAU)*SHPTOT(1,I,IGAU)*DJAC
  657. A(1,K) = A(1,K) + r_z
  658. A(2,K) = A(2,K) + r_z * QSIGAU(IGAU)
  659. A(3,K) = A(3,K) + r_z * ETAGAU(IGAU)
  660. ENDDO
  661. 4031 CONTINUE
  662. * write(*,*) 'Ecriture de la matrice t[ ]'
  663. * write(*,4) (('t(',i,',',j,')=',a(i,j),i=1,3),j=1,24)
  664. *------------résolution matricielle de [p]{a}={t}-----------------
  665. N = 3
  666. NP = 4
  667. M = 24
  668. MP = 60
  669. * write (*,*) 'appel a GAUSSJ'
  670. CALL GAUSSJ(P,N,NP,A,M,MP)
  671. * write(*,*) 'Ecriture de la matrice a[ ]'
  672. * write(*,4) (('a(',i,',',j,')=',a(i,j),i=1,3),j=1,24)
  673. GOTO 666
  674.  
  675. C=======================================================================
  676. C========== Elements MASSIFS INCOMPRESSIBLES TRIDIMENSIONNELS ==========
  677. C=======================================================================
  678.  
  679. 100 CONTINUE
  680. IF (ifour.NE.2) GOTO 999
  681.  
  682. DIM3=1.d0
  683. FACAR=1.D0
  684.  
  685. C= 3.12.3 - Boucle sur les points d'integration
  686. ESTEL=XZero
  687. DO iGau=1,NBPGAU
  688. * IBMN=MIN(IB,MELVA1.VELCHE(/2))
  689. * IGMN=MIN(iGau,MELVA1.VELCHE(/1))
  690. * FACSCA=MELVA1.VELCHE(IGMN,IBMN)
  691. FACSCA = 1.D0
  692.  
  693. * BGENE au point IGAU
  694. * if (I2B.EQ.1.and.I2GAU.EQ.1) then
  695. * write(6,*) IGAU,NBPGAU,POIGAU(1),QSIGAU(1),ETAGAU(1),DZEGAU(1)
  696. * endif
  697. CALL BMATST(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  698. 1 MELE,MFR,NBNN,LRE,IFOUR,NSTRS,NHRM,DIM3,
  699. 2 XE,SHPTOT,SHP,BGENE,DJAC,XDPGE,YDPGE)
  700.  
  701. * non actif en 3D
  702. if(.false.) then
  703. C= 3.12.3.6 - Caracteristiques pour les elements MASSIFS (EPAISSEUR en contrainte plane)
  704. IF (IVACAR.EQ.0) GOTO 80
  705. MPTVAL=IVACAR
  706. MELVAL=IVAL(1)
  707. IF (MELVAL.NE.0) THEN
  708. IGMN=MIN(iGau,VELCHE(/1))
  709. IBMN=MIN(IB,VELCHE(/2))
  710. FACAR=VELCHE(IGMN,IBMN)
  711. ENDIF
  712. GOTO 80
  713. endif
  714.  
  715. C= 3.12.3.9 - Calcul de la composante integree en ce point de Gauss
  716. 80 DJAC=ABS(DJAC)*POIGAU(iGau)*FACAR*DIM3
  717. ESTEL=ESTEL+FACSCA*DJAC
  718.  
  719. if (MELE.EQ.72.or.(MELE.ge.76.and.MELE.le.78).or.
  720. &MELE.EQ.274) then
  721. * second membre T
  722. K = 1
  723. jj= 0
  724. DO i = 1,NBNN
  725. jj = jj + 1
  726. A(1,jj) = A(1,jj) + BGENE(1,K)*DJAC
  727. A(2,jj) = A(2,jj) + BGENE(1,K)*qsigau(igau)*DJAC
  728. A(3,jj) = A(3,jj) + BGENE(1,K)*etagau(igau)*DJAC
  729. if(mele.ge.76)
  730. & A(4,jj) = A(4,jj) + BGENE(1,K)*dzegau(igau)*DJAC
  731. jj = jj + 1
  732. A(1,jj) = A(1,jj) + BGENE(2,K+1)*DJAC
  733. A(2,jj) = A(2,jj) + BGENE(2,K+1)*qsigau(igau)*DJAC
  734. A(3,jj) = A(3,jj) + BGENE(2,K+1)*etagau(igau)*DJAC
  735. if(mele.ge.76) then
  736. A(4,jj) = A(4,jj) + BGENE(2,K+1)*dzegau(igau)*DJAC
  737. jj = jj + 1
  738. A(1,jj) = A(1,jj) + BGENE(3,K+2)*DJAC
  739. A(2,jj) = A(2,jj) + BGENE(3,K+2)*qsigau(igau)*DJAC
  740. A(3,jj) = A(3,jj) + BGENE(3,K+2)*etagau(igau)*DJAC
  741. A(4,jj) = A(4,jj) + BGENE(3,K+2)*dzegau(igau)*DJAC
  742. endif
  743. K = K + KINC
  744. ENDDO
  745. * matrice M
  746. P(1,1) = P(1,1) + DJAC
  747. P(1,2) = P(1,2) + qsigau(igau)*DJAC
  748. P(1,3) = P(1,3) + etagau(igau)*DJAC
  749. P(2,2) = P(2,2) + qsigau(igau)*qsigau(igau)*DJAC
  750. P(2,3) = P(2,3) + qsigau(igau)*etagau(igau)*DJAC
  751.  
  752. P(3,3) = P(3,3) + etagau(igau)*etagau(igau)*DJAC
  753.  
  754. if(mele.ge.76) then
  755. P(1,4) = P(1,4) + dzegau(igau)*DJAC
  756.  
  757. P(2,4) = P(2,4) + qsigau(igau)*dzegau(igau)*DJAC
  758.  
  759. P(3,4) = P(3,4) + etagau(igau)*dzegau(igau)*DJAC
  760.  
  761. P(4,4) = P(4,4) + dzegau(igau)*dzegau(igau)*DJAC
  762. endif
  763.  
  764. else if(mele.le.71) then
  765.  
  766. K = 1
  767. DO i = 1,NBNN
  768. bb(1,K) = bb(1,K) + BGENE(1,K)*DJAC
  769. bb(2,K+1) = bb(2,K+1) + BGENE(2,K+1)*DJAC
  770. bb(3,K) = bb(3,K) + BGENE(3,K)*DJAC
  771. K = K + KINC
  772. ENDDO
  773.  
  774. elseif((mele.ge.73.and.mele.le.75)
  775. &.or.mele.eq.273) then
  776.  
  777. K = 1
  778. DO i = 1,NBNN
  779. bb(1,K) = bb(1,K) + BGENE(1,K)*DJAC
  780. bb(2,K+1) = bb(2,K+1) + BGENE(2,K+1)*DJAC
  781. bb(3,K+2) = bb(3,K+2) + BGENE(3,K+2)*DJAC
  782. K = K + KINC
  783. ENDDO
  784.  
  785. endif
  786.  
  787.  
  788. * if (I2B.EQ.34) then
  789. * write(6,*) 'djac',igau,DJAC,POIGAU(IGAU)
  790. * endif
  791.  
  792. ENDDO
  793.  
  794. if (MELE.EQ.72) then
  795. P(1,4) = XZer
  796.  
  797. P(2,1) = P(1,2)
  798. P(2,4) = XZer
  799. P(3,1) = P(1,3)
  800. P(3,2) = P(2,3)
  801. P(3,4) = XZer
  802.  
  803. P(4,1) = XZer
  804. P(4,2) = XZer
  805. P(4,3) = XZer
  806. P(4,4) = XZer
  807.  
  808. elseif(mele.ge.76) then
  809. P(2,1) = P(1,2)
  810. P(3,1) = P(1,3)
  811. P(3,2) = P(2,3)
  812. P(4,1) = P(1,4)
  813. P(4,2) = P(2,4)
  814. P(4,3) = P(3,4)
  815.  
  816. endif
  817.  
  818.  
  819. *
  820. ** synthese
  821. *
  822.  
  823. if (MELE.EQ.72) then
  824. * if (I2B.EQ.34) then
  825. * do j = 1,16
  826. * write(*,*) (a(i,j),i=1,3)
  827. * enddo
  828. * do j = 1,4
  829. * write(6,*) (p(i,j),i=1,4)
  830. * enddo
  831. * endif
  832. N = 3
  833. NP = 4
  834. M = 16
  835. MP = 60
  836. * write (*,*) 'appel a GAUSSJ'
  837. CALL GAUSSJ(P,N,NP,A,M,MP)
  838. * write(*,*) 'Ecriture de la matrice a[ ]'
  839. * write(*,4) (('a(',i,',',j,')=',a(i,j),i=1,3),j=1,24)
  840.  
  841. elseif (mele.ge.76.and.mele.ne.273) then
  842. N = 4
  843. NP = 4
  844. M = NBNN*3
  845. MP = 60
  846. * write (*,*) 'appel a GAUSSJ'
  847. CALL GAUSSJ(P,N,NP,A,M,MP)
  848.  
  849. elseif(mele.le.71) then
  850.  
  851. K = 1
  852. DO i = 1,NBNN
  853. bb(1,K) = bb(1,K)/ESTEL
  854. bb(2,K+1) = bb(2,K+1)/ESTEL
  855. bb(3,K) = bb(3,K)/ESTEL
  856. * if (I2B.EQ.34) then
  857. * write(6,*) 'bb',i,K, bb(1,K),bb(2,K+1),bb(3,K)
  858. * endif
  859. K = K + KINC
  860. ENDDO
  861.  
  862. elseif((mele.ge.73.and.mele.le.75).or.
  863. &mele.eq.273) then
  864.  
  865. K = 1
  866. DO i = 1,NBNN
  867. bb(1,K) = bb(1,K)/ESTEL
  868. bb(2,K+1) = bb(2,K+1)/ESTEL
  869. bb(3,K+2) = bb(3,K+2)/ESTEL
  870. * if (I2B.EQ.34) then
  871. * write(6,*) 'bb',i,K, bb(1,K),bb(2,K+1),bb(3,K+2)
  872. * endif
  873. K = K + KINC
  874. ENDDO
  875.  
  876. endif
  877.  
  878. RETURN
  879.  
  880. C=======================================================================
  881. C======================= ERREUR : Cas non prevus =======================
  882. C=======================================================================
  883. 999 CONTINUE
  884. CALL ERREUR(5)
  885.  
  886. 666 CONTINUE
  887. * write (*,*) 'Sortie de BBCAL2'
  888.  
  889. RETURN
  890.  
  891. END
  892.  
  893.  
  894.  

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