Télécharger inelcu.eso

Retour à la liste

Numérotation des lignes :

inelcu
  1. C INELCU SOURCE GOUNAND 21/06/02 21:16:26 11022
  2. SUBROUTINE INELCU(MYLRFS,IMPR,IRET)
  3. IMPLICIT REAL*8 (A-H,O-Z)
  4. IMPLICIT INTEGER (I-N)
  5. C***********************************************************************
  6. C NOM : INELCU
  7. C PROJET : Noyau linéaire NLIN
  8. C DESCRIPTION : Remplit le segment des éléments de référence
  9. C avec les éléments de référence de dimension 3,
  10. C de forme géométrique cubique.
  11. C
  12. C LANGAGE : ESOPE
  13. C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
  14. C mél : gounand@semt2.smts.cea.fr
  15. C***********************************************************************
  16. C APPELES : INILRF, FILRF, PROLRF, INILAG, PROBAP, GBAPCO
  17. C APPELE PAR : INLRFS
  18. C***********************************************************************
  19. C ENTREES : -
  20. C ENTREES/SORTIES : MYLRFS
  21. C SORTIES : -
  22. C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
  23. C***********************************************************************
  24. C VERSION : v1, 23/03/00, version initiale
  25. C HISTORIQUE : v1, 23/03/00, création
  26. C HISTORIQUE : v2, 10/05/00, modif. du segment ELREF
  27. C HISTORIQUE :
  28. C HISTORIQUE :
  29. C***********************************************************************
  30. C Prière de PRENDRE LE TEMPS de compléter les commentaires
  31. C en cas de modification de ce sous-programme afin de faciliter
  32. C la maintenance !
  33. C***********************************************************************
  34.  
  35. -INC PPARAM
  36. -INC CCOPTIO
  37. -INC TNLIN
  38. *-INC SELREF
  39. POINTEUR MYLRFS.ELREFS
  40. POINTEUR ELCOUR.ELREF
  41. POINTEUR ELPRO1.ELREF
  42. POINTEUR ELPRO2.ELREF
  43. *-INC SPOLYNO
  44. POINTEUR MYBPOL.POLYNS
  45. POINTEUR MYBPO1.POLYNS
  46. POINTEUR MYBPO2.POLYNS
  47. *
  48. INTEGER IMPR,IRET
  49. * Elément de nom : L2D1CU4
  50. REAL*8 ZERO,UNS2,UN
  51. PARAMETER (ZERO=0.D0,UNS2=0.5D0,UN=1.D0)
  52. *
  53. INTEGER INDDL
  54. *
  55. * Executable statements
  56. *
  57. IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans inelcu'
  58. *
  59. * Elément de nom : L2D0CU1
  60. * Sur un cube : élément de Lagrange, fonction L2, approximation
  61. * nodale, espace de référence de dimension 3, 1 noeud, 1 degrés de
  62. * liberté, degré de l'approximation : 0
  63. *
  64. * In INILRF : SEGINI ELCOUR
  65. CALL INILRF('L2D0CU1','CUBE','LAGRANGE','L2',
  66. $ 3,1,1,0,
  67. $ ELCOUR,
  68. $ IMPR,IRET)
  69. IF (IRET.NE.0) GOTO 9999
  70. CALL FILRF('L2D0QU1',MYLRFS,ELPRO1,IMPR,IRET)
  71. IF (IRET.NE.0) GOTO 9999
  72. CALL FILRF('L2D0SE1',MYLRFS,ELPRO2,IMPR,IRET)
  73. IF (IRET.NE.0) GOTO 9999
  74. CALL PROLRF(ELPRO1,ELPRO2,ELCOUR,IMPR,IRET)
  75. IF (IRET.NE.0) GOTO 9999
  76. ELCOUR.NPQUAF(1)=27
  77. ELCOUR.NUMCMP(1)=1
  78. * Initialise la correspondance ddl-noeud+ord.der
  79. CALL INILAG(ELCOUR,IMPR,IRET)
  80. IF (IRET.NE.0) GOTO 9999
  81. SEGACT ELPRO1
  82. MYBPO1=ELPRO1.MBPOLY
  83. SEGDES ELPRO1
  84. SEGACT ELPRO2
  85. MYBPO2=ELPRO2.MBPOLY
  86. SEGDES ELPRO2
  87. * Calcule la base polynômiale produit
  88. CALL PROBAP(MYBPO1,MYBPO2,MYBPOL,IMPR,IRET)
  89. IF (IRET.NE.0) GOTO 9999
  90. ELCOUR.MBPOLY=MYBPOL
  91. SEGDES ELCOUR
  92. MYLRFS.LISEL(**)=ELCOUR
  93. *
  94. * Elément de nom : L2D1CU4
  95. * Sur un cube : élément de Lagrange, fonction L2, approximation
  96. * nodale, espace de référence de dimension 3, 4 noeuds, 4 degrés de
  97. * liberté, degré de l'approximation : 1
  98. *
  99. * In INILRF : SEGINI ELCOUR
  100. CALL INILRF('L2D1CU4','CUBE','LAGRANGE','L2',
  101. $ 3,4,4,1,
  102. $ ELCOUR,
  103. $ IMPR,IRET)
  104. IF (IRET.NE.0) GOTO 9999
  105. ELCOUR.XCONOD(1,1)= UNS2
  106. ELCOUR.XCONOD(2,1)= UNS2
  107. ELCOUR.XCONOD(3,1)=-UNS2
  108. ELCOUR.XCONOD(1,2)=-UNS2
  109. ELCOUR.XCONOD(2,2)= ZERO
  110. ELCOUR.XCONOD(3,2)=-UNS2
  111. ELCOUR.XCONOD(1,3)= ZERO
  112. ELCOUR.XCONOD(2,3)=-UNS2
  113. ELCOUR.XCONOD(3,3)=-UNS2
  114. ELCOUR.XCONOD(1,4)= ZERO
  115. ELCOUR.XCONOD(2,4)=-UNS2
  116. ELCOUR.XCONOD(3,4)= UNS2
  117. ELCOUR.NPQUAF(1)=27
  118. ELCOUR.NUMCMP(1)=1
  119. ELCOUR.NPQUAF(2)=27
  120. ELCOUR.NUMCMP(2)=2
  121. ELCOUR.NPQUAF(3)=27
  122. ELCOUR.NUMCMP(3)=3
  123. ELCOUR.NPQUAF(4)=27
  124. ELCOUR.NUMCMP(4)=4
  125. * Initialise la correspondance ddl-noeud+ord.der
  126. CALL INILAG(ELCOUR,IMPR,IRET)
  127. IF (IRET.NE.0) GOTO 9999
  128. * Génère une base polynômiale complète (dimension 3, degré 1)
  129. CALL GBAPCO(3,1,MYBPOL,IMPR,IRET)
  130. IF (IRET.NE.0) GOTO 9999
  131. ELCOUR.MBPOLY=MYBPOL
  132. SEGDES ELCOUR
  133. MYLRFS.LISEL(**)=ELCOUR
  134. *
  135. * Elément de nom : CRD1CU6
  136. * Sur un cube : élément de Lagrange, fonction continue au centre
  137. * des faces, approximation
  138. * nodale, espace de référence de dimension 3, 6 noeuds, 6 degrés de
  139. * liberté, degré de l'approximation : 1
  140. *
  141. * In INILRF : SEGINI ELCOUR
  142. CALL INILRF('CRD1CU6','CUBE','LAGRANGE','HFAC',
  143. $ 3,6,6,1,
  144. $ ELCOUR,
  145. $ IMPR,IRET)
  146. IF (IRET.NE.0) GOTO 9999
  147. ELCOUR.XCONOD(1,1)=ZERO
  148. ELCOUR.XCONOD(2,1)=-UN
  149. ELCOUR.XCONOD(3,1)=ZERO
  150. ELCOUR.XCONOD(1,2)=UN
  151. ELCOUR.XCONOD(2,2)=ZERO
  152. ELCOUR.XCONOD(3,2)=ZERO
  153. ELCOUR.XCONOD(1,3)=ZERO
  154. ELCOUR.XCONOD(2,3)=UN
  155. ELCOUR.XCONOD(3,3)=ZERO
  156. ELCOUR.XCONOD(1,4)=-UN
  157. ELCOUR.XCONOD(2,4)=ZERO
  158. ELCOUR.XCONOD(3,4)=ZERO
  159. ELCOUR.XCONOD(1,5)=ZERO
  160. ELCOUR.XCONOD(2,5)=ZERO
  161. ELCOUR.XCONOD(3,5)=-UN
  162. ELCOUR.XCONOD(1,6)=ZERO
  163. ELCOUR.XCONOD(2,6)=ZERO
  164. ELCOUR.XCONOD(3,6)=+UN
  165. * Les d.d.l. sont aux noeuds 21 à 26...
  166. DO INDDL=1,6
  167. ELCOUR.NPQUAF(INDDL)=INDDL+20
  168. ELCOUR.NUMCMP(INDDL)=1
  169. ENDDO
  170. * Initialise la correspondance ddl-noeud+ord.der
  171. CALL INILAG(ELCOUR,IMPR,IRET)
  172. IF (IRET.NE.0) GOTO 9999
  173. * Génère une base polynômiale complète (dimension 3, degré 1)
  174. CALL GBAPCO(3,1,MYBPOL,IMPR,IRET)
  175. IF (IRET.NE.0) GOTO 9999
  176. * On rajoute les polynômes spécifiques à crouzeix-raviart quadrilatère
  177. CALL GPOCRQ(3,MYBPOL,IMPR,IRET)
  178. IF (IRET.NE.0) GOTO 9999
  179. ELCOUR.MBPOLY=MYBPOL
  180. SEGDES ELCOUR
  181. MYLRFS.LISEL(**)=ELCOUR
  182. *
  183. * Elément de nom : H1D1CU8
  184. * Sur un cube : élément de Lagrange, fonction H1, approximation
  185. * nodale, espace de référence de dimension 3, 8 noeuds, 8 degrés de
  186. * liberté, degré de l'approximation : 1
  187. *
  188. * In INILRF : SEGINI ELCOUR
  189. CALL INILRF('H1D1CU8','CUBE','LAGRANGE','H1',
  190. $ 3,8,8,1,
  191. $ ELCOUR,
  192. $ IMPR,IRET)
  193. IF (IRET.NE.0) GOTO 9999
  194. CALL FILRF('H1D1QU4',MYLRFS,ELPRO1,IMPR,IRET)
  195. IF (IRET.NE.0) GOTO 9999
  196. CALL FILRF('H1D1SE2',MYLRFS,ELPRO2,IMPR,IRET)
  197. IF (IRET.NE.0) GOTO 9999
  198. CALL PROLRF(ELPRO1,ELPRO2,ELCOUR,IMPR,IRET)
  199. IF (IRET.NE.0) GOTO 9999
  200. * Les d.d.l. sont aux noeuds 1,3,5,7...
  201. DO 201 INDDL=1,4
  202. ELCOUR.NPQUAF(INDDL)=(2*INDDL)-1
  203. ELCOUR.NUMCMP(INDDL)=1
  204. 201 CONTINUE
  205. * ... et 13,15,17,19.
  206. DO 203 INDDL=5,8
  207. ELCOUR.NPQUAF(INDDL)=(2*(INDDL-5))+13
  208. ELCOUR.NUMCMP(INDDL)=1
  209. 203 CONTINUE
  210. * Initialise la correspondance ddl-noeud+ord.der
  211. CALL INILAG(ELCOUR,IMPR,IRET)
  212. IF (IRET.NE.0) GOTO 9999
  213. SEGACT ELPRO1
  214. MYBPO1=ELPRO1.MBPOLY
  215. SEGDES ELPRO1
  216. SEGACT ELPRO2
  217. MYBPO2=ELPRO2.MBPOLY
  218. SEGDES ELPRO2
  219. * Calcule la base polynômiale produit
  220. CALL PROBAP(MYBPO1,MYBPO2,MYBPOL,IMPR,IRET)
  221. IF (IRET.NE.0) GOTO 9999
  222. ELCOUR.MBPOLY=MYBPOL
  223. SEGDES ELCOUR
  224. MYLRFS.LISEL(**)=ELCOUR
  225. *
  226. * Elément de nom : H1D2CU20
  227. * Sur un cube : élément de Lagrange, fonction H1, approximation
  228. * nodale, espace de référence de dimension 3, 20 noeuds, 20 degrés de
  229. * liberté, degré de l'approximation : 2
  230. *
  231. * In INILRF : SEGINI ELCOUR
  232. CALL INILRF('H1D2CU20','CUBE','LAGRANGE','H1',
  233. $ 3,20,20,2,
  234. $ ELCOUR,
  235. $ IMPR,IRET)
  236. IF (IRET.NE.0) GOTO 9999
  237. C Inutile
  238. C ELCOUR.XCONOD(1,1)=ZERO
  239. C ELCOUR.XCONOD(2,1)=ZERO
  240. C ELCOUR.XCONOD(3,1)=ZERO
  241. * Les d.d.l sont aux noeuds 1,...,20
  242. DO IDDL=1,20
  243. ELCOUR.NPQUAF(IDDL)=IDDL
  244. ELCOUR.NUMCMP(IDDL)=1
  245. ENDDO
  246. * Initialise la correspondance ddl-noeud+ord.der
  247. CALL INILAG(ELCOUR,IMPR,IRET)
  248. IF (IRET.NE.0) GOTO 9999
  249. * Pas de base polynômiale (on recopie l'élément de castem)
  250. ELCOUR.MBPOLY=0
  251. SEGDES ELCOUR
  252. MYLRFS.LISEL(**)=ELCOUR
  253. *
  254. * Elément de nom : H1D2CU27
  255. * Sur un cube : élément de Lagrange, fonction H1, approximation
  256. * nodale, espace de référence de dimension 3, 27 noeuds, 27 degrés de
  257. * liberté, degré de l'approximation : 2
  258. *
  259. * In INILRF : SEGINI ELCOUR
  260. CALL INILRF('H1D2CU27','CUBE','LAGRANGE','H1',
  261. $ 3,27,27,2,
  262. $ ELCOUR,
  263. $ IMPR,IRET)
  264. IF (IRET.NE.0) GOTO 9999
  265. CALL FILRF('H1D2QU9',MYLRFS,ELPRO1,IMPR,IRET)
  266. IF (IRET.NE.0) GOTO 9999
  267. CALL FILRF('H1D2SE3',MYLRFS,ELPRO2,IMPR,IRET)
  268. IF (IRET.NE.0) GOTO 9999
  269. CALL PROLRF(ELPRO1,ELPRO2,ELCOUR,IMPR,IRET)
  270. IF (IRET.NE.0) GOTO 9999
  271. * Les ieme d.d.l sont aux noeuds j
  272. DO IDDL=1,8
  273. ELCOUR.NPQUAF(IDDL)=IDDL
  274. ENDDO
  275. ELCOUR.NPQUAF( 9)=25
  276. ELCOUR.NPQUAF(10)=9
  277. ELCOUR.NPQUAF(11)=21
  278. ELCOUR.NPQUAF(12)=10
  279. ELCOUR.NPQUAF(13)=22
  280. ELCOUR.NPQUAF(14)=11
  281. ELCOUR.NPQUAF(15)=23
  282. ELCOUR.NPQUAF(16)=12
  283. ELCOUR.NPQUAF(17)=24
  284. ELCOUR.NPQUAF(18)=27
  285. DO IDDL=19,26
  286. ELCOUR.NPQUAF(IDDL)=IDDL-6
  287. ENDDO
  288. ELCOUR.NPQUAF(27)=26
  289. DO IDDL=1,27
  290. ELCOUR.NUMCMP(IDDL)=1
  291. ENDDO
  292. * Initialise la correspondance ddl-noeud+ord.der
  293. CALL INILAG(ELCOUR,IMPR,IRET)
  294. IF (IRET.NE.0) GOTO 9999
  295. SEGACT ELPRO1
  296. MYBPO1=ELPRO1.MBPOLY
  297. SEGDES ELPRO1
  298. SEGACT ELPRO2
  299. MYBPO2=ELPRO2.MBPOLY
  300. SEGDES ELPRO2
  301. * Calcule la base polynômiale produit
  302. CALL PROBAP(MYBPO1,MYBPO2,MYBPOL,IMPR,IRET)
  303. IF (IRET.NE.0) GOTO 9999
  304. ELCOUR.MBPOLY=MYBPOL
  305. SEGDES ELCOUR
  306. MYLRFS.LISEL(**)=ELCOUR
  307. *
  308. * Normal termination
  309. *
  310. IRET=0
  311. RETURN
  312. *
  313. * Format handling
  314. *
  315. *
  316. * Error handling
  317. *
  318. 9999 CONTINUE
  319. IRET=1
  320. WRITE(IOIMP,*) 'An error was detected in subroutine inelcu'
  321. RETURN
  322. *
  323. * End of subroutine INELCU
  324. *
  325. END
  326.  
  327.  
  328.  
  329.  
  330.  
  331.  

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