1 : C ZNS SOURCE MAGN 10/05/18 21:17:56 6675 2 : SUBROUTINE ZNS(MTABX,MTAB1) 3 : PARAMETER (NTB=1) 4 : IMPLICIT INTEGER(I-N) 5 : IMPLICIT REAL*8 (A-H,O-Z) 6 : C*********************************************************************** 7 : C 8 : C CET OPERATEUR DISCRETISE LES EQUATIONS DE NAVIER STOKES 9 : C EN 2D SUR LES ELEMENTS QUA4 ET TRI3 PLAN OU AXI 10 : C EN 3D SUR LES ELEMENTS CUB8 ET PRI6 11 : C LES OPERATEURS SONT "SOUS-INTEGRES" 12 : C 13 : C SYNTAXE : 14 : C --------- 15 : C 1/ Cas incompréssible 16 : C 17 : C du/dt + u Grad u = nu Lapl u - 1/ro Grad p < + S > 18 : C < + g beta (T-Tref) > 19 : C 20 : C 'OPER' 'NS' nu 'INCO' UN : 21 : C 'OPER' 'NS' nu s 'INCO' UN : 22 : C 'OPER' 'NS' nu gb tn tref 'INCO' UN : 23 : C 24 : C 25 : C 2/ Cas compréssible 26 : C 27 : C dG/dt + u Grad G + G Div u = mu Lapl u - Grad p < + S > 28 : C 29 : C 'OPER' 'NS' mu un 'INCO' GN : 30 : C 'OPER' 'NS' mu un s 'INCO' GN : 31 : C 32 : C 33 : C nu,mu viscosité cinématique resp. dynamique 34 : C FLOTTANT ou CHPOINT SCAL CENTRE 35 : C s source volumique de qdm 36 : C POINT ou CHPOINT VECT CENTRE 37 : C gb coéfficient de flottabilité (g*beta où g est l'accéllération 38 : C de la pesanteur et beta le coéfficient de dilatabilité) 39 : C POINT ou CHPOINT VECT CENTRE 40 : C tn Champ de température CHPOINT SCAL SOMMET 41 : C tref température de référence 42 : C FLOTTANT ou CHPOINT SCAL SOMMET 43 : C 44 : C Champ de vitesse -> VITESS 45 : C un Champ de vitesse transportant -> UTRANS 46 : C CHPOINT VECT SOMMET 47 : C gn Champ de vitesse massique (transporté) -> IZTU1 (Inconnue) 48 : C CHPOINT VECT SOMMET 49 : C 50 : C************************************************************************ 51 : 52 : C?-INC CCVQUA4 53 : -INC CCOPTIO 54 : -INC CCGEOME 55 : -INC SIZFFB 56 : POINTEUR IZF1.IZFFM 57 : -INC SMCHAML 58 : -INC SMCOORD 59 : -INC SMLENTI 60 : POINTEUR IPADI.MLENTI,IPADU.MLENTI,IPADF.MLENTI,IPADS.MLENTI 61 : POINTEUR IPADQ.MLENTI 62 : -INC SMELEME 63 : POINTEUR MELEM1.MELEME,MELEMC.MELEME,IGEOM0.MELEME,MELEMS.MELEME 64 : POINTEUR MELEMI.MELEME,MELEP1.MELEME 65 : -INC SMCHPOI 66 : POINTEUR IZTU1.MPOVAL,IZGG1.MPOVAL 67 : POINTEUR VITESS.MPOVAL,UTRANS.MPOVAL 68 : POINTEUR VISCO.MPOVAL,IZTGG2.MPOVAL 69 : POINTEUR IZTGG3.MPOVAL,IZTGG4.MPOVAL 70 : POINTEUR IZVOL.MPOVAL,IZTCO.MPOVAL 71 : 72 : -INC SMMATRIK 73 : POINTEUR IPM.IZAFM 74 : 75 : SEGMENT IMATRS 76 : INTEGER LIZAFS(NBSOUS,NBME) 77 : ENDSEGMENT 78 : POINTEUR IPMS.IZAFM,IPS1.IZAFM,IPS2.IZAFM,IPS3.IZAFM 79 : 80 : -INC SMLMOTS 81 : POINTEUR LINCO.MLMOTS 82 : CHARACTER*8 NOMZ,NOMI,TYPE,TYPC,NOM0,NOMA,NOM 83 : CHARACTER*8 LTAB(NTB) 84 : DIMENSION KTAB(NTB),IXV(3),RO(1) 85 : SAVE IPAS 86 : DATA LTAB/'KIZX '/,IPAS/0/,RO/1.D0/ 87 : C***************************************************************************** 88 : CNS 89 : C write(6,*)' DEBUT NS ' 90 : 91 : C 92 : C- Récupération de la table INCO (pointeur KINC) 93 : C 93 : CALL LEKTAB(MTAB1,'INCO',KINC) 94 : IF(KINC.EQ.0)THEN 95 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 96 : MOTERR( 1: 8) = ' INCO ' 97 : MOTERR( 9:16) = ' INCO ' 98 : MOTERR(17:24) = ' EQEX ' 98 : CALL ERREUR(786) 99 : RETURN 100 : ENDIF 101 : 102 : C***************************************************************************** 103 : C OPTIONS 104 : C KIMPL = 0 -> EXPL 1 -> IMPL 2 -> SEMI 105 : C KFORM = 0 -> SI 1 -> EF 2 -> VF 3 -> EFMC 106 : C IDCEN = 0-> rien 1-> CENTREE 2-> SUPGDC 3-> SUPG 4-> TVISQUEU 5-> CNG 107 : 108 : IAXI=0 109 : IF(IFOMOD.EQ.0)IAXI=2 110 : C 111 : C- Récupération de la table des options KOPT (pointeur KOPTI) 112 : C 112 : CALL LEKTAB(MTABX,'KOPT',KOPTI) 113 : IF (KOPTI.EQ.0) THEN 114 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 115 : MOTERR( 1: 8) = ' KOPT ' 116 : MOTERR( 9:16) = ' KOPT ' 117 : MOTERR(17:24) = ' KIZX ' 117 : CALL ERREUR(786) 118 : RETURN 119 : ENDIF 120 : 120 : CALL ACME(KOPTI,'MTRMASS ',MMPG) 121 : IPG=0 122 : IF(MMPG.EQ.3)IPG=1 122 : CALL ACME(KOPTI,'IDCEN',IDCEN) 122 : CALL ACME(KOPTI,'IKOMP',IKOMP) 122 : CALL ACME(KOPTI,'ALE',KKALE) 122 : CALL ACME(KOPTI,'KIMPL',KIMPL) 122 : CALL ACME(KOPTI,'KPOIN',KPRE) 122 : CALL ACME(KOPTI,'KFORM',KFORM) 123 : 124 : IF(KFORM.NE.0.AND.KFORM.NE.1)THEN 125 : C Option %m1:8 incompatible avec les données 126 : MOTERR( 1: 8) = 'EF/EFM1 ' 126 : CALL ERREUR(803) 127 : RETURN 128 : ENDIF 128 : CALL ACME(KOPTI,'KMACO',KMACO) 128 : CALL ACMF(KOPTI,'AIMPL',AIMPL) 129 : IF (IERR.NE.0) RETURN 130 : 131 : C write(6,*)' Apres les options ' 132 : C***************************************************************************** 133 : C 134 : C- Récupération de la table DOMAINE associée au domaine local 135 : C 135 : CALL ACMM(MTABX,'NOMZONE',NOMZ) 135 : CALL LEKTAB(MTABX,'DOMZ',MTABZ) 136 : IF(MTABZ.EQ.0)THEN 137 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 138 : MOTERR( 1: 8) = ' DOMZ ' 139 : MOTERR( 9:16) = ' DOMZ ' 140 : MOTERR(17:24) = ' KIZX ' 140 : CALL ERREUR(786) 141 : RETURN 142 : ENDIF 143 : 143 : CALL LEKTAB(MTABZ,'MAILLAGE',MELEME) 143 : CALL LEKTAB(MTABZ,'SOMMET',MELEMS) 143 : CALL LEKTAB(MTABZ,'CENTRE',MELEMC) 143 : CALL LEKTAB(MTABZ,'MACRO',MACRO) 143 : CALL LEKTAB(MTABZ,'MACRO1',MELEMI) 143 : CALL LEKTAB(MTABZ,'QUADRATI',MQUAD) 144 : IF (IERR.NE.0) RETURN 145 : 146 : MELEMQ=MELEMC 147 : MELEP1=MELEMC 148 : C write(6,*)' KPRE=',kpre,' MACRO=',macro,' QUADR=',MQUAD 149 : C? IF(KPRE.LE.3.AND.MACRO.EQ.0.AND.MQUAD.EQ.0.AND.KFORM.NE.0)THEN 150 : C? WRITE(6,*)' Operateur NS ' 151 : C? WRITE(6,*)' Incompatibilité du terme source et des éléments' 152 : C? WRITE(6,*)' MACRO ou QUADRATIQUE attendu ' 153 : C? WRITE(6,*)' Interuption anormale de NS' 154 : C Option %m1:8 incompatible avec les données 155 : C? MOTERR( 1: 8) = ' EF ' 155 : C? CALL ERREUR(803) 156 : C? RETURN 157 : IF(KPRE.EQ.3)THEN 157 : CALL LEKTAB(MTABZ,'CENTREP0',MELEMQ) 158 : MELEP1=MELEMQ 159 : ELSEIF(KPRE.EQ.4)THEN 159 : CALL LEKTAB(MTABZ,'CENTREP1',MELEMQ) 159 : CALL LEKTAB(MTABZ,'ELTP1NC ',MELEP1) 160 : ENDIF 161 : 162 : C************************************************************************* 163 : C VERIFICATIONS SUR LES INCONNUES 164 : C 165 : C- Récupération du nombre d'inconnues et du nom de l'inconnue NOMI 166 : C 167 : TYPE='LISTMOTS' 167 : CALL ACMO(MTABX,'LISTINCO',TYPE,LINCO) 168 : IF (IERR.NE.0) RETURN 169 : SEGACT LINCO 170 : NBINC=LINCO.MOTS(/2) 171 : IF(NBINC.NE.1)THEN 172 : C Indice %m1:8 : contient plus de %i1 %m9:16 173 : MOTERR( 1:8) = 'LISTINCO' 174 : INTERR(1) = 1 175 : MOTERR(9:16) = ' MOTS ' 175 : CALL ERREUR(799) 176 : RETURN 177 : ENDIF 178 : 179 : NOMI=LINCO.MOTS(1) 180 : NOMA=NOMI 181 : C 182 : C- Récupération de l'inconnue 183 : C 184 : TYPE=' ' 184 : CALL ACMO(KINC,NOMI,TYPE,MCHPOI) 185 : IF(TYPE.NE.'CHPOINT ')THEN 186 : C Indice %m1:8 : ne contient pas un objet de type %m9:16 187 : MOTERR( 1: 8) = 'INC '//NOMI 188 : MOTERR( 9:16) = 'CHPOINT ' 188 : CALL ERREUR(800) 189 : RETURN 190 : ELSE 190 : CALL LICHT(MCHPOI,IZTU1,TYPC,MELEM1) 191 : NINKO = IZTU1.VPOCHA(/2) 192 : IF (NINKO.NE.IDIM) THEN 193 : C Indice %m1:8 : Le %m9:16 n'a pas le bon nombre de composantes 194 : MOTERR( 1: 8) = 'INC '//NOMI 195 : MOTERR( 9:16) = 'CHPOINT ' 195 : CALL ERREUR(784) 196 : RETURN 197 : ENDIF 198 : C write(6,*)' MCHPOI,MELEM1=',MCHPOI,MELEM1 199 : C On fait pointer ces deux tableaux sur le champ U inconu (tjs présent) pour 200 : C eviter de les enlever lors de l'appel FORTRAN si les options sont absentes 201 : UTRANS=IZTU1 202 : IKW=0 203 : VITESS=IZTU1 204 : ENDIF 205 : 206 : C***************************************************************************** 207 : C Le domaine de definition est donne par le SPG de la premiere inconnue 208 : C Les inconnues suivantes devront posseder ce meme pointeur 209 : C On verifie que les points de la zone sont tous inclus dans ce SPG 210 : 210 : CALL KRIPAD(MELEM1,IPADI) 211 : IPADU=IPADI 212 : IPADF=IPADI 213 : IPADS=IPADI 213 : IF(MELEM1.NE.MELEMS)CALL KRIPAD(MELEMS,IPADS) 214 : 215 : IF(IPAS.EQ.0)THEN 216 : IPAS=1 216 : CALL VERPAD(IPADI,MELEME,IRET) 217 : IF(IRET.NE.0)THEN 218 : C Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique 219 : MOTERR(1: 8) = 'INC '//NOMI 220 : MOTERR(9:16) = 'CHPOINT ' 220 : CALL ERREUR(788) 220 : CALL ERREUR(788) 221 : RETURN 222 : ENDIF 223 : ENDIF 224 : 225 : C***************************************************************************** 226 : 227 : 228 : 229 : C************************************************************************* 230 : C Lecture des coefficient 231 : C Type du coefficient : 232 : C IK1=0 CHPOINT IK1=1 scalaire IK1=2 vecteur 233 : 233 : CALL ACME(MTABX,'IARG',IARG) 234 : IF(IKOMP.EQ.0)THEN 235 : IF(IARG.NE.1.AND.IARG.NE.2.AND.IARG.NE.3.AND.IARG.NE.4. 236 : & AND.IARG.NE.5)THEN 237 : WRITE(6,*)' Operateur NS : option incompréssible ' 238 : WRITE(6,*)' Nombre d''arguments incorrect : ',IARG 239 : WRITE(6,*)' On attend 1 , 2 ou 4 ' 240 : C Indice %m1:8 : nombre d'arguments incorrect 241 : MOTERR(1:8) = 'IARG ' 241 : CALL ERREUR(804) 242 : RETURN 243 : ENDIF 244 : ELSEIF(IKOMP.EQ.1)THEN 245 : IF(IARG.NE.2.AND.IARG.NE.3)THEN 246 : WRITE(6,*)' Operateur NS : option compréssible ' 247 : WRITE(6,*)' Nombre d''arguments incorrect : ',IARG 248 : WRITE(6,*)' On attend 2 ou 3 ' 249 : C Indice %m1:8 : nombre d'arguments incorrect 250 : MOTERR(1:8) = 'IARG ' 250 : CALL ERREUR(804) 251 : RETURN 252 : ENDIF 253 : ENDIF 254 : 255 : IF(KKALE.NE.0)IARG=IARG-1 256 : C--Cas incompréssible 257 : IF(IKOMP.EQ.0)THEN 258 : IXV(1)=MELEMC 259 : IXV(2)=1 260 : IXV(3)=0 260 : CALL LEKCOF('Opérateur NS :', 261 : & MTABX,KINC,1,IXV,IZTG1,VISCO,NPT1,NC1,IK1,IRET) 262 : IF(IRET.EQ.0)RETURN 263 : 264 : IK2=-1 265 : IKG=-1 266 : IK3=-1 267 : IK4=-1 268 : IZTGG2=VISCO 269 : IZTGG3=VISCO 270 : IZTGG4=VISCO 271 : 272 : C 2ème coef -> source S 273 : IF(IARG.GE.2)THEN 274 : IXV(1)=-MELEMQ 275 : IXV(2)=0 276 : IXV(3)=1 276 : CALL LEKCOF('Opérateur NS :', 277 : & MTABX,KINC,2,IXV,IZTG2,IZTGG2,NELG,NC2,IKG,IRET) 278 : IF(IRET.EQ.0)RETURN 279 : 280 : IF(IKG.EQ.2)IKG=1 281 : 282 : 283 : IF(IARG.EQ.4)THEN 284 : 285 : IXV(1)=MELEMS 286 : IXV(2)=0 287 : IXV(3)=0 287 : CALL LEKCOF('Opérateur NS :', 288 : & MTABX,KINC,3,IXV,IZTG3,IZTGG3,NPT3,NC3,IK3,IRET) 289 : IF(IRET.EQ.0)RETURN 290 : 291 : IXV(1)=MELEMS 292 : IXV(2)=1 293 : IXV(3)=0 293 : CALL LEKCOF('Opérateur NS :', 294 : & MTABX,KINC,4,IXV,IZTG4,IZTGG4,NPT4,NC4,IK4,IRET) 295 : IF(IRET.EQ.0)RETURN 296 : 297 : ENDIF 298 : ENDIF 299 : 300 : C--Cas compréssible 301 : ELSEIF(IKOMP.EQ.1)THEN 302 : C 1er coef : mu , viscosité dynamique 303 : IXV(1)=MELEMC 304 : IXV(2)=1 305 : IXV(3)=0 305 : CALL LEKCOF('Opérateur NS :', 306 : & MTABX,KINC,1,IXV,IZTG1,VISCO,NPT1,NC1,IK1,IRET) 307 : IF(IRET.EQ.0)RETURN 308 : 309 : IK2=-1 310 : IKG=-1 311 : IK3=-1 312 : IK4=-1 313 : IZTGG2=VISCO 314 : IZTGG3=VISCO 315 : IZTGG4=VISCO 316 : 317 : IF(IARG.EQ.3)THEN 318 : C 3ème coef : rog 319 : IXV(1)=-MELEMC 320 : IXV(2)=0 321 : IXV(3)=1 321 : CALL LEKCOF('Opérateur NS :', 322 : & MTABX,KINC,3,IXV,IZTG2,IZTGG2,NELG,NC2,IKG,IRET) 323 : IF(IRET.EQ.0)RETURN 324 : 325 : IF(IKG.EQ.2)IKG=1 326 : 327 : ENDIF 328 : 329 : C 2ème coef : un , champ de vitesse transportant 330 : IXV(1)=-MELEMS 331 : IXV(2)=0 332 : IXV(3)=0 332 : CALL LEKCOF('Opérateur NS :', 333 : & MTABX,KINC,2,IXV,MUTRAN,UTRANS,NPTU,NC3,IKW,IRET) 334 : IF(IRET.EQ.0)RETURN 335 : VITESS=UTRANS 336 : IPADU=IPADS 337 : IPADF=IPADS 338 : 339 : ENDIF 340 : 341 : IF(KKALE.NE.0)THEN 342 : 343 : IF(IKOMP.EQ.1)THEN 344 : WRITE(6,*)' Operateur NS : option compréssible en ALE' 345 : WRITE(6,*)' Cas non prévu pour l''instant: ' 346 : C Option %m1:8 incompatible avec les données 347 : MOTERR( 1: 8) = ' EF ' 347 : CALL ERREUR(803) 348 : RETURN 349 : ENDIF 350 : 351 : C n+ unième coefficient champ de vitesse transportant pour l'option ALE 352 : IXV(1)=-MELEMS 353 : IXV(2)=0 354 : IXV(3)=0 354 : CALL LEKCOF('Opérateur NS :', 355 : & MTABX,KINC,(IARG+1),IXV,MUTRAN,UTRANS,NPTW,NCW,IKW,IRET) 356 : IF(IRET.EQ.0)RETURN 357 : IPADU=IPADS 358 : 359 : ENDIF 360 : 361 : C write(6,*)' Operateur NS : Fin lecture Arguments ' 362 : C Fin lecture Arguments ************************************************ 363 : 363 : IF(IPADS.EQ.0)CALL KRIPAD(MELEMS,IPADS) 364 : 364 : CALL LEKTAB(MTABZ,'XXDXDY',MCHPOI) 365 : IF(MCHPOI.EQ.0)GO TO 90 365 : CALL LICHT(MCHPOI,IZTCO,TYPC,IGEOM) 366 : NELZ=IZTCO.VPOCHA(/1) 367 : 367 : CALL LEKTAB(MTAB1,'PASDETPS',MTABT) 368 : IF(MTABT.EQ.0)THEN 368 : CALL CRTABL(MTABT) 368 : CALL ECMM(MTABT,'SOUSTYPE','PASDETPS') 368 : CALL ECMO(MTAB1,'PASDETPS','TABLE ',MTABT) 369 : DT=1.E30 370 : DTP=1.E30+DT 371 : IPAT=1 371 : CALL ECME(MTABT,'NUPASDT',IPAT) 372 : ELSE 372 : CALL ACMF(MTABT,'DELTAT',DTP) 373 : ENDIF 374 : 375 : C*********** FORMULATIONS ********** 376 : 377 : IF(KFORM.EQ.0)THEN 378 : C Formulation EFM1 379 : IF(KIMPL.NE.0)THEN 380 : WRITE(6,*)' Operateur NS ' 381 : C Option %m1:8 incompatible avec les données 382 : MOTERR( 1: 8) = ' EF ' 382 : CALL ERREUR(803) 383 : RETURN 384 : ENDIF 385 : 385 : CALL LEKTAB(MTABZ,'MATESI',MATRIK) 386 : IF(MATRIK.EQ.0)GO TO 90 387 : SEGACT MATRIK 388 : 389 : IMATRI=IRIGEL(4,1) 390 : SEGACT IMATRI 391 : 391 : CALL LEKTAB(MTABZ,'XXPSOML',MCHELM) 392 : IF(MCHELM.EQ.0)GO TO 90 393 : SEGACT MCHELM 394 : 394 : CALL LEKTAB(MTABZ,'XXVOLUM',MCHVOL) 395 : IF(MCHVOL.EQ.0)GO TO 90 395 : CALL LICHT(MCHVOL,IZVOL,TYPC,IGEOM) 396 : 396 : CALL LEKTAB(MTAB1,'KIZG',KIZG) 397 : IF(KIZG.EQ.0)THEN 397 : CALL CRTABL(KIZG) 397 : CALL ECMM(KIZG,'SOUSTYPE','KIZG') 397 : CALL ECMO(MTAB1,'KIZG','TABLE ',KIZG) 398 : ENDIF 399 : 400 : TYPE=' ' 400 : CALL ACMO(KIZG,NOMI,TYPE,IZG1) 401 : IF(TYPE.NE.'CHPOINT ')THEN 402 : NC=IZTU1.VPOCHA(/2) 403 : TYPE='SOMMET' 403 : CALL CRCHPT(TYPE,MELEM1,NC,IZG1) 403 : CALL ECMO(KIZG,NOMI,'CHPOINT ',IZG1) 404 : ENDIF 405 : 405 : CALL LICHT(IZG1,IZGG1,TYPC,IGEOM) 406 : 407 : IF(IGEOM.NE.MELEM1)THEN 408 : WRITE(6,*)' Opérateur NS' 409 : WRITE(6,*)'Incompatibilite de SPG entre 1eres inconnues' 410 : C Option %m1:8 incompatible avec les données 411 : MOTERR( 1: 8) = ' EF ' 411 : CALL ERREUR(803) 412 : RETURN 413 : ENDIF 414 : 415 : NPT=IZGG1.VPOCHA(/1) 416 : 417 : SEGACT MELEME 418 : NBSOUS=LISOUS(/1) 419 : IF(NBSOUS.EQ.0)NBSOUS=1 420 : NUTOEL=0 421 : DT=1.E30 422 : 423 : 424 : DO 1 L=1,NBSOUS 425 : IPT1=MELEME 426 : IF(NBSOUS.NE.1)IPT1=LISOUS(L) 427 : SEGACT IPT1 428 : IZAFM=LIZAFM(L,1) 429 : IPM1=IZAFM 430 : SEGACT IZAFM 431 : IF(IAXI.NE.0)THEN 432 : IPM1=LIZAFM(L,2) 433 : SEGACT IPM1 434 : ENDIF 435 : 436 : NOM0=NOMS(IPT1.ITYPEL)//' ' 436 : CALL KALPBG(NOM0,'FONFORM0',IZFFM) 437 : SEGACT IZFFM*MOD 438 : 439 : MCHAML=ICHAML(L) 440 : SEGACT MCHAML 441 : MELVAL=IELVAL(1) 442 : SEGACT MELVAL 443 : 444 : IF(IMACHE(L).NE.IPT1)THEN 445 : write(*,*)'IPT1,IMACHE ',IPT1,IMACHE(L) 446 : goto 90 447 : ENDIF 448 : 449 : NP =IPT1.NUM(/1) 450 : NBEL=IPT1.NUM(/2) 451 : IES=IDIM 452 : 453 : NPTF=VITESS.VPOCHA(/1) 454 : NPTU=UTRANS.VPOCHA(/1) 454 : CALL ZCVI(AM,IPM1.AM,VELCHE, 455 : & IPT1.NUM,NBEL,NUTOEL,NPT,IES,NP,IAXI,IKOMP,IARG, 456 : & VISCO.VPOCHA,IK1,IZTGG2.VPOCHA,IKG,NELG,IZTGG3.VPOCHA,IK3, 457 : & IZTGG4.VPOCHA,IK4,IPADS.LECT, 458 : & UTRANS.VPOCHA,IPADU.LECT,NPTU, 459 : & IZTU1.VPOCHA,IZGG1.VPOCHA,IPADI.LECT, 460 : & VITESS.VPOCHA,IPADF.LECT,NPTF, 461 : & IZVOL.VPOCHA,IZTCO.VPOCHA,NELZ,IDCEN,IPG, 462 : & DTM1,DT,DTT1,DTT2,NUEL,DIAEL,IZFFM.FN) 463 : 464 : SEGDES IZAFM,IPT1 465 : IF(IAXI.NE.0)SEGDES IPM1 466 : NUTOEL=NUTOEL+NBEL 467 : 468 : 1 CONTINUE 469 : SEGDES MELEME 470 : 471 : 472 : IF(DT.LT.DTP)THEN 472 : CALL ECMF(MTABT,'DELTAT',DT) 472 : CALL ECMM(MTABT,'OPER','NS') 472 : CALL ECMM(MTABT,'ZONE',NOMZ) 472 : CALL ECMF(MTABT,'DTCONV',DTT1) 472 : CALL ECMF(MTABT,'DTDIFU',DTT2) 472 : CALL ECMF(MTABT,'DIAEL',DIAEL) 472 : CALL ECME(MTABT,'NUEL',NUEL) 473 : ENDIF 474 : 475 : SEGDES VISCO,IZTGG2,IZTGG3,IZTGG4 476 : SEGDES UTRANS,IZTU1 477 : SEGDES IZGG1 478 : SEGDES IZVOL,IZTCO 479 : SEGDES LINCO 480 : SEGDES MATRIK,IMATRI 481 : 482 : 483 : C************************************************************************* 484 : ELSE IF(KFORM.EQ.1)THEN 485 : C CAS FORMULATION EF 486 : 487 : IF(KIMPL.EQ.0)THEN 488 : C Option %m1:8 incompatible avec les données 489 : MOTERR( 1: 8) = ' EF ' 489 : CALL ERREUR(803) 490 : RETURN 491 : ENDIF 492 : 493 : DT=0.D0 494 : IF(IDCEN.EQ.4)THEN 495 : TYPE=' ' 495 : CALL ACMO(MTABT,'DELTAT',TYPE,IENT) 496 : IF(TYPE.NE.'ENTIER')THEN 497 : WRITE(6,*)' Opérateur NS ' 498 : WRITE(6,*)' On reclame un pas de temps' 499 : C Option %m1:8 incompatible avec les données 500 : MOTERR( 1: 8) = ' EF ' 500 : CALL ERREUR(803) 501 : RETURN 502 : ELSE 502 : CALL ACMF(MTABT,'DELTAT',DT) 503 : ENDIF 504 : ENDIF 505 : 506 : IHV=1 507 : NUTOEL=0 508 : NINKO=VITESS.VPOCHA(/2) 509 : SEGACT MELEME 510 : NBSOUS=LISOUS(/1) 511 : IF(NBSOUS.EQ.0)NBSOUS=1 512 : 513 : TYPE=' ' 513 : CALL ACMO(MTABX,'MATELM',TYPE,MATRIK) 514 : IF(TYPE.EQ.'MATRIK'.AND.KMACO.NE.0)THEN 515 : SEGACT MATRIK 516 : NMATRI=IRIGEL(/2) 517 : MELEME=IRIGEL(1,1) 518 : SEGACT MELEME 519 : IMATRI=IRIGEL(4,1) 520 : SEGACT IMATRI 521 : NBME=LIZAFM(/2) 522 : NINKO=NBME 523 : MELEMS=KSPGP 524 : 525 : 526 : ELSE 527 : 528 : NRIGE=7 529 : NKID =9 530 : NKMT =7 531 : NMATRI=1 532 : SEGINI MATRIK 533 : IRIGEL(1,1)=MELEME 534 : IRIGEL(2,1)=MELEME 535 : 536 : IRIGEL(7,1)=2 537 : NBOP=0 538 : NBME=NINKO 539 : NBELC=0 540 : SEGINI IMATRI,IMATRS 541 : IRIGEL(4,1)=IMATRI 542 : KSPGP=MELEMS 543 : KSPGD=MELEMS 544 : 545 : IF(NBME.EQ.1)THEN 546 : LISPRI(1)=NOMI(1:4)//' ' 547 : LISDUA(1)=NOMA(1:4)//' ' 548 : ELSE 549 : DO 102 I=1,NBME 550 : WRITE(NOM,FMT='(I1,A7)')I,NOMI(1:7) 551 : LISPRI(I)=NOM(1:4)//' ' 552 : WRITE(NOM,FMT='(I1,A7)')I,NOMA(1:7) 553 : LISDUA(I)=NOM(1:4)//' ' 554 : 102 CONTINUE 555 : ENDIF 556 : 557 : NUTOEL=0 558 : DO 101 L=1,NBSOUS 559 : IPT1=MELEME 560 : IF(NBSOUS.NE.1)IPT1=LISOUS(L) 561 : SEGACT IPT1 562 : NOM0=NOMS(IPT1.ITYPEL)//' ' 562 : CALL KALPBG(NOM0,'FONFORM ',IZFFM) 563 : SEGACT IZFFM*MOD 564 : IZHR=KZHR(1) 565 : SEGACT IZHR*MOD 566 : NES=GR(/1) 567 : NPG=GR(/3) 568 : 569 : NP = IPT1.NUM(/1) 570 : MP = NP 571 : NBEL=IPT1.NUM(/2) 572 : SEGINI IPM1,IPS1 573 : LIZAFM(L,1)=IPM1 574 : LIZAFS(L,1)=IPS1 575 : IPM2=IPM1 576 : IPM3=IPM1 577 : IPS2=IPS1 578 : IPS3=IPS1 579 : IF(NBME.GE.2)THEN 580 : SEGINI IPM2,IPS2 581 : LIZAFM(L,2)=IPM2 582 : LIZAFS(L,2)=IPS2 583 : ENDIF 584 : IF(NBME.GE.3)THEN 585 : SEGINI IPM3,IPS3 586 : LIZAFM(L,3)=IPM3 587 : LIZAFS(L,3)=IPS3 588 : ENDIF 589 : 590 : 591 : KITT=2 592 : KJTT=IK1 593 : NPT=UTRANS.VPOCHA(/1) 594 : NPTG=IZTU1.VPOCHA(/1) 595 : 595 : CALL ZCONV(FN,GR,PG,XYZ,HR,PGSQ,RPG, 596 : & NES,IDIM,NP,NPG,IAXI,AIMPL,IKOMP, 597 : & RO,1,UTRANS.VPOCHA,IKW,NPT,IPADU.LECT,VISCO.VPOCHA,IK1, 598 : & IPT1.NUM,NBEL,NUTOEL,XCOOR, 599 : & IPM1.AM,IPM2.AM,IPM3.AM, 600 : & IPS1.AM,IPS2.AM,IPS3.AM, 601 : & NINKO,IDCEN,DT, 602 : & IZTU1.VPOCHA,0,NPTG,IPADI.LECT,IZTCO.VPOCHA,NELZ) 603 : 603 : CALL XLAPL(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI, 604 : & VISCO.VPOCHA,VISCO.VPOCHA,VISCO.VPOCHA,KITT,KJTT,IK1, 605 : & IPT1.NUM,NBEL,NUTOEL,XCOOR,AIMPL,IKOMP, 606 : & IPM1.AM,IPM2.AM,IPM3.AM, 607 : & IPS1.AM,IPS2.AM,IPS3.AM, 608 : & NINKO,IHV,IARG,VISCO.VPOCHA) 609 : 610 : NUTOEL=NUTOEL+NBEL 611 : 101 CONTINUE 612 : 613 : ENDIF 614 : 615 : 616 : IF(KIMPL.EQ.2.OR.KIMPL.EQ.0.OR.IARG.GT.1)THEN 617 : NAT=2 618 : NSOUPO=1 619 : SEGACT MELEMS 620 : N=MELEMS.NUM(/2) 621 : NC=NINKO 622 : SEGINI MCHPO1,MSOUP1,MPOVA1 623 : MCHPO1.IFOPOI=IFOMOD 624 : MCHPO1.MOCHDE=TITREE 625 : MCHPO1.MTYPOI='SMBR' 626 : MCHPO1.JATTRI(1)=2 627 : MCHPO1.IPCHP(1)=MSOUP1 628 : DO 177 N=1,NINKO 629 : MSOUP1.NOCOMP(N)=LISDUA(N)(1:4) 630 : 177 CONTINUE 631 : MSOUP1.IGEOC=MELEMS 632 : MSOUP1.IPOVAL=MPOVA1 633 : ENDIF 634 : 635 : IF(IARG.EQ.2.OR.IARG.EQ.4)THEN 635 : CALL KRIPAD(MELEMQ,IPADQ) 636 : IF(IARG.EQ.2)THEN 637 : IKAS=2 638 : ELSEIF(IARG.EQ.4)THEN 639 : IKAS=3 640 : ENDIF 641 : 642 : 643 : 644 : IF(MACRO.NE.0.AND.KPRE.EQ.2)MELEMI=MELEME 645 : SEGACT MELEMI 646 : NBSOUS=MELEMI.LISOUS(/1) 647 : IF(NBSOUS.EQ.0)NBSOUS=1 648 : 649 : NUTOEL=0 650 : DO 1102 L=1,NBSOUS 651 : IPT1=MELEMI 652 : IF(NBSOUS.NE.1)IPT1=MELEMI.LISOUS(L) 653 : SEGACT IPT1 654 : 655 : IF(MQUAD.NE.0)THEN 656 : IF(KPRE.EQ.2)NOM0=NOMS(IPT1.ITYPEL)//'PRP0' 657 : IF(KPRE.EQ.3)NOM0=NOMS(IPT1.ITYPEL)//'PRP0' 658 : IF(KPRE.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'PRP1' 659 : ELSEIF(MACRO.NE.0)THEN 660 : IF(KPRE.EQ.2)NOM0=NOMS(IPT1.ITYPEL)//' ' 661 : IF(KPRE.EQ.3)NOM0=NOMS(IPT1.ITYPEL)//'MCP0' 662 : IF(KPRE.EQ.4)NOM0=NOMS(IPT1.ITYPEL)//'MCP1' 663 : ELSE 664 : IF(KPRE.EQ.2)NOM0=NOMS(IPT1.ITYPEL)//' ' 665 : ENDIF 665 : CALL KALPBG(NOM0,'FONFORM ',IZFFM) 666 : 667 : 668 : SEGACT IZFFM*MOD 669 : IZHR=KZHR(1) 670 : SEGACT IZHR*MOD 671 : NES=GR(/1) 672 : NPG=GR(/3) 673 : IZF1=KTP(1) 674 : SEGACT IZF1*MOD 675 : MP1=IZF1.FN(/1) 676 : NP = IPT1.NUM(/1) 677 : MP = NP 678 : NBEL=IPT1.NUM(/2) 679 : NELG=IZTGG2.VPOCHA(/1) 680 : NPT =MPOVA1.VPOCHA(/1) 681 : 682 : SEGACT MELEP1 683 : 684 : 684 : CALL XSOUR(FN,IZF1.FN,GR,PG,XYZ,HR,PGSQ,RPG, 685 : & NES,IDIM,NP,MP1,NPG,IAXI,IPT1.NUM,IKAS,KPRE, 686 : & IZTGG2.VPOCHA,IKG,NELG,IPADQ.LECT,MELEP1.NUM, 687 : & IZTGG3.VPOCHA,IK3,IZTGG4.VPOCHA,IK4,IPADS.LECT, 688 : & NBEL,NUTOEL,XCOOR,MPOVA1.VPOCHA,NPT) 689 : 690 : NUTOEL=NUTOEL+NBEL 691 : 1102 CONTINUE 692 : 693 : ENDIF 694 : 695 : IF(KIMPL.EQ.2.OR.KIMPL.EQ.0)THEN 696 : 697 : NBSOUS=LISOUS(/1) 698 : IF(NBSOUS.EQ.0)NBSOUS=1 699 : 700 : DO 1533 L=1,NBSOUS 701 : IPT1=MELEME 702 : IF(NBSOUS.NE.1)IPT1=LISOUS(L) 703 : SEGACT IPT1 704 : NP=IPT1.NUM(/1) 705 : NBEL=IPT1.NUM(/2) 706 : DO 2 N=1,NINKO 707 : IPMS=LIZAFS(L,N) 708 : SEGACT IPMS 709 : DO 12 K=1,NBEL 710 : DO 13 J=1,NP 711 : UU=0.D0 712 : IU=IPADS.LECT(IPT1.NUM(J,K)) 713 : DO 14 I=1,NP 714 : IK=IPADI.LECT(IPT1.NUM(I,K)) 715 : UU=UU+IPMS.AM(K,I,J)*IZTU1.VPOCHA(IK,N) 716 : 14 CONTINUE 717 : MPOVA1.VPOCHA(IU,N)=MPOVA1.VPOCHA(IU,N)+UU 718 : 13 CONTINUE 719 : 12 CONTINUE 720 : 721 : 2 CONTINUE 722 : 723 : 1533 CONTINUE 724 : ENDIF 725 : 726 : SEGDES IPM1,IPM2,IPM3 727 : IPS=IPS1 728 : SEGSUP IPS1 729 : IF(IPS2.NE.IPS)SEGSUP IPS2 730 : IF(IPS3.NE.IPS)SEGSUP IPS3 731 : SEGDES IZTCO 732 : 733 : IF(KIMPL.EQ.2.OR.KIMPL.EQ.0.OR.IARG.GT.1)THEN 734 : 735 : TYPE=' ' 735 : CALL ACMO(MTAB1,'SMBR',TYPE,MCHPO2) 736 : IF(TYPE.NE.'CHPOINT')THEN 736 : CALL ECMO(MTAB1,'SMBR','CHPOINT',MCHPO1) 737 : ELSE 737 : CALL ECROBJ('CHPOINT',MCHPO2) 737 : CALL ECROBJ('CHPOINT',MCHPO1) 737 : C? CALL OPERAD 737 : CALL PRFUSE 737 : CALL LIROBJ('CHPOINT',MCHPOI,1,IRET) 737 : CALL ECMO(MTAB1,'SMBR','CHPOINT',MCHPOI) 738 : ENDIF 739 : 740 : ENDIF 741 : 742 : SEGDES IMATRI 743 : SEGDES MELEME,MATRIK 744 : IF(IK1.EQ.0)THEN 745 : SEGDES VISCO 746 : ENDIF 746 : CALL ECMO(MTABX,'MATELM','MATRIK',MATRIK) 747 : 748 : C************************************************************************* 749 : ELSE IF(KFORM.EQ.2)THEN 750 : C CAS FORMULATION VF 751 : WRITE(6,*)' FORMULATION VF NON DISPONIBLE ' 752 : ENDIF 753 : C************************************************************************* 754 : 755 : IPDI=IPADI 756 : SEGSUP IPADI 757 : IF(IPADS.NE.IPDI)SEGSUP IPADS 758 : IPAS=1 759 : C write(6,*)' FIN de NS' 760 : RETURN 761 : 90 CONTINUE 762 : WRITE(6,*)' Interuption anormale de NS' 763 : C Option %m1:8 incompatible avec les données 764 : MOTERR( 1: 8) = ' EF ' 764 : CALL ERREUR(803) 765 : RETURN 766 : 1002 FORMAT(10(1X,1PE11.4)) 767 : 1001 FORMAT(20(1X,I5)) 768 : END 769 : 770 : 771 : 772 : 773 : 774 : 775 : 776 : 777 : 778 : 779 :
© Cast3M 2003 - All rights reserved.
Disclaimer