Télécharger la plaquette Cast3M


fr en

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