Télécharger ef0kon.eso

Retour à la liste

Numérotation des lignes :

  1. C EF0KON SOURCE PV 16/11/17 21:59:16 9180
  2. C**************************************************************
  3. C==============================================================
  4. C Cas de l'equation convection diffusion C
  5. C formulation EFM (elements finis mixte) C
  6. C Auteur : M. MARIN C
  7. C C
  8. C IERRKON : indicatif d'erreur = 1 si il y a une erreur C
  9. C = 0 sinon C
  10. C Cette routine cree les matrices elementaires pour un C
  11. C terme u*teta ou : C
  12. C u : CHPOINT vectoriel au sommet C
  13. C teta : CHPOINT scalaire au centre C
  14. C C
  15. C==============================================================
  16. C**************************************************************
  17.  
  18. SUBROUTINE EF0KON(KIZX,IZTU1,IZTGG1,IZTGG2,IZTGG3,TYPC,IERRKON,
  19. & MELEME,MTABZ,NOMI,IKR,MLENTI,IAXI,NOMII,MZPHI,TYPCFI,
  20. & IZTCO,NELZ,IKU,IKM,AIMPL,IDCEN,MLENT1,DT)
  21. IMPLICIT INTEGER(I-N)
  22. IMPLICIT REAL*8 (A-H,O-Z)
  23.  
  24.  
  25. -INC CCVQUA4
  26. -INC CCOPTIO
  27. -INC CCGEOME
  28. -INC SMCOORD
  29. -INC SMLENTI
  30. -INC SMELEME
  31. POINTEUR MELEMC.MELEME,MELEMS.MELEME
  32. -INC SMCHAML
  33. -INC SMCHPOI
  34. POINTEUR IZTU1.MPOVAL,MZPHI.MPOVAL,IZTCO.MPOVAL
  35. POINTEUR IZTGG1.MPOVAL,IZTGG2.MPOVAL,IZTGG3.MPOVAL
  36. -INC SIZFFB
  37. -INC SMLMOTS
  38. -INC SMTABLE
  39. POINTEUR KIZX.MTABLE,MTABZ.MTABLE
  40.  
  41. C===========================================
  42. C Variables locale ou argument C
  43. CHARACTER*8 TYPC,NOMI,NOM,NOM0,NOMII,TYPCFI
  44. INTEGER IERRKON,IKR,IKU,IKM
  45. C===========================================
  46.  
  47. c WRITE(6,*) KIZX,IZTU1,IZTGG1,IZTGG2,IZTGG3,TYPC,IERRKON,
  48. c & MELEM,MTABZ,NOMI,IKR,MLENTI,IAXI,NOMII,MZPHI,TYPCFI,
  49. c & IZTCO,NELZ,IKU,IKM,AIMPL
  50. C Validation du passage d'argument
  51.  
  52. IERRKON=0
  53.  
  54. C Initialisation pour le nombre total final d'elements (TRI3+QUA4+...)
  55. C sur le domaine ou s'applique KONV
  56. NUTOEL=0
  57. c WRITE(6,*) 'MELEME=',MELEME
  58.  
  59. C WRITE(6,*) 'Type du Coef 1:',IKR
  60.  
  61. C On teste si la variable primale est de type CENTRE
  62. IF (TYPC.NE.'CENTRE ') THEN
  63. CALL ERRKON(1,IERRKON)
  64. RETURN
  65. END IF
  66.  
  67. C Test si la variable primale est de type scalaire
  68. IF (IZTU1.VPOCHA(/2).NE.1) THEN
  69. CALL ERRKON(5,IERRKON)
  70. RETURN
  71. END IF
  72.  
  73. C On teste si la variable duale est de type SOMMET
  74. IF (TYPCFI.NE.'SOMMET ') THEN
  75. CALL ERRKON(4,IERRKON)
  76. RETURN
  77. END IF
  78.  
  79. C On teste si la variable duale est un vecteur
  80. IF ((MZPHI.VPOCHA(/2).EQ.1).AND.(IDIM.NE.1)) THEN
  81. CALL ERRKON(6,IERRKON)
  82. RETURN
  83. END IF
  84.  
  85. C On teste pour savoir si RO est bien un scalaire
  86. IF (IKR.NE.1) THEN
  87. CALL ERRKON(7,IERRKON)
  88. RETURN
  89. END IF
  90. C On recupere le nombre de composante de la vitesse
  91. NBME=IDIM
  92. c WRITE(6,*) 'Nombre de composante de UN :',NBME
  93.  
  94. C On active le segment contenant le maillage de la zone ou
  95. C KONV s'applique et on initialise le nombre de sous-objet.
  96. C ------------------------------------------------------
  97. C NBSOUS et NBME servent a initialiser le segment IMATRI
  98. SEGACT MELEME
  99. NBSOUS=MELEME.LISOUS(/1)
  100. c WRITE(6,*) 'NBSOUS=',NBSOUS
  101. IF (NBSOUS.EQ.0) NBSOUS=1
  102.  
  103. C Initialisation de la matrice elementaire de l'operateur
  104. NRIGE=7
  105. NKID=9
  106. NKMT=7
  107. NMATRI=1
  108. SEGINI MATRIK
  109.  
  110.  
  111. C On recupere le maillage des centres
  112. CALL LEKTAB(MTABZ,'CENTRE',MELEMC)
  113.  
  114. C On recupere le maillage des sommets
  115. CALL LEKTAB(MTABZ,'SOMMET',MELEMS)
  116.  
  117. C Inconnue primale
  118. IRIGEL(1,1)=MELEMC
  119.  
  120. C Inconnue duale
  121. IRIGEL(2,1)=MELEME
  122.  
  123. C Matrice rectangulaire
  124. IRIGEL(7,1)=3
  125.  
  126. SEGINI IMATRI
  127. IRIGEL(4,1)= IMATRI
  128.  
  129. C Intialisation des supports inconnue primale et duale
  130. C pour le segment IMATRI
  131. KSPGP=MELEMC
  132. KSPGD=MELEMS
  133.  
  134. DO I=1,NBME
  135. WRITE(NOM,FMT='(I1,A3)')I,NOMII(1:3)
  136. LISPRI(I)=NOMI(1:4)//' '
  137. LISDUA(I)=NOM(1:4)//' '
  138. END DO
  139.  
  140. C On recupere le segment qui contient int(Ni)
  141. CALL LEKTAB(MTABZ,'XXPSOML',MCHELM)
  142.  
  143. IF (MCHELM.EQ.0) THEN
  144. CALL ERRKON(2,IERRKON)
  145. RETURN
  146. END IF
  147.  
  148. SEGACT MCHELM
  149.  
  150. C==========Boucle principale==============C
  151. C On effectue cette boucle pour tous les
  152. C sous-objets
  153.  
  154.  
  155. DO L=1,NBSOUS
  156. IPT1=MELEME
  157. SEGACT MELEME
  158.  
  159. C Si NBSOUS > 1 on agit sur les sous-objets LISOUS
  160. IF (NBSOUS.NE.1) IPT1=LISOUS(L)
  161. SEGACT IPT1
  162.  
  163. c WRITE(6,*) 'IPT1=',IPT1,'NBSOUS=',NBSOUS
  164. C On rempli NOM0 du veritable nom des elements traites (TRI3, QUA4,...)
  165. NOM0=NOMS(IPT1.ITYPEL)//' '
  166.  
  167. C Cette routine cree des objets contenant les fonctions
  168. C de forme lie aux elements traites
  169. C Le segment IZFFM contient les fonctions de forme et les gradients
  170. C a chaque point de Gauss
  171. CALL KALPBG(NOM0,'FONFORM ',IZFFM)
  172.  
  173. SEGACT IZFFM
  174. IZHR=KZHR(1)
  175. SEGACT IZHR*MOD
  176.  
  177. C Dimension de l'element de reference
  178. NES=GR(/1)
  179.  
  180. C Nombre de points de Gauss
  181. NPG=GR(/3)
  182. c WRITE(6,*) 'Nbre de points de Gauss :' , NPG
  183.  
  184. C NP est le nombre de points de l'element de reference
  185. C MP est la troisieme dimension du tableau contenant les matrices
  186. C elementaire de chaque elements (segment IZAFM). C'est le
  187. C nombre de point par element
  188. NP = 1
  189. MP = IPT1.NUM(/1)
  190. C NBEL est le nombre d'element de type NOM0
  191. NBEL=IPT1.NUM(/2)
  192.  
  193. C On declare un segment IZAFM par composante
  194. SEGINI IPM1
  195. LIZAFM(L,1)=IPM1
  196.  
  197. IF (IDIM.EQ.2) THEN
  198. SEGINI IPM2
  199. LIZAFM(L,2)=IPM2
  200. ELSE IF (IDIM.EQ.3) THEN
  201. SEGINI IPM2
  202. SEGINI IPM3
  203. LIZAFM(L,2)=IPM2
  204. LIZAFM(L,3)=IPM3
  205. END IF
  206.  
  207. C Nombre total de points NTP
  208. NPT=IZTGG2.VPOCHA(/1)
  209.  
  210. C On recupere le tableau contenant les integrales
  211. C de fonction test
  212. MCHAML=ICHAML(L)
  213. SEGACT MCHAML
  214. MELVAL=IELVAL(1)
  215. SEGACT MELVAL
  216. SEGACT MELEMC
  217. C ******
  218. CALL XCNEF0(NBEL,NP,MP,IPM1.AM,IZTGG1.VPOCHA,IZTGG2
  219. & .VPOCHA,IZTGG3.VPOCHA,NPT,IDIM,IDCEN,XYZ,
  220. & NUTOEL,XCOOR,IPT1.NUM,MLENTI.LECT,IPM2.AM,IPM3.AM,
  221. & FN,GR,PG,HR,PGSQ,RPG,NES,NPG,IAXI,VELCHE,
  222. & NBME,IZTGG3.VPOCHA,IZTCO.VPOCHA,NELZ,
  223. & IKR,IKU,IKM,IZTU1.VPOCHA,AIMPL,MLENT1.LECT,DT,
  224. & MELEMC)
  225. C ******
  226. SEGDES MELEMC
  227. SEGDES IPT1,IPM1,IPM2
  228. IF (IDIM.EQ.3) THEN
  229. SEGDES IPM3
  230. END IF
  231. SEGDES MCHAML,MELVAL
  232. SEGSUP IZHR,IZFFM
  233.  
  234.  
  235. C On met a jour le nombre d'element total (quelquesoit son type)
  236. NUTOEL=NUTOEL+NBEL
  237. END DO
  238.  
  239. C==========Fin Boucle principale==========C
  240.  
  241. C On desactive les segments
  242.  
  243. SEGDES MELEME
  244. SEGDES MCHELM
  245.  
  246. C Si IKR = 0 alors IZTGG1 (coef 1) est un CHPOINT
  247. C que l'on desactive.
  248. IF(IKR.EQ.0)THEN
  249. SEGDES IZTGG1
  250. ENDIF
  251.  
  252. C On ecrit dans la table KONV la matrice elementaire
  253. C de l'operateur
  254. CALL ECMO(KIZX,'MATELM','MATRIK',MATRIK)
  255.  
  256. SEGDES IMATRI,MATRIK
  257. RETURN
  258.  
  259. END
  260.  
  261.  
  262.  
  263.  
  264.  
  265.  
  266.  
  267.  
  268.  
  269.  
  270.  
  271.  
  272.  
  273.  

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