Télécharger kpcoq8.eso

Retour à la liste

Numérotation des lignes :

kpcoq8
  1. C KPCOQ8 SOURCE CHAT 05/01/13 01:04:13 5004
  2. SUBROUTINE KPCOQ8(XX, P, XKP, IANT)
  3. C
  4. C Procedure de calcul de la matrice Kppour un element COQ4
  5. C Entrees : XX(3, 8) : REAL*8 : Coordonnees des noeuds
  6. C XP(4) : REAL*8 : Pression aux points de Gauss
  7. C IANT : INTEGER : 1 si calcul asymétrique, 0 sinon
  8. C Sortie : XKP(48, 48) : REAL*8 : Matrice Kp elementaire
  9. C
  10. C (D'apres "Design variations of nonlinear elastic structures
  11. C subjected to follower forces"
  12. C M.J. Poldneff, I.S. Rai, J.S. Arora)
  13. C
  14. IMPLICIT INTEGER(I-N)
  15. IMPLICIT REAL*8(A-H,O-Z)
  16. DIMENSION XP(4), XX(3, 8), XKP(48, 48), SKRO(3, 3, 3),
  17. 1 XN(8, 30), XNB1(30), XNB2(30), XNB3(30), XNB4(30), XNB(30),
  18. 2 DNB1(30), DNB2(30), DN(2, 8, 30), DX(2, 3, 30),
  19. 3 XTIN(30), A(30), B(30), C(30), D(30), E(30), F(30),
  20. 4 XNB5(30), XNB6(30), XNB7(30), XNB8(30), XNB9(30), XNB10(30),
  21. 5 XDUM(48, 48)
  22. DATA XNB1/0.5D0, -0.5D0, 0.D0, 0.D0, 0.D0,
  23. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  24. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  25. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  26. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  27. DATA XNB2/0.5D0, 0.5D0, 0.D0, 0.D0, 0.D0,
  28. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  29. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  30. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  31. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  32. DATA XNB3/0.5D0, 0.D0, -0.5D0, 0.D0, 0.D0,
  33. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  34. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  35. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  36. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  37. DATA XNB4/0.5D0, 0.D0, 0.5D0, 0.D0, 0.D0,
  38. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  39. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  40. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  41. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  42. DATA XNB5/1.D0, 0.D0, 0.D0, 0.D0, -1.D0,
  43. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  44. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  45. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  46. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  47. DATA XNB6/1.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  48. 1 -1.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  49. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  50. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  51. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  52. DATA XNB7/-1.D0, -1.D0, -1.D0, 0.D0, 0.D0,
  53. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  54. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  55. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  56. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  57. DATA XNB8/-1, 1.D0, 1.D0, 0.D0, 0.D0, 0.D0,
  58. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  59. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  60. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  61. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  62. DATA XNB9/-1.D0, 1.D0, -1.D0, 0.D0, 0.D0,
  63. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  64. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  65. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  66. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  67. DATA XNB10/-1.D0, -1.D0, 1.D0, 0.D0, 0.D0,
  68. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  69. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  70. 1 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0, 0.D0,
  71. 2 0.D0, 0.D0, 0.D0, 0.D0, 0.D0/
  72. C
  73. C Initialisation du symbole epsilon
  74. C
  75. DO 10 I = 1, 3
  76. DO 10 J = 1, 3
  77. DO 10 K = 1, 3
  78. 10 SKRO(I, J, K) = 0.D0
  79. SKRO(1, 2, 3) = 1.D0
  80. SKRO(1, 3, 2) = -1.D0
  81. SKRO(2, 3, 1) = 1.D0
  82. SKRO(2, 1, 3) = -1.D0
  83. SKRO(3, 1, 2) = 1.D0
  84. SKRO(3, 2, 1) = -1.D0
  85. C
  86. C Calcul de la pression
  87. C
  88. * P = 0.
  89. * DO 20 I = 1, 9
  90. * 20 P = P + XP(I)
  91. * P = P/9.
  92. C
  93. C Fonctions de forme et derivees
  94. C
  95. C Les coefficients sont ranges comme suit :
  96. C indice : 1 2 3 4 5 6 7 8 9
  97. C terme : 1 T1 T2 T1*T2 T1^2 T2^2 T1*T2^2 T1^2*T2 T1^3
  98. C indice : 10 11 12 13 14 15
  99. C terme : T2^3 T1*T2^3 T1^2*T2^2 T1^3*T2^3 T1^4 T2^4
  100. C indice : 16 17 18 19 20 21
  101. C terme : T1*T2^4 T1^2*T2^3 T1^3*T2^2 T1^4*T2 T1^5 T2^5
  102. C indice : 22 23 24 25 26
  103. C terme : T1*T2^5 T1^2*T2^4 T1^3*T2^3 T1^4*T2^2 T1^5*T2
  104. C indice : 27 28 29 30
  105. C terme : T1^2*T2^5 T1^3*T2^4 T1^4*T2^3 T1^5*T2^2
  106. C
  107. CALL MULQP2(XNB1, XNB3, A)
  108. CALL MULQP2(A, XNB7, XNB)
  109. CALL DERQP2(XNB, 1, DNB1)
  110. CALL DERQP2(XNB, 2, DNB2)
  111. DO 31 I = 1, 30
  112. XN(1, I) = XNB(I)
  113. DN(1, 1, I) = DNB1(I)
  114. DN(2, 1, I) = DNB2(I)
  115. 31 CONTINUE
  116. CALL MULQP2(XNB3, XNB5, XNB)
  117. CALL DERQP2(XNB, 1, DNB1)
  118. CALL DERQP2(XNB, 2, DNB2)
  119. DO 32 I = 1, 30
  120. XN(2, I) = XNB(I)
  121. DN(1, 2, I) = DNB1(I)
  122. DN(2, 2, I) = DNB2(I)
  123. 32 CONTINUE
  124. CALL MULQP2(XNB2, XNB3, A)
  125. CALL MULQP2(A, XNB9, XNB)
  126. CALL DERQP2(XNB, 1, DNB1)
  127. CALL DERQP2(XNB, 2, DNB2)
  128. DO 33 I = 1, 30
  129. XN(3, I) = XNB(I)
  130. DN(1, 3, I) = DNB1(I)
  131. DN(2, 3, I) = DNB2(I)
  132. 33 CONTINUE
  133. CALL MULQP2(XNB2, XNB6, XNB)
  134. CALL DERQP2(XNB, 1, DNB1)
  135. CALL DERQP2(XNB, 2, DNB2)
  136. DO 34 I = 1, 30
  137. XN(4, I) = XNB(I)
  138. DN(1, 4, I) = DNB1(I)
  139. DN(2, 4, I) = DNB2(I)
  140. 34 CONTINUE
  141. CALL MULQP2(XNB2, XNB4, A)
  142. CALL MULQP2(A, XNB8, XNB)
  143. CALL DERQP2(XNB, 1, DNB1)
  144. CALL DERQP2(XNB, 2, DNB2)
  145. DO 35 I = 1, 30
  146. XN(5, I) = XNB(I)
  147. DN(1, 5, I) = DNB1(I)
  148. DN(2, 5, I) = DNB2(I)
  149. 35 CONTINUE
  150. CALL MULQP2(XNB4, XNB5, XNB)
  151. CALL DERQP2(XNB, 1, DNB1)
  152. CALL DERQP2(XNB, 2, DNB2)
  153. DO 36 I = 1, 30
  154. XN(6, I) = XNB(I)
  155. DN(1, 6, I) = DNB1(I)
  156. DN(2, 6, I) = DNB2(I)
  157. 36 CONTINUE
  158. CALL MULQP2(XNB1, XNB4, A)
  159. CALL MULQP2(A, XNB10, XNB)
  160. CALL DERQP2(XNB, 1, DNB1)
  161. CALL DERQP2(XNB, 2, DNB2)
  162. DO 37 I = 1, 30
  163. XN(7, I) = XNB(I)
  164. DN(1, 7, I) = DNB1(I)
  165. DN(2, 7, I) = DNB2(I)
  166. 37 CONTINUE
  167. CALL MULQP2(XNB1, XNB6, XNB)
  168. CALL DERQP2(XNB, 1, DNB1)
  169. CALL DERQP2(XNB, 2, DNB2)
  170. DO 38 I = 1, 30
  171. XN(8, I) = XNB(I)
  172. DN(1, 8, I) = DNB1(I)
  173. DN(2, 8, I) = DNB2(I)
  174. 38 CONTINUE
  175. C
  176. C Vecteurs tangents a la coque dans la configuration initiale
  177. C
  178. C Initialisation
  179. DO 40 I = 1, 2
  180. C Boucle sur les parametres
  181. DO 40 J = 1, 3
  182. C Boucle sur les composantes
  183. DO 40 K = 1, 30
  184. C Boucle sur les coefficients des polynomes
  185. DX(I, J, K) = 0.D0
  186. 40 CONTINUE
  187. C Calcul
  188. DO 50 I = 1, 2
  189. C Boucle sur les parametres
  190. DO 50 J = 1, 3
  191. C Boucle sur les composantes
  192. DO 50 K = 1, 30
  193. C Boucle sur les coefficients des polynomes
  194. DO 50 L = 1, 8
  195. C Boucle sur les noeuds
  196. DX(I, J, K) = DX(I, J, K) + XX(J, L)*DN(I, L, K)
  197. 50 CONTINUE
  198. C
  199. C Calcul des termes de la matrice Kp
  200. C
  201. C Initialisation
  202. DO 55 I = 1, 48
  203. DO 55 J = 1, 48
  204. XDUM(I,J) = 0.D0
  205. 55 XKP(I, J) = 0.D0
  206. C Calcul
  207. DO 60 II = 1, 3
  208. DO 60 IL = 1, 8
  209. DO 60 IS = 1, 3
  210. DO 60 IT = 1, 8
  211. XRES = 0.D0
  212. DO 70 J = 1, 3
  213. IFLAG = 0
  214. DO 71 K = 1, 30
  215. A(K) = XN(IL, K)
  216. D(K) = 0.D0
  217. 71 CONTINUE
  218. IF (SKRO(II, J, IS) .NE. 0.D0) THEN
  219. DO 81 K = 1, 30
  220. B(K) = DN(2, IT, K)
  221. C(K) = DX(1, J, K)
  222. 81 CONTINUE
  223. CALL MULQP2(B, C, E)
  224. DO 72 K = 1, 30
  225. 72 D(K) = SKRO(II, J, IS)*E(K)
  226. IFLAG = 1
  227. ENDIF
  228. IF (SKRO(II, IS, J) .NE. 0.D0) THEN
  229. DO 82 K = 1, 30
  230. B(K) = DN(1, IT, K)
  231. C(K) = DX(2, J, K)
  232. 82 CONTINUE
  233. CALL MULQP2(B, C, E)
  234. IF (IFLAG .EQ. 1) THEN
  235. DO 73 K = 1, 30
  236. 73 D(K) = D(K) + SKRO(II, IS, J)*E(K)
  237. ELSE
  238. DO 74 K = 1, 30
  239. 74 D(K) = SKRO(II, IS, J)*E(K)
  240. ENDIF
  241. IFLAG = 1
  242. ENDIF
  243. IF (IFLAG .NE. 0) THEN
  244. CALL MULQP2(A, D, XTIN)
  245. CALL INTQP2(XTIN, TING)
  246. XRES = XRES + TING
  247. ENDIF
  248. 70 CONTINUE
  249. XRES = XRES*P
  250. 60 XDUM(6*(IL-1) + II, 6*(IT-1) + IS) = XRES
  251. IF (IANT .EQ. 0) THEN
  252. DO 90 II = 1, 48
  253. DO 90 IJ = II, 48
  254. XKP(II, IJ) = (XDUM(II, IJ) + XDUM(IJ, II))*0.5D0
  255. 90 XKP(IJ, II) = XKP(II, IJ)
  256. ELSE
  257. DO 91 II = 1, 48
  258. DO 91 IJ = 1, 48
  259. 91 XKP(II, IJ) = XDUM(II, IJ)
  260. ENDIF
  261. RETURN
  262. END
  263.  
  264.  
  265.  
  266.  

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