Télécharger quahm2.eso

Retour à la liste

Numérotation des lignes :

  1. C QUAHM2 SOURCE CHAT 05/01/13 02:40:56 5004
  2. C QUAHM2 SOURCE AM1 95/11/24 22:54:38 1918
  3. SUBROUTINE QUAHM2(IGAU,ITEL,MFR,NBNO,XEL,SHPTOT,SHP,IFOU,
  4. # NHARM,VKL12,VKL23,VKL33,POIGAU,ISDJC,LRE,REL,IRET)
  5. C=======================================================================
  6. C
  7. C CALCULE LES TERMES EN P * PI ,PI * (UR,RT) ,(UR,RT) *(UR,RT)
  8. C (UT,RR) * (UT,RR) , PI * (UT,RR) DE LA MATRICE
  9. C MASSE DANS LE CAS AXISYMETRIQUE OU FOURIER POUR
  10. C LA FORMULATION (37) HOMOGENE
  11. C=======================================================================
  12. C INPUT
  13. C IGAU=NUMERO DU POINT DE GAUSS
  14. C ITEL=NUMERO DE L ELEMENT DANS NOMTP
  15. C MFR =NUMERO DE LA FORMULATION
  16. C NBNO=NOMBRE DE NOEUDS
  17. C XEL =COORDONNEES DE L ELEMENT
  18. C IFOU=IFOUR DE CCOPTIO
  19. C NHARM=NUMERO DU MODE DE FOURIER
  20. C VKL12=-((COEFPI*COEFPR)/(RHOF*C**2))*SFLU/SCEL
  21. C VKL23=(BET11+BET22)*COEFPI/(2.*SCEL)
  22. C VKL33=(RHOS*2.+RHOF*(BET11+BET22))/SCEL
  23. C POIGAU=MINTE.POIGAU(IGAU)
  24. C LRE =NOMBRE DE D.D.L DE LA MATRICE DE RIGIDITE
  25. C SHPTOT(6,NBNO,NBGAU)=FONCTIONS DE FORMES ET DERIVEES
  26. C ISDJC = INDICATEUR SUR LE SIGNE DU JACOBIEN
  27. C ZONE DE TRAVAIL
  28. C SHP(5,NBNO)=TABLEAU DE TRAVAIL
  29. C OUTPUT
  30. C ISDJC = INDICATEUR SUR LE SIGNE DU JACOBIEN
  31. C REL=MATRICE DE MASSE
  32. C IRET:INDICATEUR = 1 : SUCCES
  33. C 0 : ECHEC (ELEMENT MELE INCOMPATIBLE )
  34. C 2 : ECHEC (JACOBIEN NUL )
  35. C 3 :ECHEC (ROUTINE N EST VALABLE QU
  36. C EN AXISYMETRIQUE OU FOURIER )
  37. C=======================================================================
  38. IMPLICIT INTEGER(I-N)
  39. IMPLICIT REAL*8(A-H,O-Z)
  40. DIMENSION XEL(3,*),SHP(6,*),SHPTOT(6,NBNO,*),REL(LRE,*)
  41. IF (ITEL.EQ.126) GOTO 10
  42. C
  43. C ERREUR : TYPE D' ELEMENT INCOMPATIBLE AVEC LA FORMULATION
  44. C
  45. IRET = 0
  46. GOTO 666
  47. 10 CONTINUE
  48. IF (IFOU.EQ.0.OR.IFOU.EQ.1) GOTO 11
  49. C
  50. C MESSAGE D ERREUR : ROUTINE N EST VALABLE QU EN FOURIER
  51. C OU EN AXISYMETRIQUE
  52. C
  53. IRET = 3
  54. GOTO 666
  55. 11 CONTINUE
  56. C
  57. C ELEMENTS HOMOGENEISES QUAH EN AXISYMETRIE OU EN FOURIER
  58. C NBDL = LRE/NBNO NOMBRE DE D.D.L PAR NOEUD
  59. C
  60. NBDL = LRE /NBNO
  61. C
  62. C SHP(1,I) : FONCTION DE FORME
  63. C SHP(2,I) : DERIVEE % R DE LA FONCTION DE FORME
  64. C SHP(3,I) : DERIVEE % Z DE LA FONCTION DE FORME
  65. C
  66. DO 101 NP=1,NBNO
  67. SHP(1,NP)=SHPTOT(1,NP,IGAU)
  68. SHP(2,NP)=SHPTOT(2,NP,IGAU)
  69. SHP(3,NP)=SHPTOT(3,NP,IGAU)
  70. 101 CONTINUE
  71. C
  72. CALL DEVOLU(XEL,SHP,MFR,NBNO,IFOU,NHARM,2,1.D0,RR,DJAC)
  73. IF (DJAC.EQ.0.) GOTO 667
  74. IF (DJAC.LT.0.) ISDJC = ISDJC + 1
  75. C
  76. C FONCTIONS DE FORME
  77. C
  78. SHP(4,1)=SHP(1,1) + SHP(1,4)
  79. SHP(4,2)=SHP(1,2) + SHP(1,3)
  80. SHP(4,3)=0.D0
  81. SHP(4,4)=0.D0
  82. C
  83. XZH=SHP(1,3) + SHP(1,4)
  84. C
  85. A1=MIN(XEL(2,1),XEL(2,2))
  86. A2=MIN(XEL(2,1),XEL(2,3))
  87. A3=MIN(XEL(2,1),XEL(2,4))
  88. A4=MIN(A1,A2)
  89. A5=MIN(A1,A3)
  90. Z1=MIN(A4,A5)
  91. C
  92. B1=MAX(XEL(2,1),XEL(2,2))
  93. B2=MAX(XEL(2,1),XEL(2,3))
  94. B3=MAX(XEL(2,1),XEL(2,4))
  95. B4=MAX(B1,B2)
  96. B5=MAX(B1,B3)
  97. Z2=MAX(B4,B5)
  98. C
  99. RH=Z2-Z1
  100. C
  101. S1=SHP(1,1) + SHP(1,4)
  102. S2=SHP(1,2) + SHP(1,3)
  103. C
  104. C FONCTIONS DE FORME EN Z
  105. C
  106. XZH= SHP(1,3) + SHP(1,4)
  107. C
  108. SHP(5,1)=1.D0-(3.D0*XZH*XZH)+(2.D0*XZH*XZH*XZH)
  109. SHP(5,2)=3.D0*(XZH**2)-(2.D0*XZH*XZH*XZH)
  110. SHP(5,3)=(RH*XZH)*(1.D0-2.D0*XZH+XZH*XZH)
  111. SHP(5,4)=(RH*XZH)*(XZH*XZH-XZH)
  112. C
  113. C FONCTIONS DE FORME POUR LA FLEXION
  114. C
  115. SHP(6,1)=S1*SHP(5,1)
  116. SHP(6,2)=S2*SHP(5,1)
  117. SHP(6,3)=S2*SHP(5,2)
  118. SHP(6,4)=S1*SHP(5,2)
  119. C
  120. SHP(6,5)=S1*SHP(5,3)
  121. SHP(6,6)=S2*SHP(5,3)
  122. SHP(6,7)=S2*SHP(5,4)
  123. SHP(6,8)=S1*SHP(5,4)
  124. C
  125. C
  126. C1=MIN(XEL(1,1),XEL(1,2))
  127. C2=MIN(XEL(1,1),XEL(1,3))
  128. C3=MIN(XEL(1,1),XEL(1,4))
  129. C4=MIN(C1,C2)
  130. C5=MIN(C1,C3)
  131. R1=MIN(C4,C5)
  132. C
  133. D1=MAX(XEL(1,1),XEL(1,2))
  134. D2=MAX(XEL(1,1),XEL(1,3))
  135. D3=MAX(XEL(1,1),XEL(1,4))
  136. D4=MAX(D1,D2)
  137. D5=MAX(D1,D3)
  138. R2=MAX(D4,D5)
  139. C
  140. DELTAR = R2 - R1
  141. DV = (DELTAR*RH)/4.d0
  142. C
  143. DJAC = ABS(DJAC)
  144. C
  145. C
  146. C TERMES EN P * PI
  147. C
  148. DJAC1 = ABS(DJAC)*POIGAU
  149. IX1=0
  150. IY1=0
  151. DO 102 IX=2,LRE ,NBDL
  152. IX1=IX1 + 1
  153. DO 103 IY=1,IX ,NBDL
  154. IY1=IY1 + 1
  155. REL(IY,IX) = REL(IY,IX) + VKL12*DJAC1*SHP(1,IX1)*SHP(1,IY1)
  156. REL(IX,IY) = REL(IY,IX)
  157. 103 CONTINUE
  158. IY1=0
  159. 102 CONTINUE
  160. DO 104 IX=2+NBDL,LRE ,NBDL
  161. IX2=IX - NBDL
  162. DO 105 IY=1,IX2 ,NBDL
  163. REL(IY+1,IX-1) = REL(IY,IX)
  164. REL(IX-1,IY+1) = REL(IY+1,IX-1)
  165. 105 CONTINUE
  166. 104 CONTINUE
  167. C
  168. C TERMES EN PI * (UR , RT )
  169. C
  170. C
  171. IX1=0
  172. IY1=0
  173. DO 106 IX=3,LRE ,NBDL
  174. IX1=IX1 + 1
  175. DO 107 IY=2,IX ,NBDL
  176. IY1=IY1 + 1
  177. REL(IY,IX) = REL(IY,IX) + VKL23*DJAC1*SHP(2,IY1)
  178. #*SHP(6,IX1)
  179. REL(IY,IX+1) = REL(IY,IX+1) + VKL23*DJAC1*SHP(2,IY1)
  180. #*SHP(6,IX1+4)
  181. REL(IX,IY) = REL(IY,IX)
  182. REL(IX+1,IY) = REL(IY,IX+1)
  183. 107 CONTINUE
  184. IY1=0
  185. 106 CONTINUE
  186. IX1=1
  187. IY1=0
  188. DO 108 IX=2+NBDL,LRE ,NBDL
  189. IX1=IX1 + 1
  190. DO 109 IY=3,IX ,NBDL
  191. IY1=IY1 + 1
  192. REL(IY,IX) = REL(IY,IX) + VKL23*DJAC1*SHP(2,IX1)
  193. #*SHP(6,IY1)
  194. REL(IY+1,IX) = REL(IY+1,IX) + VKL23*DJAC1*SHP(2,IX1)
  195. #*SHP(6,IY1+4)
  196. REL(IX,IY) = REL(IY,IX)
  197. REL(IX,IY+1) = REL(IY+1,IX)
  198. 109 CONTINUE
  199. IY1=0
  200. 108 CONTINUE
  201. C
  202. IF ( IFOU.EQ.1) THEN
  203. C
  204. C TERMES EN PI * (UT , RR )
  205. C NON NULS QU EN FOURIER
  206. C
  207. DJAC2 = ABS(DJAC)*POIGAU
  208. C
  209. VKL25 = -1.D0* VKL23*NHARM / RR
  210. IX1=0
  211. IY1=0
  212. DO 110 IX=5,LRE ,NBDL
  213. IX1=IX1 + 1
  214. DO 111 IY=2,IX ,NBDL
  215. IY1=IY1 + 1
  216. REL(IY,IX) = REL(IY,IX) + VKL25*DJAC2*SHP(1,IY1)*SHP(6,IX1)
  217. REL(IY,IX+1) = REL(IY,IX+1) + VKL25*DJAC2*SHP(1,IY1)
  218. #*SHP(6,IX1+4)*(-1.D0)
  219. REL(IX,IY) = REL(IY,IX)
  220. REL(IX+1,IY) = REL(IY,IX+1)
  221. 111 CONTINUE
  222. IY1=0
  223. 110 CONTINUE
  224. IX1=1
  225. IY1=0
  226. DO 112 IX=2+NBDL,LRE ,NBDL
  227. IX1=IX1 + 1
  228. DO 113 IY=5,IX ,NBDL
  229. IY1=IY1 + 1
  230. REL(IY,IX) = REL(IY,IX) + VKL25*DJAC2*SHP(1,IX1)*SHP(6,IY1)
  231. REL(IY+1,IX) = REL(IY+1,IX) + VKL25*DJAC2*SHP(1,IX1)
  232. #*SHP(6,IY1+4)*(-1.D0)
  233. REL(IX,IY) = REL(IY,IX)
  234. REL(IX,IY+1) = REL(IY+1,IX)
  235. 113 CONTINUE
  236. IY1=0
  237. 112 CONTINUE
  238. ENDIF
  239. C
  240. C TERMES EN (UR,RT ) * (UR , RT )
  241. C
  242. C
  243. IX1=0
  244. IY1=0
  245. DO 114 IX=3,LRE ,NBDL
  246. IX1=IX1 + 1
  247. DO 115 IY=3,IX ,NBDL
  248. IY1=IY1 + 1
  249. C
  250. REL(IY,IX) = REL(IY,IX) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1)
  251. REL(IY,IX+1) = REL(IY,IX+1) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1+4)
  252. REL(IY+1,IX) = REL(IY+1,IX) + VKL33*DJAC1*SHP(6,IY1+4)*SHP(6,IX1)
  253. REL(IY+1,IX+1) = REL(IY+1,IX+1)+VKL33*DJAC1*SHP(6,IY1+4)
  254. #*SHP(6,IX1+4)
  255. REL(IX,IY) = REL(IY,IX)
  256. REL(IX+1,IY) = REL(IY,IX+1)
  257. REL(IX,IY+1) = REL(IY+1,IX)
  258. REL(IX+1,IY+1) = REL(IY+1,IX+1)
  259. 115 CONTINUE
  260. IY1=0
  261. 114 CONTINUE
  262. C
  263. IF ( IFOU.EQ.1) THEN
  264. C
  265. C TERMES EN (UT,RR ) * (UT , RR )
  266. C NON NULS QU EN FOURIER
  267. C
  268. IX1=0
  269. IY1=0
  270. DO 303 IX=5,LRE ,NBDL
  271. IX1=IX1 + 1
  272. DO 403 IY=5,IX ,NBDL
  273. IY1=IY1 + 1
  274. REL(IY,IX) = REL(IY,IX) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1)
  275. REL(IY,IX+1) = REL(IY,IX+1) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1+4)
  276. #*(-1.D0)
  277. REL(IY+1,IX) = REL(IY+1,IX) + VKL33*DJAC1*SHP(6,IY1+4)*SHP(6,IX1)
  278. #*(-1.D0)
  279. REL(IY+1,IX+1) = REL(IY+1,IX+1)+VKL33*DJAC1*SHP(6,IY1+4)
  280. #*SHP(6,IX1+4)
  281. REL(IX,IY) = REL(IY,IX)
  282. REL(IX+1,IY) = REL(IY,IX+1)
  283. REL(IX,IY+1) = REL(IY+1,IX)
  284. REL(IX+1,IY+1) = REL(IY+1,IX+1)
  285. 403 CONTINUE
  286. IY1=0
  287. 303 CONTINUE
  288. C
  289. ENDIF
  290. IRET = 1
  291. GOTO 666
  292. C
  293. C MESSAGE D ERREUR : ELEMENT A SURFACE NULLE
  294. C
  295. 667 CONTINUE
  296. IRET = 2
  297. GOTO 666
  298. C
  299. 666 CONTINUE
  300. RETURN
  301. END
  302.  
  303.  
  304.  
  305.  
  306.  

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