Télécharger divu.eso

Retour à la liste

Numérotation des lignes :

  1. C DIVU SOURCE CB215821 19/08/20 21:16:43 10287
  2. SUBROUTINE DIVU
  3. C-----------------------------------------------------------------------
  4. C Calcul la divergence d'un champ de vecteur (en général la vitesse)
  5. C discretise dans l'espace de Raviart Thomas de plus petit degre.
  6. C-----------------------------------------------------------------------
  7. C
  8. C On obtient la DIVERGENCE de la vitesse en sommant les débits sur
  9. C chaque élément, chaque débit étant de support la normale sortante.
  10. C C'est pourquoi nous avons besoin du MCHAML d'orientation des normales
  11. C de référence afin de modifier la valeur algébrique du débit.
  12. C
  13. C---------------------------
  14. C Phrase d'appel (GIBIANE) :
  15. C---------------------------
  16. C
  17. C CHP1 = DIVU MMODEL TABLE1 CHP2 CHAM1 (CHP3 (CHP4))
  18. C
  19. C------------------------
  20. C Opérandes et résultat :
  21. C------------------------
  22. C
  23. C CHP1 : CHAMPOINT centre contenant la divergence de la vitesse
  24. C MMODEL : Objet MMODEL spécifiant la formulation
  25. C TABLE1 : Objet table DOMAINE regroupant les infos géometriques
  26. C CHP2 : CHAMPOINT contenant le debit aux faces
  27. C CHAM1 : MCHAML sens de la normale à la face (1=out , -1=in)
  28. C CHP3 : CHAMPOIN contenant une fonction aux centres des elements
  29. C CHP4 : CHAMPOIN contenant les valeurs de cette fonction sur les
  30. C limites du domaine(le support est contitué par des points face)
  31. C
  32. C-----------------------------------------------------------------------
  33. C
  34. C Langage : ESOPE + FORTRAN77
  35. C
  36. C Auteurs : F.DABBENE 01/94
  37. C Complements : C. LE POTIER ET F. AURIOL 20/00
  38. C-----------------------------------------------------------------------
  39. IMPLICIT INTEGER(I-N)
  40. IMPLICIT REAL*8 (A-H,O-Z)
  41. *
  42. -INC CCOPTIO
  43. -INC SMCHAML
  44. -INC SMCHPOI
  45. -INC SMELEME
  46. -INC SMMODEL
  47. -INC SMTABLE
  48. *
  49. SEGMENT IPMAHY
  50. INTEGER MAHYBR(NSOUS)
  51. ENDSEGMENT
  52. *
  53. CHARACTER*4 NOMTOT(1)
  54. CHARACTER*8 CHARIN,LETYPE
  55. *
  56. * Initialisations
  57. *
  58. CHARIN = 'MAILLAGE'
  59. *
  60. * Lecture du MMODEL
  61. *
  62. CALL LIROBJ('MMODEL ',IPMODE,1,IRET1)
  63. CALL ACTOBJ('MMODEL ',IPMODE,1)
  64. IF (IERR.NE.0) RETURN
  65. *
  66. * Lecture de la TABLE domaine
  67. *
  68. IPTABL = 0
  69. CALL LEKMOD(IPMODE,IPTABL,IRET)
  70. IF (IERR.NE.0) RETURN
  71. CALL LEKTAB(IPTABL,CHARIN,IOBRE)
  72. IF (IERR.NE.0) RETURN
  73. IPGEOM = IOBRE
  74. CALL LEKTAB(IPTABL,'ELTFA',IOBRE)
  75. IF (IERR.NE.0) RETURN
  76. IELTFA = IOBRE
  77. CALL LEKTAB(IPTABL,'CENTRE',IOBRE)
  78. IF (IERR.NE.0) RETURN
  79. IPCENT = IOBRE
  80. CALL LEKTAB(IPTABL,'FACE',IOBRE)
  81. IF (IERR.NE.0) RETURN
  82. IPFACE = IOBRE
  83. CALL LEKTAB(IPTABL,'FACEL',IOBRE)
  84. IF (IERR.NE.0) RETURN
  85. IFACEL = IOBRE
  86. *
  87. * Lecture du CHAMPOINT de debit (nom de composantes FLUX)
  88. *
  89. CALL LIROBJ('CHPOINT ',ICHP1,1,IRET1)
  90. CALL ACTOBJ('CHPOINT ',ICHP1,1)
  91. IF (IERR.NE.0) RETURN
  92. *
  93. * Lecture du CHAMELEM
  94. *
  95. CALL LIROBJ('MCHAML ',IPIN,1,IRET1)
  96. CALL ACTOBJ('MCHAML ',IPIN,1)
  97. CALL REDUAF(IPIN,IPMODE,IPCHEL,0,IR,KER)
  98. IF(IR .NE. 1) CALL ERREUR(KER)
  99. IF(IERR .NE. 0) RETURN
  100. C CALL LEKTAB(IPTABL,'XXNORMAE',IPCHEL)
  101. C IF (IERR.NE.0) RETURN
  102. *
  103. * Lecture du CHAMPOINT de la fonction aux centres
  104. *
  105. IPFONC=0
  106. CALL LIROBJ('CHPOINT ',ICHP2,0,IRET2)
  107. IF (IERR.NE.0) RETURN
  108. IF(IRET2.EQ.1)THEN
  109. CALL ACTOBJ('CHPOINT ',ICHP2,1)
  110. INDIC = 1
  111. NBCOMP = 1
  112. NOMTOT(1) = 'SCAL'
  113. CALL QUEPOI(ICHP2,IPCENT,INDIC,NBCOMP,NOMTOT)
  114. IF (IERR.NE.0) RETURN
  115. ICLIM=0
  116. CALL LIROBJ('CHPOINT ',ICLIM,0,IRLIM)
  117. IF(IRLIM .EQ. 1) CALL ACTOBJ('CHPOINT ',ICLIM,1)
  118. CALL DIVUFN(ICHP2,ICLIM,IPFACE,IFACEL,ICHP1,IPCHEL,IPFONC)
  119.  
  120. * VERRUE
  121. * CALL ECCHPO(ipfonc)
  122. ENDIF
  123. *
  124. * Test CHAMPOINT de FLUX
  125. *
  126. INDIC = 1
  127. NBCOMP = 1
  128. NOMTOT(1) = 'FLUX'
  129. CALL QUEPOI(ICHP1,IPFACE,INDIC,NBCOMP,NOMTOT)
  130. IF (IERR.NE.0) RETURN
  131. *
  132. * Test de la formulation contenu dans MMODEL
  133. *
  134. MMODEL = IPMODE
  135. SEGACT MMODEL
  136. NSOUS = KMODEL(/1)
  137. SEGINI IPMAHY
  138. IDARCY = 0
  139. IPT1=IPGEOM
  140. IPT2=IPGEOM
  141. DO 10 ISOUS=1,NSOUS
  142. IF(NSOUS.GT.1)THEN
  143. SEGACT IPT2
  144. IPT1=IPT2.LISOUS(ISOUS)
  145. ENDIF
  146. IMODEL = KMODEL(ISOUS)
  147. SEGACT IMODEL
  148. LETYPE = FORMOD(1)
  149. IF (LETYPE.EQ.'DARCY') THEN
  150. IDARCY = IDARCY + 1
  151. MAHYBR(ISOUS) = IPT1
  152. ENDIF
  153. 10 CONTINUE
  154. IF (IDARCY.EQ.0) THEN
  155. SEGSUP IPMAHY
  156. MOTERR = LETYPE
  157. CALL ERREUR(193)
  158. RETURN
  159. ENDIF
  160. *
  161. * Récuperation des pointeurs ELTFA pour les zones ou DARCY est defini
  162. *
  163. MELEME = IELTFA
  164. SEGACT MELEME
  165. LZONES = LISOUS(/1)
  166. IF (LZONES.EQ.0) LZONES = 1
  167. IPT1 = IPGEOM
  168. SEGACT IPT1
  169. DO 30 ISOUS=1,NSOUS
  170. IMAMEL = MAHYBR(ISOUS)
  171. IF (IMAMEL.NE.0) THEN
  172. DO 20 ISZ=1,LZONES
  173. IF (LZONES.EQ.1) THEN
  174. IPT2 = IPT1
  175. IPT3 = MELEME
  176. ELSE
  177. IPT2 = IPT1.LISOUS(ISZ)
  178. IPT3 = LISOUS(ISZ)
  179. ENDIF
  180. IF (IPT2.EQ.IMAMEL) THEN
  181. MAHYBR(ISOUS) = IPT3
  182. GOTO 30
  183. ENDIF
  184. 20 CONTINUE
  185. IF (IMAMEL.EQ.MAHYBR(ISOUS)) THEN
  186. MOTERR(1:8) = ' MODELE '
  187. MOTERR(1:16) = ' ELTFA '
  188. INTERR(1) = ISOUS
  189. CALL ERREUR(698)
  190. SEGSUP IPMAHY
  191. RETURN
  192. ENDIF
  193. ENDIF
  194. 30 CONTINUE
  195. *
  196. * Test du SUPPORT du MCHAML d'orientation : On teste qu'à chaque zone
  197. * du MODELE ou Darcy est défini correspond une zone du MCHAML
  198. *
  199. MCHELM = IPCHEL
  200. SEGACT MCHELM
  201. NCH1 = IMACHE(/1)
  202. NBOK = 0
  203. DO 50 ISOUS=1,NSOUS
  204. IMAMEL = MAHYBR(ISOUS)
  205. IF (IMAMEL.NE.0) THEN
  206. DO 40 JSOUS=1,NCH1
  207. IF (IMACHE(JSOUS).EQ.IMAMEL) NBOK=NBOK+1
  208. 40 CONTINUE
  209. ENDIF
  210. 50 CONTINUE
  211. SEGSUP IPMAHY
  212. IF (NBOK.NE.NCH1) THEN
  213. MOTERR(1:8) = ' ORIENT '
  214. MOTERR(9:16)= ' ELTFA '
  215. CALL ERREUR(135)
  216. RETURN
  217. ENDIF
  218. *
  219. * Calcul de la divergence du champ de vecteur connu via ses flux
  220. *
  221. CALL DIVU1(ICHP1,IPCHEL,IPFACE,IPCENT,IPFONC,IRET)
  222.  
  223. CALL ACTOBJ('CHPOINT ',IRET,1)
  224. CALL ECROBJ('CHPOINT ',IRET)
  225.  
  226. IF(IPFONC.NE.0)THEN
  227. CALL ACTOBJ('CHPOINT ',IPFONC,1)
  228. CALL ECROBJ('CHPOINT ',IPFONC)
  229. * CALL ECCHPO(ipfonc)
  230. * CALL ECCHPO(iret)
  231. MCHPOI=IPFONC
  232. MSOUPO=IPCHP(1)
  233. MPOVAL=IPOVAL
  234. ENDIF
  235.  
  236. END
  237.  
  238.  
  239.  
  240.  

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