1 : C SMTP SOURCE CHAT 11/03/16 21:31:37 6902 2 : SUBROUTINE SMTP 3 : C----------------------------------------------------------------------- 4 : C 1) Calcul du second membre du systeme en trace de charge dans le cas 5 : C de la résolution des équations de Darcy par EFMH avec le modèle DARCY. 6 : C 2) Calcul du second membre du systeme en trace de charge et/ou d'une 7 : C contribution matricielle provenant de la convection dans le cas de la 8 : C résolution d'une équation de diffusion-convection par EFMH avec le 9 : C modèle DARCY. 10 : C----------------------------------------------------------------------- 11 : C 12 : C--------------------------- 13 : C Phrase d'appel (GIBIANE) : 14 : C--------------------------- 15 : C 16 : C (RIG2) (CHP3) = 'SMTP' MMODEL TABLE1 RIG1 (TABLE2) (CHP1) (CHP2) ; 17 : C 18 : C------------------------ 19 : C Operandes et resultat : 20 : C------------------------ 21 : C 22 : C MMODEL : Objet modele décrivant la formulation. 23 : C TABLE1 : TABLE DOMAINE contenant les maillages et les connectivités. 24 : C RIG1 : Matrices hybrides elementaires de DARCY crees par MHYB. 25 : C TABLE2 : TABLE DARCY_TRANSITOIRE contenant les infos transitoires. 26 : C CHP1 : CHPO centre de composante SOUR, source moyenne par élément. 27 : C CHP2 : CHPO face de composante FLUX, flux de vitesse convective. 28 : C CHP3 : CHPO face de composante FLUX, second membre du pb en TH. 29 : C RIG2 : RIGIDITE, contribution du terme convectif au LHS en TH. 30 : C 31 : C----------------------------------------------------------------------- 32 : C 33 : C Langage : ESOPE + FORTRAN77 34 : C 35 : C Auteurs : 09/93 F.DABBENE - Cas permanent 36 : C 09/94 X.NOUVELLON - Extension au transitoire 37 : C 12/95 F.DABBENE - Extension à la convection 38 : C 39 : C----------------------------------------------------------------------- 40 : IMPLICIT INTEGER(I-N) 41 : IMPLICIT REAL*8 (A-H,O-Z) 42 : * 43 : -INC CCOPTIO 44 : -INC SMCHAML 45 : -INC SMCHPOI 46 : -INC SMELEME 47 : -INC SMMODEL 48 : -INC SMRIGID 49 : -INC SMTABLE 50 : * 51 : SEGMENT IPMAHY 52 : INTEGER MAHYBR(NSOUS) 53 : ENDSEGMENT 54 : * 55 : LOGICAL LOGRE,LOGIN 56 : CHARACTER*4 NOMTOT(1) 57 : CHARACTER*6 CHAR6 58 : CHARACTER*8 TAPIND,TYPOBJ,CHARIN,CHARRE,LETYPE 59 : * 60 : * Initialisations 61 : * 62 : ISOU1 = 0 63 : IHN1 = 0 64 : ITPN1 = 0 65 : IVALIN = 0 66 : XVALIN = 0.D0 67 : LOGIN = .TRUE. 68 : IOBIN = 0 69 : TAPIND = 'MOT ' 70 : * 71 : * Lecture du MMODEL 72 : * 72 : CALL LIROBJ('MMODEL',IPMODE,1,IRET) 73 : IF (IERR.NE.0) RETURN 74 : MMODEL = IPMODE 75 : * 76 : * Lecture de la TABLE domaine 77 : * 78 : IPTABL = 0 78 : CALL LEKMOD(IPMODE,IPTABL,IRET) 79 : IF (IERR.NE.0) RETURN 80 : CHARIN = 'MAILLAGE' 81 : TYPOBJ = 'MAILLAGE' 81 : CALL LEKTAB(IPTABL,CHARIN,IOBRE) 82 : IF (IERR.NE.0) RETURN 83 : IPGEOM = IOBRE 83 : CALL LEKTAB(IPTABL,'ELTFA',IOBRE) 84 : IF (IERR.NE.0) RETURN 85 : IELTFA = IOBRE 85 : CALL LEKTAB(IPTABL,'CENTRE',IOBRE) 86 : IF (IERR.NE.0) RETURN 87 : ICENTR = IOBRE 87 : CALL LEKTAB(IPTABL,'FACE',IOBRE) 88 : IF (IERR.NE.0) RETURN 89 : IPFACE = IOBRE 90 : * 91 : * Lecture de RIGIDITE 92 : * 92 : CALL LIROBJ('RIGIDITE',IPRIGI,1,IRET) 93 : IF (IERR.NE.0) RETURN 94 : MRIGID = IPRIGI 95 : * 96 : * Lecture éventuelle de la TABLE Darcy_transitoire 97 : * 98 : IPTAB1 = 0 98 : CALL LIRTAB('DARCY_TRANSITOIRE',IPTAB1,0,IRET) 99 : IF (IERR.NE.0) RETURN 100 : IF (IRET.EQ.1) THEN 101 : THEMIN = -1.D-12 102 : THEMAX = 1.D0 - THEMIN 103 : TYPOBJ = 'FLOTTANT' 103 : CALL ACCTAB(IPTAB1,TAPIND,IVALIN,XVALIN,'THETA',LOGIN,IOBIN, 104 : S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE) 105 : IF (IERR.NE.0) RETURN 106 : THETA = XVALRE 107 : IF ((THETA.LT.THEMIN).OR.(THETA.GT.THEMAX)) THEN 108 : REAERR(1) = REAL(THETA) 109 : C REAERR(2) = REAL(0.D0) 110 : C REAERR(3) = REAL(1.D0) 111 : MOTERR(1:8) = ' THETA ' 111 : CALL ERREUR(42) 112 : RETURN 113 : ENDIF 114 : IF (THETA.LT.0.D0) THETA=0.D0 115 : IF (THETA.GT.1.D0) THETA=1.D0 116 : TYPOBJ = ' ' 116 : CALL ACCTAB(IPTAB1,TAPIND,IVALIN,XVALIN,'THETA_CONVECTION', 117 : S LOGIN,IOBIN, 118 : S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE) 119 : IF (IERR.NE.0) RETURN 120 : IF (TYPOBJ.EQ.'FLOTTANT') THEN 121 : THETAC = XVALRE 122 : IF ((THETAC.LT.THEMIN).OR.(THETAC.GT.THEMAX)) THEN 123 : REAERR(1) = REAL(THETAC) 124 : REAERR(2) = REAL(0.D0) 125 : REAERR(3) = REAL(1.D0) 126 : MOTERR(1:8) = 'TETACONV' 126 : CALL ERREUR(42) 127 : RETURN 128 : ENDIF 129 : IF (THETAC.LT.0.D0) THETA=0.D0 130 : IF (THETAC.GT.1.D0) THETA=1.D0 131 : ELSE 132 : THETAC = THETA 133 : ENDIF 133 : CALL ACCTAB(IPTAB1,TAPIND,IVALIN,XVALIN,'PAS',LOGIN,IOBIN, 134 : S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE) 135 : IF (IERR.NE.0) RETURN 136 : DELTAT = XVALRE 137 : IF (DELTAT.LE.0.D0) THEN 138 : REAERR(1) = REAL(DELTAT) 139 : REAERR(2) = REAL(0.D0) 140 : MOTERR(1:8) = ' DELTAT ' 140 : CALL ERREUR(41) 141 : RETURN 142 : ENDIF 143 : TYPOBJ = 'MCHAML ' 143 : CALL ACCTAB(IPTAB1,TAPIND,IVALIN,XVALIN,'SURF',LOGIN,IOBIN, 144 : S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE) 145 : IF (IERR.NE.0) RETURN 146 : IPCK = IOBRE 147 : TYPOBJ = 'CHPOINT ' 147 : CALL ACCTAB(IPTAB1,TAPIND,IVALIN,XVALIN,'TRACE',LOGIN,IOBIN, 148 : S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE) 149 : IF (IERR.NE.0) RETURN 150 : ITP = IOBRE 150 : CALL ACCTAB(IPTAB1,TAPIND,IVALIN,XVALIN,'CHARGE',LOGIN,IOBIN, 151 : S TYPOBJ,IVALRE,XVALRE,CHARRE,LOGRE,IOBRE) 152 : IF (IERR.NE.0) RETURN 153 : IH = IOBRE 154 : ELSE 155 : THETA = 0.D0 156 : DELTAT = 0.D0 157 : IPCK = 0 158 : ITP = 0 159 : IH = 0 160 : ENDIF 161 : * 162 : * Lecture des champoints facultatifs SOUR et FLUX 163 : * prise en compte des termes sources 164 : * prise en compte de la convection 165 : * 166 : ISOUR = 0 166 : CALL LIROBJ('CHPOINT',ISOUR,0,IRET) 167 : IF (IERR.NE.0) RETURN 168 : * 169 : IFLUX = 0 169 : CALL LIROBJ('CHPOINT',IFLUX,0,IRET) 170 : IF (IERR.NE.0) RETURN 171 : * 172 : * Tri des deux champoints facultatifs 173 : * 174 : IFAC = ISOUR + IFLUX 175 : IF (IFAC.EQ.0.AND.IPTAB1.EQ.0) THEN 176 : MOTERR(1:8) = 'CHAMPOIN' 176 : CALL ERREUR(37) 177 : RETURN 178 : ENDIF 179 : IF (IFAC.NE.0) THEN 180 : NBCOMP = 0 181 : INDIC = 0 181 : CALL QUEPOI(ISOUR,ICENTR,INDIC,NBCOMP,NOMTOT) 182 : IF (NOMTOT(1).NE.'SOUR'.AND.NOMTOT(1).NE.'FLUX') THEN 183 : MOTERR(1:4) = NOMTOT(1) 183 : CALL ERREUR(197) 184 : RETURN 185 : ELSE IF (NOMTOT(1).EQ.'FLUX') THEN 186 : IFAC = IFLUX 187 : IFLUX = ISOUR 188 : ISOUR = IFAC 189 : ENDIF 190 : ENDIF 191 : * 192 : * Dans le cas de la convection, récup. MCHAML orientation normale 193 : * 194 : IF (IFLUX.NE.0) THEN 194 : CALL LEKTAB(IPTABL,'XXNORMAE',IORIE) 195 : IF (IERR.NE.0) RETURN 196 : ENDIF 197 : * 198 : * Test du CHAMPOINT de composante SOUR de spg CENTRE 199 : * 200 : IF (ISOUR.NE.0) THEN 201 : NBCOMP = 1 202 : INDIC = 1 203 : NOMTOT(1) = 'SOUR' 203 : CALL QUEPOI(ISOUR,ICENTR,INDIC,NBCOMP,NOMTOT) 204 : IF (IERR.NE.0) RETURN 205 : ENDIF 206 : * 207 : * Test du CHAMPOINT de composante FLUX de spg FACE 208 : * 209 : IF (IFLUX.NE.0) THEN 210 : NBCOMP = 1 211 : INDIC = 1 212 : NOMTOT(1) = 'FLUX' 212 : CALL QUEPOI(IFLUX,IPFACE,INDIC,NBCOMP,NOMTOT) 213 : IF (IERR.NE.0) RETURN 214 : ENDIF 215 : * 216 : * Test du CHAMPOINT de composante TH de spg FACES 217 : * 218 : IF (ITP.NE.0) THEN 219 : NBCOMP = 1 220 : INDIC = 1 221 : NOMTOT(1) = 'TH' 221 : CALL QUEPOI(ITP,IPFACE,INDIC,NBCOMP,NOMTOT) 222 : IF (IERR.NE.0) RETURN 223 : ENDIF 224 : * 225 : * Test du CHAMPOINT de composante H de spg CENTRE 226 : * 227 : IF (IH.NE.0) THEN 228 : NBCOMP = 1 229 : INDIC = 1 230 : NOMTOT(1) = 'H' 230 : CALL QUEPOI(IH,ICENTR,INDIC,NBCOMP,NOMTOT) 231 : IF (IERR.NE.0) RETURN 232 : ENDIF 233 : * 234 : * Test de la formulation 235 : * 236 : SEGACT MMODEL 237 : NSOUS = KMODEL(/1) 238 : SEGINI IPMAHY 239 : IDARCY = 0 240 : IPT1=IPGEOM 241 : IPT2=IPGEOM 242 : DO 10 ISOUS=1,NSOUS 243 : IF(NSOUS.GT.1)THEN 244 : SEGACT IPT2 245 : IPT1=IPT2.LISOUS(ISOUS) 246 : SEGDES IPT2 247 : ENDIF 248 : IMODEL = KMODEL(ISOUS) 249 : SEGACT IMODEL 250 : LETYPE = FORMOD(1) 251 : IF (LETYPE.EQ.'DARCY') THEN 252 : IDARCY = IDARCY + 1 253 : MAHYBR(ISOUS) = IPT1 254 : ENDIF 255 : SEGDES IMODEL 256 : 10 CONTINUE 257 : SEGDES MMODEL 258 : IF (IDARCY.EQ.0) THEN 259 : MOTERR = LETYPE 259 : CALL ERREUR(193) 260 : GOTO 100 261 : ENDIF 262 : * 263 : * Recuperation des pointeurs ELTFA pour les zones ou DARCY est defini 264 : * 265 : MELEME = IELTFA 266 : SEGACT MELEME 267 : LZONES = LISOUS(/1) 268 : IF (LZONES.EQ.0) LZONES = 1 269 : IPT1 = IPGEOM 270 : SEGACT IPT1 271 : DO 30 ISOUS=1,NSOUS 272 : IMAMEL = MAHYBR(ISOUS) 273 : IF (IMAMEL.NE.0) THEN 274 : DO 20 ISZ=1,LZONES 275 : IF (LZONES.EQ.1) THEN 276 : IPT2 = IPT1 277 : IPT3 = MELEME 278 : ELSE 279 : IPT2 = IPT1.LISOUS(ISZ) 280 : IPT3 = LISOUS(ISZ) 281 : ENDIF 282 : IF (IPT2.EQ.IMAMEL) THEN 283 : MAHYBR(ISOUS) = IPT3 284 : GOTO 30 285 : ENDIF 286 : 20 CONTINUE 287 : IF (IMAMEL.EQ.MAHYBR(ISOUS)) THEN 288 : INTERR(1) = ISOUS 288 : CALL ERREUR(698) 289 : SEGDES IPT1 290 : SEGDES MELEME 291 : GOTO 100 292 : ENDIF 293 : ENDIF 294 : 30 CONTINUE 295 : SEGDES IPT1 296 : SEGDES MELEME 297 : * 298 : * Test du sous-type de la matrice de rigiditée récupérée 299 : * 300 : SEGACT MRIGID 301 : LETYPE = MTYMAT 302 : IF (LETYPE.NE.'DARCY') THEN 303 : MOTERR(1:8) = 'RIGIDITE' 304 : MOTERR(9:16) = 'DARCY ' 304 : CALL ERREUR(79) 305 : SEGDES MRIGID 306 : GOTO 100 307 : ENDIF 308 : * 309 : * Controle des pointeurs de MELEME de la rigidité 310 : * 311 : DO 40 ISOUS=1,NSOUS 312 : IMAMEL = MAHYBR(ISOUS) 313 : IF (IMAMEL.NE.0) THEN 314 : IPTEST = IRIGEL(1,ISOUS) 315 : IF (MAHYBR(ISOUS).NE.IRIGEL(1,ISOUS)) THEN 316 : MOTERR(1:8) = 'DARCY ' 317 : MOTERR(9:16) = 'ELTFA ' 318 : INTERR(1) = ISOUS 318 : CALL ERREUR(698) 319 : SEGDES MRIGID 320 : GOTO 100 321 : ENDIF 322 : ENDIF 323 : 40 CONTINUE 324 : SEGDES MRIGID 325 : * 326 : * Vérification du support du MCHAML % au MMODEL 327 : * Controle du nombre de composantes reelles 328 : * 329 : IF (IPCK.NE.0) THEN 330 : ISUP = 2 331 : ICOND = 1 331 : CALL QUESUP(IPMODE,IPCK,ISUP,ICOND,IRET,IRET2) 332 : IF (IRET.GT.1) GOTO 100 333 : MCHELM = IPCK 334 : SEGACT MCHELM 335 : DO 50 ISOUS=1,NSOUS 336 : IF (MAHYBR(ISOUS).NE.0) THEN 337 : MCHAML = ICHAML(ISOUS) 338 : SEGACT MCHAML 339 : N2 = TYPCHE(/2) 340 : IF (N2.GT.1) THEN 340 : CALL ERREUR(320) 341 : SEGDES MCHELM 342 : SEGDES MCHAML 343 : GOTO 100 344 : ENDIF 345 : CHAR6 = TYPCHE(1)(1:6) 346 : IF (CHAR6.NE.'REAL*8') THEN 347 : MOTERR(1:8) = NOMCHE(1) 347 : CALL ERREUR(679) 348 : SEGDES MCHELM 349 : SEGDES MCHAML 350 : GOTO 100 351 : ENDIF 352 : SEGDES MCHAML 353 : ENDIF 354 : 50 CONTINUE 355 : SEGDES MCHELM 356 : ENDIF 357 : * 358 : * Construction de MCHPOI 359 : * 360 : SEGDES IPMAHY 360 : CALL SMTP1(IPMAHY,IPFACE,IPRIGI,THETA,THETAC,DELTAT,IPCK, 361 : S ITP,IH,ISOUR,IFLUX,IORIE,MCHPOI,IRI1) 362 : * 362 : IF (MCHPOI.NE.0) CALL ECROBJ('CHPOINT',MCHPOI) 362 : IF (IRI1.NE.0) CALL ECROBJ('RIGIDITE',IRI1) 363 : * 364 : * Ménage 365 : * 366 : 100 CONTINUE 367 : SEGSUP IPMAHY 368 : RETURN 369 : END 370 : 371 : 372 : 373 : 374 : 375 : 376 : 377 : 378 : 379 : 380 : 381 : 382 : 383 : 384 : 385 :
© Cast3M 2003 - All rights reserved.
Disclaimer