1 : C ZDUDW SOURCE MAGN 10/05/18 21:17:45 6675 2 : SUBROUTINE ZDUDW(MTABX,MTAB1) 3 : IMPLICIT INTEGER(I-N) 4 : IMPLICIT REAL*8 (A-H,O-Z) 5 : C*********************************************************************** 6 : C 7 : C CET OPERATEUR DISCRETISE L'OPERATEUR DE PENALISATION DIV(U)=EPS*P 8 : C EN 2D SUR LES ELEMENTS QUA4 TRI3 TRI7 et QUA9 PLAN OU AXI 9 : C EN 3D SUR LES ELEMENTS CUB8 ET PRI6 10 : C 11 : C SYNTAXE : 12 : C --------- 13 : C 14 : C DUDW(EPS) INCO UN : 15 : C 16 : C COEFFICIENTS : 17 : C -------------- 18 : C 19 : C 20 : C EPS (SCAL DOMA) PARAMETRE DE PENALISATION 21 : C (SCAL ELEM) 22 : C 23 : C INCONNUES : 24 : C ----------- 25 : C 26 : C UN CHAMP DE VITESSE 27 : C 28 : C 29 : C OPTIONS : POROSITE DIV(PU)=EPS*P 30 : C SOURCE DIV(U)-Q=EPS*P 31 : C 32 : C************************************************************************ 33 : 34 : -INC CCOPTIO 35 : -INC CCGEOME 36 : -INC SMCOORD 37 : -INC SMLENTI 38 : -INC SMELEME 39 : POINTEUR MELEM1.MELEME,MELEMS.MELEME,MELEMI.MELEME,MELEMP.MELEME 40 : POINTEUR MELEMC.MELEME 41 : -INC SMCHPOI 42 : POINTEUR IZGG1.MPOVAL,SRCE.MPOVAL,IZTGG1.MPOVAL 43 : 44 : -INC SIZFFB 45 : POINTEUR IZF1.IZFFM 46 : -INC SMMATRIK 47 : 48 : -INC SMLMOTS 49 : POINTEUR LINCO.MLMOTS 50 : SEGMENT TRAV 51 : REAL*8 AF(NP,NP,NINC,NINC),AS(NP,NINC),CT(MP1,NP,NINC),PQ(MP1) 52 : ENDSEGMENT 53 : CHARACTER*8 NOMZ,NOMI,TYPE,TYPC 54 : CHARACTER*8 NOM0,NOM1,NOM2,NOM3 55 : DIMENSION IXV(3) 56 : PARAMETER (NTB=1) 57 : CHARACTER*8 LTAB(NTB) 58 : DIMENSION KTAB(NTB) 59 : DATA LTAB/'KIZX '/ 60 : C***************************************************************************** 61 : CDUDW 62 : C write(6,*)' Operateur DUDW ' 63 : 64 : C 65 : C- Récupération de la table INCO (pointeur KINC) 66 : C 66 : CALL LEKTAB(MTAB1,'INCO',KINC) 67 : IF(KINC.EQ.0)THEN 68 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 69 : MOTERR( 1: 8) = ' INCO ' 70 : MOTERR( 9:16) = ' INCO ' 71 : MOTERR(17:24) = ' EQEX ' 71 : CALL ERREUR(786) 72 : RETURN 73 : ENDIF 74 : 75 : C***************************************************************************** 76 : C OPTIONS 77 : C KIMPL = 0 -> EXPL 1 -> IMPL 2 -> SEMI 78 : C KFORM = 0 -> SI 1 -> EF 2 -> VF 3 -> EFMC 79 : C IDCEN = 0-> rien 1-> CENTREE 2-> SUPGDC 3-> SUPG 4-> TVISQUEU 5-> CNG 80 : 81 : IAXI=0 82 : IF(IFOMOD.EQ.0)IAXI=2 83 : C 84 : C- Récupération de la table des options KOPT (pointeur KOPTI) 85 : C 85 : CALL LEKTAB(MTABX,'KOPT',KOPTI) 86 : IF (KOPTI.EQ.0) THEN 87 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 88 : MOTERR( 1: 8) = ' KOPT ' 89 : MOTERR( 9:16) = ' KOPT ' 90 : MOTERR(17:24) = ' KIZX ' 90 : CALL ERREUR(786) 91 : RETURN 92 : ENDIF 93 : 93 : CALL ACME(KOPTI,'KIMPL',KIMPL) 93 : CALL ACME(KOPTI,'KPOIN',KPRE) 93 : CALL ACME(KOPTI,'KFORM',KFORM) 94 : 95 : IF(KFORM.NE.0.AND.KFORM.NE.1)THEN 96 : C Option %m1:8 incompatible avec les données 97 : MOTERR( 1: 8) = 'EF/EFM1 ' 97 : CALL ERREUR(803) 98 : RETURN 99 : ENDIF 100 : 101 : C write(6,*)' Apres les options ' 102 : C***************************************************************************** 103 : C 104 : C- Récupération de la table DOMAINE associée au domaine local 105 : C 105 : CALL ACMM(MTABX,'NOMZONE',NOMZ) 105 : CALL LEKTAB(MTABX,'DOMZ',MTABZ) 106 : IF(MTABZ.EQ.0)THEN 107 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 108 : MOTERR( 1: 8) = ' DOMZ ' 109 : MOTERR( 9:16) = ' DOMZ ' 110 : MOTERR(17:24) = ' KIZX ' 110 : CALL ERREUR(786) 111 : RETURN 112 : ENDIF 113 : 113 : CALL LEKTAB(MTABZ,'SOMMET',MELEMS) 113 : CALL LEKTAB(MTABZ,'MAILLAGE',MELEME) 114 : MACRO1=0 114 : CALL LEKTAB(MTABZ,'MACRO1',MELEMI) 115 : IF(MELEMI.NE.0)MACRO1=1 115 : CALL LEKTAB(MTABZ,'QUADRATI',MQUAD) 116 : IF (IERR.NE.0) RETURN 117 : 118 : IF(MQUAD.EQ.0.AND.MACRO1.EQ.0.AND.KPRE.NE.2)THEN 119 : WRITE(6,*)'Operateur DUDW ' 120 : WRITE(6,*)'Type d''éléments non prévu' 121 : RETURN 122 : ENDIF 123 : 124 : IF(KPRE.EQ.2)THEN 124 : CALL LEKTAB(MTABZ,'CENTRE',MELEMC) 125 : MELEMP=MELEMC 126 : MELEMI=MELEME 127 : ELSEIF(KPRE.EQ.3)THEN 127 : CALL LEKTAB(MTABZ,'CENTREP0',MELEMC) 128 : MELEMP=MELEMC 129 : ELSEIF(KPRE.EQ.4)THEN 129 : CALL LEKTAB(MTABZ,'CENTREP1',MELEMC) 129 : CALL LEKTAB(MTABZ,'ELTP1NC ',MELEMP) 130 : ENDIF 131 : 132 : SEGACT MELEME 133 : SEGACT MELEMP 134 : C************************************************************************* 135 : C VERIFICATIONS SUR LES INCONNUES 136 : C write(6,*)' Verification des inconnues ' 137 : 138 : TYPE='LISTMOTS' 138 : CALL ACMO(MTABX,'LISTINCO',TYPE,LINCO) 139 : SEGACT LINCO 140 : 141 : NBINC=LINCO.MOTS(/2) 142 : IF(NBINC.NE.1)THEN 143 : C Indice %m1:8 : contient plus de %i1 %m9:16 144 : MOTERR( 1:8) = 'LISTINCO' 145 : INTERR(1) = 1 146 : MOTERR(9:16) = ' MOTS ' 146 : CALL ERREUR(799) 147 : RETURN 148 : ENDIF 149 : 150 : NOMI=LINCO.MOTS(1) 151 : 152 : TYPE=' ' 152 : CALL ACMO(KINC,NOMI,TYPE,MCHPOI) 153 : IF(TYPE.NE.'CHPOINT ')THEN 154 : C Indice %m1:8 : ne contient pas un objet de type %m9:16 155 : MOTERR( 1: 8) = 'INC '//NOMI 156 : MOTERR( 9:16) = 'CHPOINT ' 156 : CALL ERREUR(800) 157 : RETURN 158 : ELSE 158 : CALL LICHT(MCHPOI,IZTU1,TYPC,MELEM1) 159 : ENDIF 160 : 161 : C************************************************************************* 162 : C Le domaine de definition est donne par le SPG de la premiere inconnue 163 : C Les inconnues suivantes devront posseder ce meme pointeur 164 : C On verifie que les points de la zone sont tous inclus dans ce SPG 165 : 165 : CALL KRIPAD(MELEM1,MLENTI) 166 : 166 : CALL VERPAD(MLENTI,MELEME,IRET) 167 : IF(IRET.NE.0)THEN 168 : C Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique 169 : MOTERR(1: 8) = 'INC '//NOMI 170 : MOTERR(9:16) = 'CHPOINT ' 170 : CALL ERREUR(788) 171 : RETURN 172 : ENDIF 173 : 174 : NINC=IDIM 175 : 176 : MLENT1=MLENTI 177 : 178 : C***************************************************************************** 179 : C Lecture du ou des coefficients 180 : C write(6,*)' Lecture des coefficients' 181 : 181 : CALL ACME(MTABX,'IARG',IARG) 182 : 183 : C? IXV(1)=MELEMC 184 : IXV(1)=0 185 : IXV(2)=1 186 : IXV(3)=0 186 : CALL LEKCOF('Opérateur DUDW :', 187 : & MTABX,KINC,1,IXV,IZTG1,IZTGG1,NPT1,NC1,IK1,IRET) 188 : IF(IRET.EQ.0)RETURN 189 : KITT=2 190 : 191 : INDGS=0 192 : SRCE=IZTGG1 193 : MPOVA1 =IZTGG1 194 : MPOVA3 =IZTGG1 195 : 196 : IF(IARG.EQ.2)THEN 197 : IXV(1)=MELEMC 198 : IXV(2)=1 199 : IXV(3)=0 199 : CALL LEKCOF('Opérateur DUDW :', 200 : & MTABX,KINC,2,IXV,IZS,SRCE,NPTS,NCS,IKS,IRET) 201 : IF(IRET.EQ.0)RETURN 202 : INDGS=1 202 : CALL KRIPAD(MELEMC,MLENT1) 203 : ENDIF 204 : 205 : C************************************************************************* 206 : IF(KFORM.EQ.0)THEN 207 : C CAS FORMULATION EF SI (GRESHO) 208 : 209 : IF(KIMPL.NE.0)THEN 210 : GO TO 90 211 : ENDIF 212 : 213 : WRITE(6,*)' Operateur DUDW ' 214 : WRITE(6,*)' Cas Formulation EF SI ' 215 : WRITE(6,*)' Cas invalide ' 216 : GO TO 90 217 : 218 : C************************************************************************* 219 : ELSE IF(KFORM.EQ.1)THEN 220 : C CAS FORMULATION EF 221 : 222 : NUTOEL=0 223 : 224 : SEGACT MELEMI 225 : NBSOUS=MELEMI.LISOUS(/1) 226 : IF(NBSOUS.EQ.0)NBSOUS=1 227 : 228 : NRIGE=7 229 : NKID =9 230 : NKMT =7 231 : NMATRI=1 232 : SEGINI MATRIK 233 : IRIGEL(1,1)=MELEMI 234 : IRIGEL(2,1)=MELEMI 235 : NBOP=0 236 : NBME=NINC*NINC 237 : C NBME=2 238 : NBELC=0 239 : SEGINI IMATRI 240 : IRIGEL(4,1)=IMATRI 241 : KSPGP=MELEMS 242 : KSPGD=MELEMS 243 : 244 : IF(NBME.EQ.4)THEN 245 : NOM1='1'//NOMI(1:7) 246 : NOM2='2'//NOMI(1:7) 247 : LISPRI(1)=NOM1 248 : LISDUA(1)=NOM1 249 : LISPRI(2)=NOM2 250 : LISDUA(2)=NOM1 251 : LISPRI(3)=NOM1 252 : LISDUA(3)=NOM2 253 : LISPRI(4)=NOM2 254 : LISDUA(4)=NOM2 255 : ELSEIF(NBME.EQ.9)THEN 256 : NOM1='1'//NOMI(1:7) 257 : NOM2='2'//NOMI(1:7) 258 : NOM3='3'//NOMI(1:7) 259 : LISPRI(1)=NOM1 260 : LISDUA(1)=NOM1 261 : LISPRI(2)=NOM2 262 : LISDUA(2)=NOM1 263 : LISPRI(3)=NOM3 264 : LISDUA(3)=NOM1 265 : 266 : LISPRI(4)=NOM1 267 : LISDUA(4)=NOM2 268 : LISPRI(5)=NOM2 269 : LISDUA(5)=NOM2 270 : LISPRI(6)=NOM3 271 : LISDUA(6)=NOM2 272 : 273 : LISPRI(7)=NOM1 274 : LISDUA(7)=NOM3 275 : LISPRI(8)=NOM2 276 : LISDUA(8)=NOM3 277 : LISPRI(9)=NOM3 278 : LISDUA(9)=NOM3 279 : 280 : ELSEIF(NBME.EQ.2)THEN 281 : 282 : LISPRI(1)='1'//NOMI(1:7) 283 : LISDUA(1)='1'//NOMI(1:7) 284 : LISPRI(2)='2'//NOMI(1:7) 285 : LISDUA(2)='2'//NOMI(1:7) 286 : 287 : ELSE 288 : WRITE(6,*)' Operateur DUDW ' 289 : WRITE(6,*)' Cas invalide' 290 : GO TO 90 291 : ENDIF 292 : 293 : 294 : IF(INDGS.NE.0)THEN 295 : NAT=2 296 : NSOUPO=1 297 : SEGACT MELEMS 298 : N=MELEMS.NUM(/2) 299 : NC=NINC 300 : SEGINI MCHPO1,MSOUP1,MPOVA1 301 : MCHPO1.IFOPOI=IFOMOD 302 : MCHPO1.MOCHDE=TITREE 303 : MCHPO1.MTYPOI='SMBR' 304 : MCHPO1.JATTRI(1)=2 305 : MCHPO1.IPCHP(1)=MSOUP1 306 : MSOUP1.NOCOMP(1)=NOM1 307 : MSOUP1.NOCOMP(2)=NOM2 308 : IF(NINC.EQ.3)MSOUP1.NOCOMP(3)=NOM3 309 : MSOUP1.IGEOC=MELEMS 310 : MSOUP1.IPOVAL=MPOVA1 311 : 312 : SEGACT MELEMC 313 : N=MELEMC.NUM(/2) 314 : NC=1 315 : SEGINI MCHPO3,MSOUP3,MPOVA3 316 : MCHPO3.IFOPOI=IFOMOD 317 : MCHPO3.MOCHDE=TITREE 318 : MCHPO3.MTYPOI='SMBR' 319 : MCHPO3.JATTRI(1)=2 320 : MCHPO3.IPCHP(1)=MSOUP3 321 : MSOUP3.NOCOMP(1)='SCAL' 322 : MSOUP3.IGEOC=MELEMC 323 : MSOUP3.IPOVAL=MPOVA3 324 : 325 : ENDIF 326 : 327 : DO 101 L=1,NBSOUS 328 : IPT1=MELEMI 329 : IPT2=MELEMI 330 : IF(NBSOUS.NE.1)IPT1=MELEMI.LISOUS(L) 331 : SEGACT IPT1 332 : IF(INDGS.NE.0)THEN 333 : IPT2=MELEMP 334 : IF(NBSOUS.NE.1)IPT2=MELEMP.LISOUS(L) 335 : SEGACT IPT2 336 : ENDIF 337 : 338 : IF(MQUAD.NE.0)THEN 339 : IF(KPRE.EQ.2)NOM0=NOMS(IPT1.ITYPEL)//'PRP0' 340 : IF(KPRE.EQ.3)NOM0=NOMS(IPT1.ITYPEL)//'PRP0' 341 : IF(KPRE.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'PRP1' 342 : ELSEIF(MACRO1.NE.0)THEN 343 : IF(KPRE.EQ.2)NOM0=NOMS(IPT1.ITYPEL)//' ' 344 : IF(KPRE.EQ.3)NOM0=NOMS(IPT1.ITYPEL)//'MCP0' 345 : IF(KPRE.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'MCP1' 346 : ELSE 347 : IF(KPRE.EQ.2)NOM0=NOMS(IPT1.ITYPEL)//' ' 348 : ENDIF 349 : 349 : CALL KALPBG(NOM0,'FONFORM ',IZFFM) 350 : 351 : 352 : SEGACT IZFFM*MOD 353 : IZHR=KZHR(1) 354 : SEGACT IZHR*MOD 355 : NES=GR(/1) 356 : NPG=GR(/3) 357 : IZF1=KTP(1) 358 : SEGACT IZF1*MOD 359 : MP1=IZF1.FN(/1) 360 : NP = IPT1.NUM(/1) 361 : 362 : MP = NP 363 : NBEL=IPT1.NUM(/2) 364 : SEGINI IPM1,IPM2,IPM3,IPM4 365 : LIZAFM(L,1)=IPM1 366 : C LIZAFM(L,2)=IPM4 367 : LIZAFM(L,2)=IPM2 368 : LIZAFM(L,3)=IPM3 369 : LIZAFM(L,4)=IPM4 370 : IPM5=IPM1 371 : IPM6=IPM1 372 : IPM7=IPM1 373 : IPM8=IPM1 374 : IPM9=IPM1 375 : IF(NBME.EQ.9)THEN 376 : SEGINI IPM5,IPM6,IPM7,IPM8,IPM9 377 : LIZAFM(L,5)=IPM5 378 : LIZAFM(L,6)=IPM6 379 : LIZAFM(L,7)=IPM7 380 : LIZAFM(L,8)=IPM8 381 : LIZAFM(L,9)=IPM9 382 : ENDIF 383 : 384 : C Pour l'instant les sources et puits de masse ne sont pas actifs 385 : INDG2=0 386 : KI2=0 387 : KJ2=0 388 : 389 : SEGINI TRAV 390 : 391 : NPT=MPOVA1.VPOCHA(/1) 391 : CALL XDUDW(FN,IZF1.FN,GR,PG,XYZ,HR,PGSQ,RPG,AJ, 392 : & NES,IDIM,NP,MP1,NPG,IAXI,NINC, 393 : & IZTGG1.VPOCHA,IK1,SRCE.VPOCHA,INDGS,IKS, 394 : & IPT1.NUM,NBEL,NUTOEL,XCOOR,AF,AS,CT,PQ, 395 : & IPM1.AM,IPM2.AM,IPM3.AM,IPM4.AM,IPM5.AM,IPM6.AM,IPM7.AM, 396 : & IPM8.AM,IPM9.AM,MPOVA1.VPOCHA,NPT,LECT,IPT2.NUM,MLENT1.LECT, 397 : & MPOVA3.VPOCHA) 398 : 399 : SEGSUP TRAV 400 : SEGDES IPT1 401 : NUTOEL=NUTOEL+NBEL 402 : 101 CONTINUE 403 : SEGDES IMATRI 404 : SEGDES MELEME,MATRIK 405 : IF(IK1.EQ.1)THEN 406 : SEGSUP IZTGG1 407 : ELSE 408 : SEGDES IZTGG1 409 : ENDIF 409 : CALL ECMO(MTABX,'MATELM','MATRIK',MATRIK) 410 : 411 : 412 : IF(INDGS.NE.0)THEN 413 : 414 : TYPE=' ' 414 : CALL ACMO(MTAB1,'SMBR',TYPE,MCHPO2) 415 : IF(TYPE.NE.'CHPOINT')THEN 415 : CALL ECMO(MTAB1,'SMBR','CHPOINT',MCHPO1) 416 : ELSE 416 : CALL ECROBJ('CHPOINT',MCHPO2) 416 : CALL ECROBJ('CHPOINT',MCHPO1) 416 : CALL PRFUSE 416 : CALL LIROBJ('CHPOINT',MCHPOI,1,IRET) 416 : CALL ECMO(MTAB1,'SMBR','CHPOINT',MCHPOI) 417 : ENDIF 418 : 419 : TYPE=' ' 419 : C CALL ACMO(MTAB1,'DUDWSRCE',TYPE,MCHPO4) 420 : C IF(TYPE.NE.'CHPOINT')THEN 420 : CALL ECMO(MTAB1,'DUDWSRCE','CHPOINT',MCHPO3) 421 : C ELSE 421 : C CALL ECROBJ('CHPOINT',MCHPO3) 421 : C CALL ECROBJ('CHPOINT',MCHPO4) 421 : C CALL PRFUSE 421 : C CALL LIROBJ('CHPOINT',MCHPOI,1,IRET) 421 : C CALL ECMO(MTAB1,'DUDWSRCE','CHPOINT',MCHPOI) 422 : C ENDIF 423 : 424 : ENDIF 425 : 426 : SEGSUP MLENTI 427 : IF(INDGS.NE.0)SEGSUP MLENT1 428 : 429 : RETURN 430 : C************************************************************************* 431 : ELSE IF(KFORM.EQ.2)THEN 432 : C CAS FORMULATION VF 433 : 434 : RETURN 435 : ENDIF 436 : C************************************************************************* 437 : 1002 FORMAT(10(1X,1PE11.4)) 438 : 90 CONTINUE 439 : WRITE(6,*)' Interuption anormale de DUDW' 440 : RETURN 441 : END 442 : 443 : 444 : 445 : 446 : 447 : 448 : 449 : 450 :
© Cast3M 2003 - All rights reserved.
Disclaimer