1 : C HDEBI2 SOURCE CHAT 09/10/09 21:18:44 6519 2 : C HDEBI1 SOURCE CHAT 97/12/23 22:20:50 3021 3 : SUBROUTINE HDEBI2(IPMAHY,IPFACE,IPDARC,ICHP0,IPCHEL, 4 : S ICHP1,ICHP2,ITTH,INORM,IORIE,ISURF,IPFORC,IFORC,IRET) 5 : C----------------------------------------------------------------------- 6 : C Calcul le débit aux faces et orientation suivant le sens de la normale 7 : C lorsqu'on connait la concentrations au centre 8 : C----------------------------------------------------------------------- 9 : C 10 : C--------------------------- 11 : C Parametres Entree/Sortie : 12 : C--------------------------- 13 : C 14 : C E/ IPMAHY : MELEME des connectivités éléments/faces pour Darcy 15 : C E/ IPFACE : MELEME des points FACEs -1 16 : C E/ IPDARC : RIGIDITE de sous type DARCY (contient RE ). 17 : C E/ ICHP0 : CHPO face des traces de concentration au temps n 18 : C E/ IPCHEL : MCHAML des orientations de normale (1=out,-1=in) 19 : C E/ ICHP1 : CHPO concentration au centre au temps n 20 : C 21 : C Parametre optionnel 22 : C E/ ICHP2 : CHPO face flux de vitesse 23 : C 24 : C /S IRET : CHPO face des flux les noms des composantes sont 25 : C ceux des composantes de ITPN et IPCH1. 26 : C Si ICHP2 est donné on ajoute le flux convectif 27 : C 28 : C---------------------- 29 : C Variables en COMMON : 30 : C---------------------- 31 : C 32 : C E/ IFOUR : cf CCOPTIO 33 : C E/ IDIM : cf CCOPTIO 34 : C 35 : C---------------------- 36 : C Tableaux de travail : 37 : C---------------------- 38 : C 39 : C ICPR(I)=J : Le noeud I a le numero J dans le MELEME des faces 40 : C Correspondance numérotation globale/locale 41 : C ITES : Nombre de noeuds FACE 42 : C NNGOT : Nombre de noeuds total du domaine 43 : C IVAA(I) : indice du CHAMPOINT au Ieme noeud global 44 : C 45 : C----------------------------------------------------------------------- 46 : C 47 : C Langage : ESOPE + FORTRAN77 48 : C 49 : C Auteurs : 09/93 F.DABBENE - Cas permanent 50 : C 09/94 X.NOUVELLON - Extension au cas transitoire 51 : C 02/96 L.V.BENET - Prise en compte de forces de volume 52 : C 05/98 F.AURIOL expression en fonction des concentrations 53 : C (charges) et traces de concentrations (traces 54 : C de concentrations) possibilités de champs 55 : C à plusieurs composantes, en transitoire. 56 : C 57 : C----------------------------------------------------------------------- 58 : IMPLICIT INTEGER(I-N) 59 : IMPLICIT REAL*8 (A-H,O-Z) 60 : * 61 : -INC CCOPTIO 62 : -INC SMELEME 63 : -INC SMCHPOI 64 : POINTEUR MCHPO5.MCHPOI,MCHPO6.MCHPOI,MSOUP6.MSOUPO 65 : -INC SMCHAML 66 : -INC SMRIGID 67 : -INC SMCOORD 68 : * 69 : CHARACTER*4 MOREFD,MOREFP 70 : SEGMENT IPMAHY 71 : INTEGER MAHYBR(NSOUS) 72 : ENDSEGMENT 73 : SEGMENT ICCPR 74 : INTEGER ICPR(NNGOT) 75 : ENDSEGMENT 76 : SEGMENT IORGA 77 : INTEGER IVAA(ITES), IVBB(ITES) 78 : ENDSEGMENT 79 : SEGMENT ITRAV 80 : REAL*8 RLIGNE(NBDDL) 81 : REAL*8 FLFOR(NBDDL),RFOR(NBDDL) 82 : ENDSEGMENT 83 : C 84 : C 85 : *' Initialisations' 86 : C 87 : MRIGID = IPDARC 88 : RI1 = IPFORC 89 : MCHELM = IPCHEL 90 : IPT1 = IPFACE 91 : MCHPO1 = ICHP0 92 : MCHPO2 = ICHP1 93 : MCHPO3 = ICHP2 94 : IRET = 0 95 : MCHEL2 = IORIE 96 : MCHPO4 = IFORC 97 : MCHPO5 = INORM 98 : MCHPO6 = ISURF 99 : C 100 : SEGACT IPMAHY 101 : NBMAIL = MAHYBR(/1) 102 : SEGACT MRIGID 103 : * 104 : *' Creation du tableau ICPR pour le maillage IPT1' 105 : * 106 : NNGOT = XCOOR(/1)/(IDIM+1) 107 : SEGINI ICCPR 108 : SEGACT IPT1 109 : N2 = IPT1.NUM(/2) 110 : IK = 0 111 : DO 15 I2=1,N2 112 : K = IPT1.NUM(1,I2) 113 : IF (ICPR(K).EQ.0) THEN 114 : IK = IK + 1 115 : ICPR(K) = IK 116 : ENDIF 117 : 15 CONTINUE 118 : SEGDES IPT1 119 : ITES = IK 120 : * 121 : *' Activation du MPOVAL du CHPO de traces de concentrations 122 : * 123 : SEGACT MCHPO1 124 : MSOUP1 = MCHPO1.IPCHP(1) 125 : SEGDES MCHPO1 126 : SEGACT MSOUP1 127 : MPOVA1 = MSOUP1.IPOVAL 128 : SEGACT MPOVA1 129 : NBCOMP= MPOVA1.VPOCHA(/2) 130 : * 131 : * Activation du MPOVAL du CHPO des concentrations au centre 132 : * 133 : SEGACT MCHPO2 134 : MSOUP2 = MCHPO2.IPCHP(1) 135 : SEGDES MCHPO2 136 : SEGACT MSOUP2 137 : MPOVA2 = MSOUP2.IPOVAL 138 : SEGACT MPOVA2 139 : * 140 : * Activation du MPOVAL du CHPO flux de vitesses aux faces 141 : * 142 : IF (ICHP2.NE.0) THEN 143 : SEGACT MCHPO3 144 : MSOUP3 = MCHPO3.IPCHP(1) 145 : SEGDES MCHPO3 146 : SEGACT MSOUP3 147 : MPOVA3 = MSOUP3.IPOVAL 148 : SEGACT MPOVA3 149 : ENDIF 150 : * 151 : * activation des objets liés à la présence d'une force volumique 152 : * 153 : IF (IFORC.NE.0) THEN 154 : * 155 : * Activation du MPOVAL du CHPO force appuyé au centre des éléments volumiques 156 : * 157 : SEGACT MCHPO4 158 : MSOUP4 = MCHPO4.IPCHP(1) 159 : SEGDES MCHPO4 160 : SEGACT MSOUP4 161 : MPOVA4 = MSOUP4.IPOVAL 162 : SEGDES MSOUP4 163 : SEGACT MPOVA4 164 : * 165 : * Activation du MPOVAL du CHPO des vecteurs normales appuyé au centre des faces 166 : * 167 : SEGACT MCHPO5 168 : MSOUP5 = MCHPO5.IPCHP(1) 169 : SEGDES MCHPO5 170 : SEGACT MSOUP5 171 : MPOVA5 = MSOUP5.IPOVAL 172 : SEGDES MSOUP5 173 : SEGACT MPOVA5 174 : * 175 : * Activation du MPOVAL du CHPO des surfaces appuyé au centre des faces 176 : * 177 : SEGACT MCHPO6 178 : MSOUP6 = MCHPO6.IPCHP(1) 179 : SEGDES MCHPO6 180 : SEGACT MSOUP6 181 : MPOVA6 = MSOUP6.IPOVAL 182 : SEGDES MSOUP6 183 : SEGACT MPOVA6 184 : * 185 : * Activation du MCHEL des orientations des faces 186 : * 187 : SEGACT MCHEL2 188 : * 189 : * Activation du MRIGI de la matrice masse hybride 190 : * 191 : SEGACT RI1 192 : ENDIF 193 : * 194 : * On recherche l ordre des traces de concentrations par rapport à IPT1 195 : * 196 : SEGINI IORGA 197 : MELEME = MSOUP1.IGEOC 198 : SEGACT MELEME 199 : N2 = NUM(/2) 200 : DO 25 I2=1,N2 201 : K = NUM(1,I2) 202 : IF (ICPR(K).EQ.0) THEN 203 : INTERR(1) = K 204 : MOTERR(1:8) = 'FACE ' 204 : CALL ERREUR(64) 205 : SEGDES MELEME, MSOUP1 206 : SEGDES MCHPO2, MRIGID, IPMAHY 207 : SEGSUP ICCPR, IORGA 208 : RETURN 209 : ELSE 210 : IVAA(ICPR(K)) = I2 211 : ENDIF 212 : 25 CONTINUE 213 : SEGDES MELEME 213 : CALL INITI( IVBB,ITES,0) 214 : * 215 : * Construction de CHPOIN resultat les composantes ont les noms 216 : * de celles des concentrations au centre ( ou aux faces) 217 : * 218 : SEGACT IPT1 219 : NPN=IPT1.NUM(/2) 220 : SEGDES IPT1 221 : NSOUPO=1 222 : NAT=1 223 : SEGINI MCHPOI 224 : MTYPOI=' ' 225 : MOCHDE=' CHPOIN CREE PAR HDEBI1 ' 226 : IFOPOI=IFOUR 227 : JATTRI(1)=2 228 : NC=NBCOMP 229 : SEGINI MSOUPO 230 : IPCHP(1)=MSOUPO 231 : DO 5 L=1,NBCOMP 232 : NOCOMP(L)=MSOUP1.NOCOMP(L) 233 : NOHARM(L)=MSOUP1.NOHARM(L) 234 : 5 CONTINUE 235 : IGEOC=IPFACE 236 : N=NPN 237 : SEGINI MPOVAL 238 : IPOVAL=MPOVAL 239 : NB=N*NC 239 : CALL INITD(VPOCHA,NB,0.D0) 240 : IF(ITTH.EQ.1) THEN 241 : C 242 : C cas des traces de charges récupération du nom des composantes 243 : C 244 : NBMAIL = MAHYBR(/1) 245 : DO 27 IMAIL=1,NBMAIL 246 : IF (MAHYBR(IMAIL).NE.0) THEN 247 : DESCR = IRIGEL(3,IMAIL) 248 : SEGACT DESCR 249 : MOREFD = LISDUA(1) 250 : MOREFP = LISINC(1) 251 : SEGDES DESCR 252 : GOTO 30 253 : ENDIF 254 : 27 CONTINUE 255 : 30 CONTINUE 256 : NOCOMP(1)=MOREFD 257 : ENDIF 258 : * 259 : * 260 : * 261 : C 262 : C-------------------------------------------------- 263 : *' Boucle 310 sur les OBJETS RIGIDITES ELEMENTAIRES' 264 : C-------------------------------------------------- 265 : C 266 : ITELEM = 0 267 : SEGACT MCHELM 268 : DO 310 IRI=1,NBMAIL 269 : C 270 : C Récupération MELEME ou Darcy est défini 271 : C 272 : MELEME = MAHYBR(IRI) 273 : IF (MELEME.EQ.0) GOTO 310 274 : SEGACT MELEME 275 : N1 = NUM(/1) 276 : N2 = NUM(/2) 277 : C 278 : C Récupération des infos pour la zone IRI dans le chapeau MRIGID 279 : C 280 : DESCR = IRIGEL(3,IRI) 281 : SEGACT DESCR 282 : NBDDL = LISDUA(/2) 283 : SEGDES DESCR 284 : SEGINI ITRAV 285 : xMATRI = IRIGEL(4,IRI) 286 : SEGACT xMATRI 287 : C 288 : C Activation du MELVAL du MCHAML d'orientation 289 : C 290 : MCHAML = ICHAML(IRI) 291 : SEGACT MCHAML 292 : MELVAL = IELVAL(1) 293 : SEGDES MCHAML 294 : SEGACT MELVAL 295 : * 296 : * Activation des objets necessaires à la prise en compte des forces de volumes 297 : * 298 : IF (IFORC.NE.0) THEN 299 : MCHAM2 = MCHEL2.ICHAML(IRI) 300 : SEGACT MCHAM2 301 : MELVA2 = MCHAM2.IELVAL(1) 302 : SEGDES MCHAM2 303 : SEGACT MELVA2 304 : xMATR1 = RI1.IRIGEL(4,IRI) 305 : SEGACT xMATR1 306 : ELSE 307 : DO 35 I=1,NBDDL 308 : RFOR(I)=0.D0 309 : 35 CONTINUE 310 : ENDIF 311 : C 312 : C------------------------------------------ 313 : *' Boucle 300 sur les MATRICES ELEMENTAIRES.' 314 : C------------------------------------------ 315 : C 316 : DO 300 I2=1,N2 317 : ITELEM = ITELEM + 1 318 : 319 : IF (IFORC.NE.0) THEN 320 : * 321 : *- calcul des flux de forces aux faces de l'element 322 : * 323 : DO 55 IDDL=1,NBDDL 324 : FLFOR(IDDL)= 0.D0 325 : IPOPTS = ICPR(NUM(IDDL,I2)) 326 : DO 50 I=1,IDIM 327 : FLFOR(IDDL) = FLFOR(IDDL) + MPOVA5.VPOCHA(IPOPTS,I) * 328 : S MELVA2.VELCHE(IDDL,I2) * MPOVA4.VPOCHA(ITELEM,I) * 329 : S MPOVA6.VPOCHA(IPOPTS,1) 330 : 50 CONTINUE 331 : 55 CONTINUE 332 : * 333 : *- Construction du tableau aux faces M.FORCE 334 : * 335 : * XMATR1 = IMATR1.IMATTT(I2) 336 : * SEGACT XMATR1 337 : DO 65 I=1,NBDDL 338 : RFOR(I)=0.D0 339 : DO 60 J=1,NBDDL 340 : RFOR(I) = RFOR(I) + XMATR1.RE(I,J,i2)*FLFOR(J) 341 : 60 CONTINUE 342 : 65 CONTINUE 343 : * SEGDES XMATR1 344 : ENDIF 345 : * 346 : * Recuperation de la matrice elementaire 347 : * 348 : * XMATRI = IMATTT(I2) 349 : * SEGACT XMATRI 350 : * 351 : *- De la somme des coefs pour une ligne 352 : *- -1 t 353 : *- LIGNE = RE * DIV 354 : *- -1 t 355 : DO 140 I=1,NBDDL 356 : RLIGNE(I) = 0.D0 357 : DO 130 J=1,NBDDL 358 : RLIGNE(I) = RLIGNE(I) + RE(I,J,i2) 359 : 130 CONTINUE 360 : 140 CONTINUE 361 : C 362 : C Calcul du flux aux faces 363 : C 364 : DO 200 IN=1,NBDDL 365 : NUMFA = ICPR(NUM(IN,I2)) 366 : IF (IVBB(NUMFA).EQ.0) THEN 367 : VVV= 0.D0 368 : IF(ICHP2.NE.0)THEN 369 : VVV=MPOVA3.VPOCHA(NUMFA,1) 370 : ENDIF 371 : DO 180 K=1,NBCOMP 372 : VA1 = 0.D0 373 : VA2 = 0.D0 374 : DO 190 JN=1,NBDDL 375 : VA1 = VA1+RE(IN,JN,i2)*(MPOVA1.VPOCHA( 376 : S IVAA(ICPR(NUM(JN,I2))),K)-RFOR(JN)) 377 : 190 CONTINUE 378 : VA2=RLIGNE(IN)*MPOVA2.VPOCHA(ITELEM,K) 379 : VA3= VVV*MPOVA1.VPOCHA(NUMFA,K) 380 : VPOCHA(NUMFA,K) = (VA2 -VA1+ VA3 ) * VELCHE(IN,I2) 381 : 180 CONTINUE 382 : IVBB(NUMFA)=1 383 : ENDIF 384 : 200 CONTINUE 385 : * SEGDES XMATRI 386 : 300 CONTINUE 387 : 388 : SEGDES MELVAL, xMATRI, MELEME 389 : SEGSUP ITRAV 390 : IF (IFORC.NE.0) THEN 391 : SEGDES MELVA2, xMATR1 392 : ENDIF 393 : 310 CONTINUE 394 : C 395 : C Nettoyage final 396 : C 397 : SEGDES MCHELM, MRIGID, IPMAHY, MSOUPO, MPOVAL, MCHEL2 398 : SEGDES MSOUP1, MPOVA1 399 : SEGDES MSOUP2, MPOVA2 400 : IF (MCHPO3.NE.0) SEGDES MPOVA3,MSOUP3 401 : IF (IFORC.NE.0) THEN 402 : SEGDES RI1 403 : SEGDES MSOUP4, MPOVA4 404 : SEGDES MSOUP5, MPOVA5 405 : SEGDES MSOUP6, MPOVA6 406 : ENDIF 407 : C 408 : SEGDES MCHPOI 409 : IRET = MCHPOI 410 : C 411 : SEGSUP IORGA, ICCPR 412 : C 413 : RETURN 414 : END 415 : 416 : 417 : 418 : 419 : 420 : 421 : 422 : 423 : 424 : 425 : 426 :
© Cast3M 2003 - All rights reserved.
Disclaimer