1 : C DIVC SOURCE MAGN 10/05/18 21:15:17 6675 2 : SUBROUTINE DIVC(IKAS,MTABX) 3 : IMPLICIT INTEGER(I-N) 4 : IMPLICIT REAL*8 (A-H,O-Z) 5 : C************************************************************************ 6 : C Operateur KMAC 7 : C 8 : C OBJET : Cree un objet de type MATRIK 9 : C 10 : C 11 : C IKAS=1 KMAB calcul de B (DIV U) 12 : C IKAS=2 KMBT calcul de Bt uniquement (GRAD P) 13 : C IKAS=3 KBBT calcul de B assemblage pour B et Bt 14 : C 15 : C*********************************************************************** 16 : -INC CCOPTIO 17 : -INC CCGEOME 18 : -INC CCREEL 19 : -INC SMCOORD 20 : -INC SIZFFB 21 : -INC SMCHPOI 22 : POINTEUR IZCH2.MCHPOI,IZCCH2.MPOVAL 23 : POINTEUR IZDV.MCHPOI,IZDDV.MPOVAL,IZTU1.MPOVAL,TETAN.MPOVAL 24 : POINTEUR IZTG1.MCHPOI,IZTGG1.MPOVAL,IZBETA.MPOVAL,VISCO.MPOVAL 25 : -INC SMMATRIK 26 : POINTEUR IPM.IZAFM 27 : SEGMENT IMATRS 28 : ENDSEGMENT 29 : POINTEUR IPMS.IZAFM,IPS1.IZAFM,IPS2.IZAFM,IPS3.IZAFM 30 : 31 : PARAMETER (LRV=64) 32 : DIMENSION VISC(3) 33 : 34 : -INC SMLENTI 35 : POINTEUR IZIPAD.MLENTI,MLENTI1.MLENTI,MLENTI2.MLENTI 36 : -INC SMLMOTS 37 : POINTEUR LINCO.MLMOTS 38 : -INC SMELEME 39 : POINTEUR MELEMZ.MELEME,MELEMB.MELEME,MELSTB.MELEME 40 : POINTEUR MELEM1.MELEME,MELES1.MELEME,MCTREI.MELEME 41 : POINTEUR IGEOM.MELEME,MELEMA.MELEME 42 : POINTEUR IZLEMC.MELEME,MELEMS.MELEME,MELEMC.MELEME 43 : POINTEUR MELEMI.MELEME,MELEMP.MELEME 44 : POINTEUR MACRO1.MELEME 45 : 46 : CHARACTER*8 TYPE,TYPC,NOMZ,NOMP,NOMD,NOM0 47 : CHARACTER*8 NOMPP,NOMDD 48 : CHARACTER*4 NOM 49 : INTEGER IPAD,IPAD2,IK 50 : REAL*8 KAUX,TETA1 51 : DIMENSION IXV(5) 52 : C 53 : C************************************************************************* 54 : CKMIC 55 : C write(6,*)' Operateur KMIC MTABX=',MTABX 56 : C 57 : C- Récupération de la table EQEX (pointeur MTAB1) 58 : C 58 : CALL LEKTAB(MTABX,'EQEX',MTAB1) 59 : IF(MTAB1.EQ.0)THEN 60 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 61 : MOTERR( 1: 8) = ' EQEX ' 62 : MOTERR( 9:16) = ' EQEX ' 63 : MOTERR(17:24) = ' KIZX ' 63 : CALL ERREUR(786) 64 : RETURN 65 : ENDIF 66 : 66 : CALL ACME(MTAB1,'DISCPRES',IDISCP) 67 : 67 : CALL ACME(MTAB1,'NAVISTOK',NASTOK) 68 : IF(NASTOK.EQ.0)THEN 68 : CALL ZKMIC(IKAS,MTABX,MTAB1) 69 : RETURN 70 : ENDIF 71 : C 72 : C- Récupération de la table INCO (pointeur KINC) 73 : C 73 : CALL LEKTAB(MTAB1,'INCO',KINC) 74 : IF(KINC.EQ.0)THEN 75 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 76 : MOTERR( 1: 8) = ' INCO ' 77 : MOTERR( 9:16) = ' INCO ' 78 : MOTERR(17:24) = ' EQEX ' 78 : CALL ERREUR(786) 79 : RETURN 80 : ENDIF 81 : 82 : C************************************************************************* 83 : C OPTIONS 84 : C KIMPL = 0 -> EXPL 1 -> IMPL 2 -> CN 85 : C KFORM = 0 -> EFSI 1 -> EF 2 -> VF 3 -> EFMC 86 : C KPRE=3 pression P0 KPRE=4 pression P1 KPRE=2 cas macro 1ère génération 87 : C KPRE=5 pression MSOMMET 88 : 89 : IAXI=0 90 : IK=0 91 : IF(IFOMOD.EQ.0)IAXI=2 92 : C 93 : C- Récupération de la table des options KOPT (pointeur KOPTI) 94 : C 94 : CALL LEKTAB(MTABX,'KOPT',KOPTI) 95 : IF (KOPTI.EQ.0) THEN 96 : C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24 97 : MOTERR( 1: 8) = ' KOPT ' 98 : MOTERR( 9:16) = ' KOPT ' 99 : MOTERR(17:24) = ' KIZX ' 99 : CALL ERREUR(786) 100 : RETURN 101 : ENDIF 102 : 102 : CALL ACME(KOPTI,'KIMPL',KIMPL) 102 : CALL ACME(KOPTI,'KPOIN',KPRE) 103 : C write(6,*)'KBBT KPRE=',KPRE 103 : CALL ACMF(KOPTI,'AIMPL',TETA1) 103 : CALL ACME(KOPTI,'KFORM',KFORM) 103 : C CALL ACME(KOPTI,'IKOMP',IKOMP) 104 : 105 : IF (IERR.NE.0) RETURN 106 : C write(6,*)' Apres les options ' 107 : C************************************************************************* 108 : C 109 : C- Récupération de la table DOMAINE associée au domaine local 110 : C 110 : CALL ACMM(MTABX,'NOMZONE',NOMZ) 111 : TYPE=' ' 111 : CALL ACMO(MTABX,'DOMZ',TYPE,MMODEL) 112 : IF(TYPE.NE.'MMODEL')THEN 113 : C On attend un des objets : %m1:8 %m9:16 %m17:24 %m25:32 %m33:40 114 : MOTERR( 1: 8) = ' MMODEL ' 115 : MOTERR( 8:16) = ' MMODEL ' 116 : MOTERR(17:24) = ' MMODEL ' 117 : MOTERR(25:32) = ' MMODEL ' 118 : MOTERR(33:40) = ' MMODEL ' 118 : CALL ERREUR(471) 119 : RETURN 120 : ENDIF 121 : 121 : CALL LEKMOD(MMODEL,MTABZ,INEFMD) 122 : 122 : CALL LEKTAB(MTABZ,'MAILLAGE',MELEME) 122 : CALL LEKTAB(MTABZ,'SOMMET',MELEMS) 122 : CALL LEKTAB(MTABZ,'CENTRE',MELEMC) 122 : C?? CALL LEKTAB(MTABZ,'MACRO',MACRO) 123 : MACRO1=0 123 : IF(INEFMD.EQ.2)CALL LEKTAB(MTABZ,'MACRO1',MACRO1) 123 : C?? IF(INEFMD.EQ.3)CALL LEKTAB(MTABZ,'QUADRATI',MQUAD) 124 : IF (IERR.NE.0) RETURN 125 : IF(INEFMD.EQ.4.AND.KPRE.NE.5)THEN 126 : C% Données incompatibles 126 : CALL ERREUR(21) 127 : RETURN 128 : ENDIF 129 : 130 : MELEMI=MELEME 131 : IF(MACRO1.NE.0.AND.KPRE.NE.2)THEN 132 : MELEMI=MACRO1 133 : ENDIF 134 : 135 : IF(KPRE.EQ.2.AND.INEFMD.EQ.3)KPRE=3 136 : IF(INEFMD.EQ.1.AND.KPRE.NE.5)KPRE=2 137 : 138 : IF(KPRE.EQ.3)THEN 138 : CALL LEKTAB(MTABZ,'CENTREP0',MELEMC) 139 : MELEMP=MELEMC 140 : ELSEIF(KPRE.EQ.4)THEN 140 : CALL LEKTAB(MTABZ,'CENTREP1',MELEMC) 140 : CALL LEKTAB(MTABZ,'ELTP1NC ',MELEMP) 141 : ELSEIF(KPRE.EQ.2)THEN 141 : CALL LEKTAB(MTABZ,'CENTRE',MELEMC) 142 : MELEMP=MELEMC 143 : ELSEIF(KPRE.EQ.5)THEN 143 : CALL LEKTAB(MTABZ,'MSOMMET',MELEMC) 143 : CALL LEKTAB(MTABZ,'MMAIL ',MELEMP) 144 : ENDIF 145 : C write(6,*)' KMIC KPRE=',kpre 146 : 147 : IF(IDISCP.EQ.0)THEN 147 : CALL ECME(MTAB1,'DISCPRES',KPRE) 148 : ELSEIF(IDISCP.NE.KPRE)THEN 149 : C% Données incompatibles 149 : CALL ERREUR(21) 150 : RETURN 151 : ENDIF 152 : 153 : IKOMP=1 154 : IF(KPRE.EQ.5)IKOMP=0 155 : 156 : C************************************************************************* 157 : C VERIFICATIONS SUR LES INCONNUES 158 : 159 : C write(6,*)' Verification sur les inconnues ' 160 : TYPE='LISTMOTS' 160 : CALL ACMO(MTABX,'LISTINCO',TYPE,LINCO) 161 : IF (IERR.NE.0) RETURN 162 : SEGACT LINCO 163 : NBINC=LINCO.MOTS(/2) 164 : IF(NBINC.NE.2)THEN 165 : WRITE(6,*)'Operateur KMAB ' 166 : WRITE(6,*)'Nombre d''inconnues incorrecte : ',NBINC,' On attend 2' 167 : C Indice %m1:8 : contient plus de %i1 %m9:16 168 : MOTERR( 1:8) = 'LISTINCO' 169 : INTERR(1) = 2 170 : MOTERR(9:16) = ' MOTS ' 170 : CALL ERREUR(799) 171 : RETURN 172 : ENDIF 173 : 174 : C On recupere PHI n et TETA n pour Cranck-Nicholson 175 : NOMP=LINCO.MOTS(1) 176 : TYPE=' ' 176 : CALL ACMO(KINC,NOMP,TYPE,MCHPOI) 177 : IF(TYPE.NE.'CHPOINT ')THEN 178 : WRITE(6,*)' Opérateur KMAB :' 179 : WRITE(6,*)' L objet CHPOINT ',NOMP, 180 : & ' n existe pas dans la table' 181 : C Indice %m1:8 : ne contient pas un objet de type %m9:16 182 : MOTERR( 1: 8) = 'INC '//NOMP 183 : MOTERR( 9:16) = 'CHPOINT ' 183 : CALL ERREUR(800) 184 : RETURN 185 : ELSE 185 : CALL LICHT(MCHPOI,IZTU1,TYPC,IGEOM0) 186 : ENDIF 187 : C************************************************************************* 188 : C Le domaine de definition est donne par le SPG de la premiere inconnue 189 : C Les inconnues suivantes devront posseder ce meme pointeur 190 : C On verifie que les points de la zone sont tous inclus dans ce SPG 191 : C Inconnue Primale 192 : 193 : C write(6,*)' Verification inconnue primale ' 193 : CALL KRIPAD(IGEOM0,MLENTI) 194 : IF(IKAS.EQ.1.OR.IKAS.EQ.3)THEN 195 : MELEMK=MELEMS 196 : ELSE 197 : MELEMK=MELEMC 198 : ENDIF 199 : 200 : C write(6,*)' MELEMK=',melemk 200 : CALL VERPAD(MLENTI,MELEMK,IRET) 201 : IF(IRET.NE.0)THEN 202 : WRITE(6,*)' Opérateur KMAB ' 203 : WRITE(6,*)' La zone ',NOMZ,' n''est pas incluse dans le domaine' 204 : & , ' de définition de l''inconnue ',NOMP 205 : WRITE(6,*)' MELEMK=',melemk,' IGEOM0=',IGEOM0 206 : C Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique 207 : MOTERR(1: 8) = 'INC '//NOMP 208 : MOTERR(9:16) = 'CHPOINT ' 208 : CALL ERREUR(788) 209 : IPAS=0 210 : RETURN 211 : ENDIF 212 : 213 : SEGSUP MLENTI 214 : 215 : C************************************************************************* 216 : 217 : NOMD=LINCO.MOTS(2) 218 : TYPE=' ' 218 : CALL ACMO(KINC,NOMD,TYPE,MCHPOI) 219 : IF(TYPE.NE.'CHPOINT ')THEN 220 : WRITE(6,*)' Opérateur KMAB :' 221 : WRITE(6,*)' L objet CHPOINT ',NOMD, 222 : & ' n existe pas dans la table' 223 : C Indice %m1:8 : ne contient pas un objet de type %m9:16 224 : MOTERR( 1: 8) = 'INC '//NOMD 225 : MOTERR( 9:16) = 'CHPOINT ' 225 : CALL ERREUR(800) 226 : RETURN 227 : ELSE 227 : CALL LICHT(MCHPOI,TETAN,TYPC,IGEOM0) 228 : ENDIF 229 : 230 : NC=TETAN.VPOCHA(/2) 231 : C************************************************************************* 232 : C Le domaine de definition est donne par le SPG de la premiere inconnue 233 : C Les inconnues suivantes devront posseder ce meme pointeur 234 : C On verifie que les points de la zone sont tous inclus dans ce SPG 235 : C Inconnue Duale 236 : 237 : C write(6,*)' IGEOM0=',igeom0 237 : CALL KRIPAD(IGEOM0,MLENTI) 238 : IF(IKAS.EQ.1.OR.IKAS.EQ.3)THEN 239 : MELEMK=MELEMC 240 : ELSE 241 : MELEMK=MELEMS 242 : ENDIF 243 : 244 : C write(6,*)' Verification inconnue duale ',MELEMK 244 : CALL VERPAD(MLENTI,MELEMK,IRET) 245 : IF(IRET.NE.0)THEN 246 : WRITE(6,*)' Opérateur KMAB ' 247 : WRITE(6,*)' La zone ',NOMZ,' n''est pas incluse dans le domaine' 248 : & ,' de définition de l''inconnue ',NOMD 249 : WRITE(6,*)' MELEMK=',melemk,' IGEOM0=',IGEOM0 250 : C Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique 251 : MOTERR(1: 8) = 'INC '//NOMD 252 : MOTERR(9:16) = 'CHPOINT ' 252 : CALL ERREUR(788) 253 : IPAS=0 254 : RETURN 255 : ENDIF 256 : 257 : SEGSUP MLENTI 258 : 259 : C************************************************************************* 260 : C Lecture du ou des coefficients 261 : C Type du coefficient : 262 : C IK1=0 CHPOINT IK1=1 scalaire IK1=2 vecteur 263 : 264 : C write(6,*)' Verification sur les coefficients ' 264 : CALL ACME(MTABX,'IARG',IARG) 265 : C write(6,*)' IARG=',iarg 266 : 267 : IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)THEN 268 : TYPE=' ' 268 : CALL ACMO(MTABZ,'MELSTB',TYPE,MELSTB) 269 : SEGACT MELSTB 270 : IF(IDIM.EQ.2)NBELEM=MELSTB.NUM(/2)/4 271 : IF(IDIM.EQ.3)NBELEM=MELSTB.NUM(/2)/8 272 : NBNN=MELSTB.NUM(/1) 273 : NBSOUS=0 274 : NBREF=0 275 : SEGINI MELEMA 276 : MELEMA.ITYPEL=MELSTB.ITYPEL 277 : 278 : NKPE=4 279 : IF(IDIM.EQ.3)NKPE=8 280 : do 4878 k=1,nbelem 281 : mi=(k-1)*NKPE+1 282 : do 4879 i=1,nbnn 283 : MELEMA.num(i,k)=melstb.num(i,mi) 284 : 4879 continue 285 : C write(6,*)k,(MELEMA.num(i,k),i=1,nbnn) 286 : 4878 continue 287 : 288 : TYPE=' ' 288 : CALL ACMO(MTABZ,'MCHPOC',TYPE,MCHPOC) 289 : TYPE=' ' 289 : CALL ACMO(MTABZ,'CENTRE',TYPE,MCTREI) 290 : ENDIF 291 : 292 : 293 : C 1er COEF 294 : 295 : IXV(1)=MELEMC 296 : IXV(2)=1 297 : IXV(3)=0 298 : IXV(4)=MELEMS 299 : IXV(5)=-MELEMS 300 : IRET=5 300 : CALL LEKCOF('Opérateur KMAB :', 301 : & MTABX,KINC,1,IXV,IZTG1,IZTGG1,NPT1,NC1,IK1,IRET) 302 : IF(IRET.EQ.0)RETURN 302 : IF(IK1.GE.4)CALL KRIPAD(MELEMS,MLENTI) 303 : 304 : 305 : BETA0=1.D-1 306 : IF(IARG.EQ.2)THEN 307 : C 2ème COEF 308 : IXV(1)=0 309 : IXV(2)=1 310 : IXV(3)=0 310 : CALL LEKCOF('Opérateur KMAB :', 311 : & MTABX,KINC,2,IXV,IZTG2,IZBETA,NPT2,NC2,IK2,IRET) 312 : IF(IRET.EQ.0)RETURN 313 : BETA0=IZBETA.VPOCHA(1,1) 314 : C? IF(IZBETA.VPOCHA(1,1).LT.0.D0)THEN 315 : C% Données incompatibles 315 : C? CALL ERREUR(21) 316 : C? RETURN 317 : C? ENDIF 318 : ENDIF 319 : 320 : 321 : NOMP=LINCO.MOTS(1) 322 : NOMD=LINCO.MOTS(2)(1:4) 323 : 324 : NRIGE=7 325 : NKID =9 326 : NKMT =7 327 : NMATRI=1 328 : IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)NMATRI=2 329 : IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.4.AND.IDIM.EQ.2)NMATRI=2 330 : IF(KPRE.EQ.5)NMATRI=NMATRI+1 331 : C write(6,*)' NMATRI=',nmatri,' KPRE=',kpre 332 : SEGINI MATRIK 333 : 334 : C CAS Stabilisation via MACRO CENTRE 335 : IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)THEN 336 : IF(IKAS.EQ.3)BETA0=-BETA0 337 : C write(6,*)' CAS Stabilisation via MACRO CENTRE ',BETA0,IKAS 338 : NK=0 339 : I2=2 340 : NBME=1 341 : NBSOUS=1 342 : SEGINI IMATRI 343 : IRIGEL(4,i2)=IMATRI 344 : KSPGP=MCTREI 345 : KSPGD=MCTREI 346 : IRIGEL(1,i2)=MELEMA 347 : IRIGEL(2,i2)=MELEMA 348 : IRIGEL(7,i2)=0 348 : CALL LICHT(MCHPOC,MPOVAL,TYPC,IGEOM) 349 : 350 : C write(6,*)' MELSTB=',melstb,' MCTREI=',MCTREI 351 : 352 : SEGACT MELSTB 353 : 354 : NBSOUS=MELSTB.LISOUS(/1) 355 : IF(NBSOUS.NE.0)THEN 355 : CALL ERREUR(5) 356 : ENDIF 357 : 358 : NBEL=MELEMA.NUM(/2) 359 : NBCI=MELSTB.NUM(/2) 360 : NP =MELSTB.NUM(/1) 361 : MP =NP 362 : 363 : C write(6,*)' nbel,nbci,np,mp=',nbel,nbci,np,mp 364 : SEGINI IZAFM 365 : LIZAFM(1,1)=IZAFM 366 : LISPRI(1)=NOMD 367 : LISDUA(1)=NOMD 368 : 368 : CALL KRIPAD(MCTREI,MLENT1) 369 : 370 : C write(6,*)' IK1=',ik1 371 : IF(IK1.LT.4)THEN 372 : DO 33 K=1,NBEL 373 : 374 : NK=NK+1 375 : KC=1+(1-IK1)*(NK-1) 376 : C write(6,*)' K=',K,' NK=',nk,np,kc 377 : 378 : DO 32 J=1,NP 379 : K1=MLENT1.LECT(MELEMA.NUM(J,K)) 380 : II=J 381 : DO 321 I=1,NP 382 : U=VPOCHA(K1,I)*BETA0*IZTGG1.VPOCHA(KC,1) 383 : IF(I.EQ.1)U=ABS(VPOCHA(K1,I))*BETA0*IZTGG1.VPOCHA(KC,1) 384 : IF(II.LE.NP)THEN 385 : AM(K,II,J)=U 386 : ELSE 387 : AM(K,II-NP,J)=U 388 : ENDIF 389 : II=II+1 390 : 321 CONTINUE 391 : 32 CONTINUE 392 : 33 CONTINUE 393 : C write(6,*)' FIN BCL 33 ' 394 : ELSE 395 : 396 : DO 34 K=1,NBEL 397 : 398 : NK=NK+1 399 : UA=1.D0 400 : 401 : DO 35 J=1,NP 402 : K1=MLENT1.LECT(MELEMA.NUM(J,K)) 403 : II=J 404 : DO 351 I=1,NP 405 : U=VPOCHA(K1,I)*BETA0*UA 406 : IF(I.EQ.1)U=ABS(VPOCHA(K1,I))*BETA0*UA 407 : IF(II.LE.NP)THEN 408 : AM(K,II,J)=U 409 : ELSE 410 : AM(K,II-NP,J)=U 411 : ENDIF 412 : II=II+1 413 : 351 CONTINUE 414 : 35 CONTINUE 415 : 34 CONTINUE 416 : 417 : 418 : ENDIF 419 : SEGSUP MLENT1 420 : ENDIF 421 : C write(6,*)' Apres 34 ' 422 : 423 : C CAS Stabilisation via MACRO CENTREP1 424 : DEUPI=1.D0 425 : IF(IAXI.NE.0)DEUPI=2.D0*XPI 426 : IF(IKAS.EQ.3)DEUPI=-DEUPI 427 : 428 : IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.4.AND.IDIM.EQ.2)THEN 429 : I2=2 430 : NBME=1 431 : NBSOUS=1 432 : SEGINI IMATRI 433 : IRIGEL(4,i2)=IMATRI 434 : KSPGP=MELEMC 435 : KSPGD=MELEMC 436 : IRIGEL(1,i2)=MELEMP 437 : IRIGEL(2,i2)=MELEMP 438 : IRIGEL(7,i2)=0 439 : 440 : SEGACT MELEMP 441 : NBSOUS=MELEMP.LISOUS(/1) 442 : IF(NBSOUS.NE.0)THEN 442 : CALL ERREUR(5) 443 : ENDIF 444 : 445 : SEGACT MELEMI 446 : NBSOUM=MELEMI.LISOUS(/1) 447 : IF(NBSOUM.EQ.0)NBSOUM=1 448 : 449 : NBEL=MELEMP.NUM(/2) 450 : NP =MELEMP.NUM(/1) 451 : MP =NP 452 : NOM0=NOMS(MELEMP.ITYPEL)//' ' 453 : 454 : SEGINI IZAFM 455 : LIZAFM(1,1)=IZAFM 456 : LISPRI(1)=NOMD 457 : LISDUA(1)=NOMD 458 : 458 : CALL KALPBG(NOM0,'FONFORM ',IZFFM) 459 : 460 : SEGACT IZFFM*MOD 461 : IZHR=KZHR(1) 462 : SEGACT IZHR*MOD 463 : NES=GR(/1) 464 : NPG=GR(/3) 465 : 466 : NK=0 467 : DO 430 LM=1,NBSOUM 468 : IPT1=MELEMI 469 : IF(NBSOUM.NE.1)IPT1=MELEMI.LISOUS(LM) 470 : SEGACT IPT1 471 : ITYP=IPT1.ITYPEL 472 : BETA=0.D0 473 : IF(NOMS(ITYP).EQ.'TRI6')BETA=BETA0 474 : C write(6,*)' NOMS(ITYP)=',NOMS(ITYP),BETA 475 : NBEL=IPT1.NUM(/2) 476 : 477 : C Cas coefficient FLOTTANT ou CENTRE 478 : IF(IK1.LT.4)THEN 479 : DO 43 KK=1,NBEL 480 : NK=NK+1 481 : DO 42 I=1,NP 482 : J=MELEMP.NUM(I,NK) 483 : DO 42 N=1,IDIM 484 : XYZ(N,I)=XCOOR((J-1)*(IDIM+1) +N) 485 : 42 CONTINUE 486 : 486 : CALL CALJBR 487 : &(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE,AJ,SGN) 488 : KC=1+(1-IK1)*(NK-1) 489 : 490 : DO 44 I=1,NP 491 : DO 44 J=1,NP 492 : U=0.D0 493 : DO 433 L=1,NPG 494 : UH=0.D0 495 : DO 4333 KH=1,IDIM 496 : UH=UH+HR(KH,J,L)*HR(KH,I,L) 497 : 4333 CONTINUE 498 : U=U-UH*PGSQ(L)*DEUPI*BETA*IZTGG1.VPOCHA(KC,1) 499 : 433 CONTINUE 500 : AM(NK,I,J)=U 501 : 44 CONTINUE 502 : 503 : 43 CONTINUE 504 : 505 : ELSEIF(IK1.EQ.4.OR.IK1.EQ.5)THEN 506 : C Cas coefficient SOMMET 507 : 508 : DO 45 KK=1,NBEL 509 : NK=NK+1 510 : DO 46 I=1,NP 511 : J=MELEMP.NUM(I,NK) 512 : DO 46 N=1,IDIM 513 : XYZ(N,I)=XCOOR((J-1)*(IDIM+1) +N) 514 : 46 CONTINUE 515 : 515 : CALL CALJBR 516 : &(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE,AJ,SGN) 517 : 518 : DO 47 I=1,NP 519 : DO 47 J=1,NP 520 : U=0.D0 521 : DO 453 L=1,NPG 522 : UH=0.D0 523 : DO 4533 KH=1,IDIM 524 : UH=UH+HR(KH,J,L)*HR(KH,I,L) 525 : 4533 CONTINUE 526 : 527 : UA=0.D0 528 : DO 454 JJ=1,NP 529 : J1=LECT(IPT1.NUM(JJ,NK)) 530 : DO 454 KH=1,IDIM 531 : KG=(IK1-4)*(KH-1)+1 532 : UA=UA+FN(JJ,L)*IZTGG1.VPOCHA(J1,KG) 533 : 454 CONTINUE 534 : UA=UA/KG 535 : 536 : U=U-UH*PGSQ(L)*DEUPI*BETA*UA 537 : 453 CONTINUE 538 : AM(NK,I,J)=U 539 : 47 CONTINUE 540 : 541 : 45 CONTINUE 542 : ENDIF 543 : 544 : SEGDES IPT1 545 : 430 CONTINUE 546 : 547 : ENDIF 548 : C Fin Stabilisation de toutes sortes 549 : C write(6,*)'C Fin Stabilisation de toutes sortes' 550 : 551 : NBME=IDIM 551 : CALL KRIPAD(MELEMI,MLENTI1) 552 : SEGACT MELEMI 553 : NBSOUS=MELEMI.LISOUS(/1) 554 : IF(NBSOUS.EQ.0)NBSOUS=1 555 : SEGINI IMATRI 556 : 557 : IF(IKAS.EQ.2)THEN 558 : KSPGD=MELEMS 559 : KSPGP=MELEMC 560 : IRIGEL(2,1)=MELEMI 561 : IRIGEL(1,1)=MELEMP 562 : ELSE 563 : KSPGP=MELEMS 564 : KSPGD=MELEMC 565 : IRIGEL(1,1)=MELEMI 566 : IRIGEL(2,1)=MELEMP 567 : ENDIF 568 : SEGACT MELEMP 569 : 570 : C write(6,*)' ds kmac melemp=',IRIGEL(1,1) 571 : C write(6,*)' ds kmac melemd=',IRIGEL(2,1) 572 : 573 : IRIGEL(4,1)=IMATRI 574 : * GBM modifie le 22 nov 99 pour pas de multipl de lagrange 575 : IF(IKAS.EQ.1)IRIGEL(7,1)=3 576 : IF(IKAS.EQ.2)IRIGEL(7,1)=3 577 : IF(IKAS.EQ.3)IRIGEL(7,1)=4 578 : 579 : SEGACT MELEMP 580 : NBSOUC=MELEMP.LISOUS(/1) 581 : IF(NBSOUC.EQ.0)NBSOUC=1 582 : C write(6,*)'NBSOUS,NBSOUC=',NBSOUS,NBSOUC 583 : 584 : K0=0 585 : DO 11 L=1,NBSOUS 586 : IPT1=MELEMI 587 : IPT2=MELEMP 588 : IF(NBSOUS.NE.1)IPT1=MELEMI.LISOUS(L) 589 : IF(NBSOUC.NE.1)IPT2=MELEMP.LISOUS(L) 590 : SEGACT IPT1,IPT2 591 : NBEL=IPT1.NUM(/2) 592 : 593 : IF(IKAS.EQ.2)THEN 594 : MP=IPT1.NUM(/1) 595 : C avt NP=MELEMP.NUM(/1) 596 : NP=IPT2.NUM(/1) 597 : ELSE 598 : NP=IPT1.NUM(/1) 599 : C avt MP=MELEMP.NUM(/1) 600 : MP=IPT2.NUM(/1) 601 : ENDIF 602 : 603 : DO 12 I=1,NBME 604 : SEGINI IZAFM 605 : LIZAFM(L,I)=IZAFM 606 : IF(IKAS.EQ.2)THEN 607 : WRITE(NOM,FMT='(I1,A3)')I,NOMD(1:3) 608 : LISDUA(I)=NOM//' ' 609 : LISPRI(I)=NOMP 610 : ELSE 611 : WRITE(NOM,FMT='(I1,A3)')I,NOMP(1:3) 612 : LISPRI(I)=NOM//' ' 613 : LISDUA(I)=NOMD 614 : ENDIF 615 : 12 CONTINUE 616 : IPM1=LIZAFM(L,1) 617 : IPM2=LIZAFM(L,2) 618 : IPM3=LIZAFM(L,2) 619 : IF(IDIM.EQ.3)IPM3=LIZAFM(L,3) 620 : 621 : IF(IK1.LT.4)THEN 622 : C write(6,*)' Appel KPRUSS' 622 : CALL KPRUSS(IPT1,IPT2, 623 : &IPM1,IPM2,IPM3,IAXI,IKAS,INEFMD,KPRE,IZTGG1,IK1,K0,IKOMP) 624 : C write(6,*)'Retour KPRUSS' 625 : ELSE 626 : C write(6,*)' Appel KPRJSS' 626 : CALL KPRJSS(IPT1,IPT2,IPM1,IPM2,IPM3, 627 : & IAXI,IKAS,INEFMD,KPRE,IZTGG1,MLENTI,IK1,IKOMP) 628 : ENDIF 629 : 630 : K0=K0+NBEL 631 : SEGDES IPM1,IPM2,IPM3 632 : SEGDES IPT1 633 : 11 CONTINUE 634 : SEGDES MELEMI 635 : IF(IK1.EQ.4)SEGSUP MLENTI 636 : 637 : 638 : C%%%%%%%%%%%%%%%%% CAS DES PRESSIONS CONTINUES %%%%%%%%%%%%%%%%%%%%%%% 639 : 640 : 641 : IF(KPRE.EQ.5)THEN 641 : CALL LEKTAB(MTABZ,'MMAIL',MELEME) 641 : CALL LEKTAB(MTABZ,'MSOMMET',MELEMS) 642 : 643 : SEGACT MELEME 644 : NBSOU1=LISOUS(/1) 645 : IF(NBSOU1.EQ.0)NBSOU1=1 646 : 647 : IRIGEL(1,2)=MELEME 648 : IRIGEL(2,2)=MELEME 649 : IRIGEL(7,2)=0 650 : 651 : c%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 652 : 653 : IHV=0 654 : AIMPL=1.D0 655 : 656 : C --- CALCUL DU NOMBRE DE PAQUETS DE LRV ELEMENTS --- 656 : CALL KPLRVM(MELEME,LRV,NBSOUS) 657 : NBME=1 658 : 659 : SEGINI IMATRI 660 : IRIGEL(4,2)=IMATRI 661 : KSPGP=MELEMS 662 : KSPGD=MELEMS 663 : LISPRI(1)=NOMD//' ' 664 : LISDUA(1)=NOMD//' ' 665 : 666 : LL=0 667 : NUTOEL=0 668 : DO 101 L=1,NBSOU1 669 : IPT1=MELEME 670 : IF(NBSOU1.NE.1)IPT1=LISOUS(L) 671 : SEGACT IPT1 672 : NOM0=NOMS(IPT1.ITYPEL)//' ' 672 : CALL KALPBG(NOM0,'FONFORM ',IZFFM) 673 : SEGACT IZFFM*MOD 674 : IZHR=KZHR(1) 675 : SEGACT IZHR*MOD 676 : NES=GR(/1) 677 : NPG=GR(/3) 678 : 679 : NP = IPT1.NUM(/1) 680 : MP = NP 681 : NBELEM=IPT1.NUM(/2) 682 : 683 : NNN=MOD(NBELEM,LRV) 684 : IF(NNN.EQ.0) NPACK=NBELEM/LRV 685 : IF(NNN.NE.0) NPACK=1+(NBELEM-NNN)/LRV 686 : 687 : DO 701 KPACK=1,NPACK 688 : 689 : C --- POUR CHAQUE PAQUET DE LRV ELEMENTS ou moins 690 : KDEB=1+(KPACK-1)*LRV 691 : KFIN=MIN(NBELEM,KDEB+LRV-1) 692 : 693 : NBEL=KFIN-KDEB+1 694 : LL=LL+1 695 : 696 : SEGINI IPM1,IPS1 697 : C write(6,*)' IPM1,LL=',IPM1,LL 698 : LIZAFM(LL,1)=IPM1 699 : 700 : IPM2=IPM1 701 : IPM3=IPM1 702 : IPS2=IPS1 703 : IPS3=IPS1 704 : 705 : KITT=2 706 : KJTT=IK1 707 : VISCO=IZTGG1 707 : CALL INITD(VISC,3,1.D0) 708 : 708 : CALL XLAPL(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI, 709 : & VISC,VISC,VISC,KITT,KJTT,IK1, 710 : & IPT1.NUM(1,KDEB),NBEL,NUTOEL,XCOOR,AIMPL,IKOMP, 711 : & IPM1.AM,IPM2.AM,IPM3.AM, 712 : & IPS1.AM,IPS2.AM,IPS3.AM, 713 : & NINKO,IHV,IARG,VISC) 714 : 715 : NUTOEL=NUTOEL+NBEL 716 : SEGDES IPM1 717 : 718 : 701 CONTINUE 719 : SEGDES IPT1 720 : 101 CONTINUE 721 : 722 : 723 : 724 : 725 : 726 : 727 : c%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 728 : 729 : ENDIF 730 : 730 : CALL ECROBJ('MATRIK',MATRIK) 731 : 732 : NAT=2 733 : NSOUPO=0 734 : SEGINI MCHPOI 735 : JATTRI(1)=2 735 : CALL ECROBJ('CHPOINT',MCHPOI) 736 : 737 : C write(6,*)' Fin operateur KMIC' 738 : SEGDES IMATRI,MATRIK 739 : RETURN 740 : 1001 FORMAT(20(1X,I5)) 741 : 1002 FORMAT(10(1X,1PE11.4)) 742 : END 743 : 744 : 745 : 746 : 747 : 748 : 749 : 750 : 751 : 752 : 753 : 754 :
© Cast3M 2003 - All rights reserved.
Disclaimer