1 : C PJMODE SOURCE FANDEUR 13/01/29 21:16:09 7683 2 : SUBROUTINE PJMODE(ipmode) 3 : IMPLICIT INTEGER(I-N) 4 : IMPLICIT REAL*8 (A-H,O-Z) 5 : C======================================================================= 6 : C OPERATEUR PJBA : 7 : C PROJECTION D'UN CHPOINT, D'UN CHARGEMENT OU D'UNE RIGIDITE 8 : C SUR LES ELEMENTS D'UNE BASE MODALE B. 9 : C LE RESULTAT EST DU MEME TYPE. 10 : C 11 : C SYNTAXE : 12 : C * FN = PJBA B OBJET ; SI BASE ELEMENTAIRE 13 : C * FN = PJBA B STR1 (N) OBJET ; SI BASE COMPLEXE 14 : C 15 : C OBJET POUVANT ETRE UNE FORCE OU UN CHARGEMENT, 16 : C OU UNE RIGIDITE DANS LE PREMIER CAS. 17 : C======================================================================= 18 : *********************************************************** 19 : * PROJECTION D'UNE MATRICE SUR UNE BASE DE MODES * 20 : * _______________________________________________________ * 21 : * * 22 : * DATE : le 11 Avril 1995 * 23 : * AUTEUR : Nicolas BENECH * 24 : * _______________________________________________________ * 25 : * * 26 : * MODULE(S) APPELANT(S) : PJBA * 27 : * * 28 : * MODULE(S) APPELE(S) : ACCTAB, YTMX * 29 : * _______________________________________________________ * 30 : * * 31 : * EN ENTREE : * 32 : * MRIGID : Matrice a projeter * 33 : * MTAB1 : Base de modes, reels ou complexes * 34 : * 'REEL' : indique que l'on utilise le produit * 35 : * scalaire reel (pas de conjugaison) * 36 : * * 37 : * EN SORTIE : * 38 : * RI1 : Matrice projetee (partie reelle) * 39 : * RI2 : Matrice projetee (partie imaginaire) * 40 : * _______________________________________________________ * 41 : * * 42 : * REMARQUE : * 43 : * L'operation realisee est : * 44 : * (MTAB1)t . MRIGID. MTAB1 * 45 : * Dans le cas complexe, la transposition est accompagnee * 46 : * de la conjugaison (si REEL n'est pas mentionne). * 47 : * 48 : * voir aussi PROJTA 49 : *********************************************************** 50 : * 51 : -INC SMCHPOI 52 : -INC SMCHARG 53 : -INC SMLCHPO 54 : -INC CCOPTIO 55 : -INC CCGEOME 56 : -INC CCREEL 57 : -INC SMELEME 58 : -INC CCHAMP 59 : -INC SMCHAML 60 : -INC SMMODEL 61 : -INC SMRIGID 62 : -INC SMLMOTS 63 : -INC SMLENTI 64 : 65 : C 66 : * Declarations 67 : * 68 : REAL*8 XVAL, RMAX 69 : CHARACTER*8 LETYPE 70 : CHARACTER*8 TYPMOD,cmate 71 : LOGICAL MODCOM,dedans,dchpo,l3,lr2,lirl 72 : INTEGER I, J, NBMOD, POS, IREEL, IVALRE, IOBRE 73 : REAL*8 XVALRE 74 : LOGICAL LOGRE 75 : segment plcf 76 : integer lpref(ldepl),ldefo(ldepl),lmade(ldepl) 77 : real*8 prmas(ldepl) 78 : endsegment 79 : segment pmod 80 : integer kdefo(nbmod) 81 : endsegment 82 : segment prigmat 83 : integer lrigmat(nrigmat,2+9) 84 : endsegment 85 : segment pmapmo 86 : integer lmapmo(nmapmo),defpmo(nmapmo),dimpmo(nmapmo) 87 : character*4 compmo(nmapmo) 88 : real*8 coepmo(nmapmo) 89 : endsegment 90 : segment pcompo 91 : character*4 mcol 92 : real*8 valmod(nipmod) 93 : endsegment 94 : LOGICAL L0,L1,lcf 95 : PARAMETER (ncod=8) 96 : CHARACTER*4 IDDL,lcod(ncod),lcof(ncod),motinc 97 : CHARACTER*8 TYPRET,CHARRE 98 : data xlopre/1.d-11/ 99 : DATA KZERO/0/ 100 : data lcod/'UX','UY','UZ','RX','RY','RZ','UR','UT'/ 101 : data lcof/'FX','FY','FZ','MX','MY','MZ','FR','FT'/ 102 : 103 : plcf = 0 104 : jgn = 4 105 : jgm = ncod 106 : segini mlmot5 107 : segini mlmot6 108 : do io = 1,ncod 109 : mlmot5.mots(io) = lcod(io) 110 : mlmot6.mots(io) = lcof(io) 111 : enddo 112 : 113 : modcom = .false. 114 : dchpo = .false. 115 : iriout = 0 116 : iriout1 = 0 117 : iriout2 = 0 118 : mmodel = ipmode 119 : segact mmodel 120 : n1 = kmodel(/1) 121 : segini mmode1 122 : jn = 0 123 : do im = 1, n1 124 : imodel = kmodel(im) 125 : segact imodel 126 : if (formod(1).eq.'MECANIQUE'.and.MATMOD(1).eq.'ELASTIQUE' 127 : &.and.(MATMOD(2).eq.'MODAL'.or.MATMOD(2).eq.'STATIQUE')) then 128 : jn = jn + 1 129 : mmode1.kmodel(jn) = imodel 130 : endif 131 : segdes imodel 132 : enddo 133 : segdes mmodel 134 : if (jn.ne.0) then 135 : n1 = jn 136 : segadj mmode1 137 : ipmode = mmode1 138 : else 139 : segsup mmode1 140 : * cas de projection non pr�vue 140 : call erreur(5) 141 : return 142 : endif 143 : 143 : call lirobj('MCHAML',IPCARA,1,iretou) 144 : if (ierr.ne.0) return 145 : 146 : ipchpo = 0 147 : iprigi = 0 147 : call lirobj('CHARGEME',IPCHAR,0,iretou) 147 : if (iretou.eq.0) call lirobj('CHPOINT ',IPCHPO,0,iretou) 147 : if (iretou.eq.0) call lirobj('RIGIDITE',IPRIGI,0,iretou) 148 : 149 : if (iretou.eq.0) then 150 : * manque un op�rande 150 : call erreur(5) 151 : return 152 : endif 153 : 153 : call ecrobj('MCHAML',ipcara) 153 : call ecrobj('MMODEL',ipmode) 153 : call redu 153 : call lirobj('MCHAML',IPCARA,1,iretou) 154 : if (ierr.ne.0) return 155 : 156 : lcf = .false. 157 : mmodel = ipmode 158 : segact mmodel 159 : mchelm = ipcara 160 : segact mchelm 161 : if (ipchar.ne.0) goto 100 162 : if (iprigi.ne.0) goto 200 163 : if (ipchpo.ne.0) then 164 : n = 1 165 : segini mcharg 166 : ipchar = mcharg 167 : segini icharg 168 : kcharg(1) = icharg 169 : ichpo1 = ipchpo 170 : goto 100 171 : endif 172 : 173 : 174 : 100 continue 175 : MCHAR1=IPCHAR 176 : SEGINI,MCHARG=MCHAR1 177 : NBCHG=KCHARG(/1) 178 : DO 10 INCHA=1,NBCHG 179 : ICHAR1=KCHARG(INCHA) 180 : SEGINI,ICHARG=ICHAR1 181 : KCHARG(INCHA)=ICHARG 182 : IP1=ICHPO1 183 : c 184 : IRET = 0 185 : c 186 : c deplacement impose => idepi=1 187 : c force imposee => idepi=0 188 : c 189 : IDEPI = 0 190 : c idepi = -1 191 : KDEPI = 0 192 : MCHPOI = IP1 193 : SEGACT MCHPOI 194 : IF (MTYPOI.EQ.'FLX ') IDEPI = 1 195 : SEGDES MCHPOI 196 : c if (idepi.lt.0) then 197 : c moterr(1:8) = 'chpoint' 197 : c call erreur(302) 198 : c return 199 : c endif 200 : c 201 : NBNN = 1 202 : NBREF = 0 203 : NBSOUS = 0 204 : * 205 : LDEPL = kmodel(/1) 206 : if (.not.lcf) segini plcf 207 : c 208 : c 209 : c **** on initialise le chpoint 210 : c 211 : NSOUPO = 1 212 : NAT=1 213 : SEGINI,MCHPOI 214 : IRET = MCHPOI 215 : MTYPOI = ' ' 216 : MOCHDE=' J''AI ETE FABRIQUE PAR L''OPERATEUR PJBA' 217 : IFOPOI = IFOUR 218 : * champ de force nodal: nature discrete 219 : JATTRI(1)=2 220 : NC = 1 221 : SEGINI,MSOUPO 222 : IPCHP(1) = MSOUPO 223 : NOHARM(1) = NIFOUR 224 : NOCOMP(1) = 'FALF' 225 : 226 : do 101 inocomp=1,2 227 : 228 : N = LDEPL 229 : SEGINI MPOVAL 230 : IPOVAL = MPOVAL 231 : * 232 : NBNN = 1 233 : NBELEM = LDEPL 234 : NBSOUS = 0 235 : NBREF = 0 236 : SEGINI MELEME 237 : IGEOC = MELEME 238 : ITYPEL = 1 239 : 240 : knum = 0 241 : c 242 : c ****boucle sur les chpoints de depl 243 : c 244 : DO 11 IM = 1,kmodel(/1) 245 : 246 : imodel = kmodel(im) 247 : segact imodel 248 : nomid = lnomid(2) 249 : segact nomid 250 : if (.not.lcf) then 251 : ipt1 = imamod 252 : segact ipt1 253 : iptr = ipt1.num(1,1) 254 : lpref(im) = iptr 255 : 256 : indc = 1 257 : 34 if (imache(indc).ne.imamod.or.conche(indc).ne.conmod) then 258 : indc = indc + 1 259 : if (indc.gt.imache(/1)) then 260 : * champ de caracteristiques incomplet 261 : goto 99 262 : endif 263 : goto 34 264 : endif 265 : 266 : mchaml = ichaml(indc) 267 : segact mchaml 268 : do iij = 1, nomche(/2) 269 : if (nomche(iij).eq.'DEFO') then 270 : melval = ielval(iij) 271 : segact melval 272 : ipp1 = ielche(1,1) 273 : ldefo(im) = ipp1 274 : segdes melval 275 : endif 276 : if (nomche(iij).eq.'MADE') then 277 : melval = ielval(iij) 278 : segact melval 279 : ipp2 = ielche(1,1) 280 : lmade(im) = ipp2 281 : segdes melval 282 : endif 283 : if (nomche(iij).eq.'MASS') then 284 : melval = ielval(iij) 285 : segact melval 286 : ymass = velche(1,1) 287 : prmas(im) = ymass 288 : segdes melval 289 : endif 290 : if(ldefo(im).gt.0.and.lmade(im).gt.0.and. 291 : &prmas(im).gt.0) goto 35 292 : enddo 293 : 35 continue 294 : segdes mchaml 295 : if (ldefo(im).eq.0) goto 99 296 : if (prmas(im).le.0.and.cmatee(1:5).eq.'MODAL') goto 99 297 : if (lmade(im).eq.0.and.cmatee(1:8).eq.'STATIQUE') goto 99 298 : 299 : endif 300 : 301 : if (NOCOMP(1).ne.lesobl(1)(1:4)) goto 11 302 : knum = knum + 1 303 : 304 : iptr = lpref(im) 305 : ipp1 = ldefo(im) 306 : NUM(1,knum) = IPTR 307 : ICOLOR(knum) = IDCOUL 308 : XRET = 0.D0 308 : call xty1(ipp1,ip1,mlmot5,mlmot6,xret) 309 : IF (IDEPI.NE.1) THEN 310 : ELSE 311 : * ?? 312 : indn = 1 313 : 45 if (nomche(indn).ne.'FREQ') then 314 : indn = indn + 1 315 : if (indn.gt.nomche(/2)) then 316 : * pas la composante FREQ 317 : goto 99 318 : endif 319 : goto 45 320 : endif 321 : 322 : melval = ielval(indn) 323 : segact melval 324 : x1 = velche(1,1) 325 : OM = X1 326 : OM = 2.D0 * XPI * OM 327 : OM = OM * OM 328 : XRET = -XRET / OM 329 : ENDIF 330 : IF (IFOUR .EQ. 1) THEN 331 : IF (NIFOUR .NE. 0) THEN 332 : XRET = XRET*XPI 333 : ELSE 334 : XRET = XRET*2.D0*XPI 335 : ENDIF 336 : ENDIF 337 : VPOCHA(knum,1) = XRET 338 : * 339 : if (cmatee(1:5).eq.'MODAL') then 340 : ymass = prmas(im) 341 : elseif (cmatee(1:8).eq.'STATIQUE') then 342 : ipp2 = lmade(im) 342 : call xty1(ipp1,ipp2,mlmot5,mlmot6,ymass) 343 : else 344 : endif 345 : if (lmade(im).gt.0.and.dabs(XRET).gt.(1.d-10*ymass).and. 346 : & ymass.gt.0.and.cmatee(1:5).eq.'MODAL') then 347 : * kich : on enleve la projection sur base modale - a creuser pour statique ! 347 : CALL ADCHPO(IP1,IPP2,IP2,1.d0,(XRET/ymass*(-1.))) 348 : IP1 = IP2 349 : endif 350 : * 351 : 11 CONTINUE 352 : * 353 : lcf = .true. 354 : * 355 : * 356 : if (knum.eq.kmodel(/1)) goto 102 357 : if (inocomp.eq.1) then 358 : if (knum.eq.0) then 359 : NOCOMP(1) = 'FBET' 360 : else 361 : N = knum 362 : NBELEM = knum 363 : segadj MPOVAL,MELEME 364 : SEGDES MPOVAL,MELEME 365 : NSOUPO = 2 366 : segadj MCHPOI 367 : SEGINI,MSOUPO 368 : IPCHP(2) = MSOUPO 369 : NOCOMP(1) = 'FBET' 370 : endif 371 : endif 372 : 101 continue 373 : 374 : 102 continue 375 : N = knum 376 : NBELEM = knum 377 : segadj MPOVAL,MELEME 378 : SEGDES MPOVAL,MELEME 379 : SEGDES MSOUPO 380 : 381 : SEGDES MCHPOI 382 : IF(IERR.NE.0) RETURN 383 : ICHPO1=IRET 384 : SEGDES,ICHARG 385 : 10 CONTINUE 386 : segsup mlmot5,mlmot6,plcf 387 : if (ipchpo.gt.0) then 388 : segsup icharg,mcharg 388 : call ecrobj('CHPOINT ',iret) 389 : goto 999 390 : endif 391 : SEGDES,MCHARG 391 : CALL ECROBJ('CHARGEME',MCHARG) 392 : 393 : goto 999 394 : 99 segsup mpoval,msoupo,mchpoi 394 : call erreur(26) 395 : return 396 : 397 : 398 : 200 continue 399 : ipri1 = iprigi 399 : call SEPA(ipri1,1) 400 : ipri2 = iprigi 400 : call SEPA(ipri2,2) 401 : * 402 : * 403 : * 404 : * 405 : nmapmo = 100 406 : kpmo = 0 407 : segini pmapmo 408 : do isous = 1,kmodel(/1) 409 : imodel = kmodel(isous) 410 : segact imodel 411 : cmate = cmatee 412 : meleme = imamod 413 : segact meleme 413 : if (itypel.ne.1) call erreur(5) 413 : if (num(/1).ne.1) call erreur(5) 414 : if (cmate.eq.'STATIQUE'.or.cmate.EQ.'MODAL') then 415 : do ilp = 1,num(/2) 416 : kpmo = kpmo + 1 417 : if (kpmo.gt.nmapmo) then 418 : nmapmo = nmapmo + 100 419 : segadj pmapmo 420 : endif 421 : lmapmo(kpmo) = num(1,ilp) 422 : if (cmate.eq.'STATIQUE') then 423 : compmo(kpmo) = 'BETA' 424 : elseif (cmate.eq.'MODAL') then 425 : compmo(kpmo) = 'ALFA' 426 : endif 427 : do im = 1 , imache(/1) 428 : if (imache(im).eq.imamod) then 429 : if (conche(im).eq.conmod) then 430 : mchaml = ichaml(im) 431 : segact mchaml 432 : do iv = 1,ielval(/1) 433 : if (nomche(iv).eq.'DEFO') then 434 : melval = ielval(iv) 435 : segact melval 436 : ibmn = min(ilp,ielche(/2)) 437 : defpmo(kpmo) = ielche(1,ibmn) 438 : endif 439 : if (nomche(iv).eq.'IDEF') then 440 : melval = ielval(iv) 441 : segact melval 442 : ibmn = min(ilp,ielche(/2)) 443 : dimpmo(kpmo) = ielche(1,ibmn) 444 : endif 445 : enddo 446 : endif 447 : endif 448 : enddo 449 : 450 : enddo 451 : endif 452 : segdes meleme,imodel 453 : enddo 454 : 455 : nmapmo = kpmo 456 : segadj pmapmo 457 : nbmod = nmapmo 458 : * 459 : N1 = NBMOD 460 : nbcod = 8 461 : segini pmod 462 : SEGINI, MLCHP1 463 : SEGINI, MLCHP2 464 : jgm = 1 465 : jgn = 4 466 : segini mlmot4 467 : * 468 : * Constitution du maillage support et du segment descriptif 469 : * 470 : NBNN = NBMOD 471 : NBELEM = 1 472 : NBSOUS = 0 473 : NBREF = 0 474 : SEGINI, MELEME 475 : ITYPEL = 1 476 : * 477 : NLIGRD=NBMOD 478 : NLIGRP=NBMOD 479 : SEGINI, DESCR 480 : * 481 : mrigid = ipri1 482 : segact mrigid 483 : nrigel = coerig(/1) 484 : if (nrigel.lt.1) goto 250 485 : typmod = ' ' 486 : IREEL = -1 487 : C* POS ? IF (POS.EQ.1) IREEL = 1 488 : * 489 : LETYPE = ' ' 490 : DO 210 IM=1,NBMOD 491 : IPT1 = 0 492 : * 493 : imodel = kmodel(im) 494 : segact imodel 495 : ipt1 = imamod 496 : segact ipt1 497 : iptr = ipt1.num(1,1) 498 : * Cas reel ou cas complexe ? 499 : * 500 : if (dimpmo(im).gt.0) TYPMOD = 'MODE_COM' 501 : 502 : IF (TYPMOD .EQ. 'MODE_COM') THEN 503 : MODCOM=.TRUE. 504 : mchpoi = defpmo(im) 505 : MLCHP1.ICHPOI(IM) = MCHPOI 506 : mchpoi = dimpmo(im) 507 : MLCHP2.ICHPOI(IM) = MCHPOI 508 : ELSE 509 : MODCOM = .FALSE. 510 : mchpoi = defpmo(im) 511 : MLCHP1.ICHPOI(IM) = MCHPOI 512 : ENDIF 513 : * 514 : ipt1 = iptr 515 : * 516 : MELEME.NUM(IM,1)=IPT1 517 : * 518 : if (cmatee.eq.'MODAL') then 519 : DESCR.LISINC(IM) = 'ALFA' 520 : DESCR.LISDUA(IM) = 'FALF' 521 : else if (cmatee.eq.'STATIQUE') then 522 : DESCR.LISINC(IM) = 'BETA' 523 : DESCR.LISDUA(IM) = 'FBET' 524 : endif 525 : DESCR.NOELEP(IM) = IM 526 : DESCR.NOELED(IM) = IM 527 : * 528 : 210 CONTINUE 529 : SEGDES, DESCR 530 : SEGDES, MELEME 531 : * 532 : * Constitution des segments XMATRI 533 : * 534 : NLIGRD=NBMOD 535 : NLIGRP=NBMOD 536 : nelrig=1 537 : * 538 : IF (LETYPE .EQ. 'ANNULE') THEN 539 : SEGINI, XMATR1 540 : IF (MODCOM) THEN 541 : SEGINI, XMATR2 542 : SEGDES, XMATR2 543 : ENDIF 544 : SEGDES, XMATR1 545 : GOTO 55 546 : ENDIF 547 : * 548 : * Cas reel 549 : * 550 : SEGINI, XMATR1 551 : DO 20, I=1, NBMOD 552 : MCHPO1 = MLCHP1.ICHPOI(I) 553 : DO 20, J=1, NBMOD 554 : MCHPO2 = MLCHP1.ICHPOI(J) 554 : CALL YTMX (MCHPO2, MCHPO1, MRIGID, XVAL) 555 : XMATR1.RE(I,J,1)=XVAL 556 : 20 CONTINUE 557 : SEGDES, XMATR1 558 : * 559 : * Cas complexe : calcul de termes complementaires 560 : * 561 : IF (MODCOM) THEN 562 : SEGACT, XMATR1*mod 563 : DO 30, I=1, NBMOD 564 : MCHPO1 = MLCHP2.ICHPOI(I) 565 : DO 30, J=1, NBMOD 566 : MCHPO2 = MLCHP2.ICHPOI(J) 566 : CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL) 567 : XMATR1.RE(I,J,1)=XMATR1.RE(I,J,1)-IREEL*XVAL 568 : 30 CONTINUE 569 : SEGDES, XMATR1 570 : * 571 : SEGINI, XMATR2 572 : DO 40, I=1, NBMOD 573 : MCHPO1 = MLCHP1.ICHPOI(I) 574 : DO 40, J=1, NBMOD 575 : MCHPO2 = MLCHP2.ICHPOI(J) 575 : CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL) 576 : XMATR2.RE(I,J,1)=XVAL 577 : 40 CONTINUE 578 : DO 50, I=1, NBMOD 579 : MCHPO1 = MLCHP2.ICHPOI(I) 580 : DO 50, J=1, NBMOD 581 : MCHPO2 = MLCHP1.ICHPOI(J) 581 : CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL) 582 : XMATR2.RE(I,J,1)=XMATR2.RE(I,J,1)+IREEL*XVAL 583 : 50 CONTINUE 584 : SEGDES, XMATR2 585 : ENDIF 586 : SEGDES, MLCHP1 587 : SEGDES, MLCHP2 588 : * 589 : SEGACT, MRIGID 590 : LETYPE = MRIGID.MTYMAT 591 : SEGDES, MRIGID 592 : * 593 : * Creation des segments IMATRI 594 : * 595 : 55 NELRIG = 1 596 : * SEGINI, IMATR1 597 : * IMATR1.IMATTT(1) = XMATR1 598 : SEGDES, xMATR1 599 : IF (MODCOM) THEN 600 : * SEGINI, IMATR2 601 : * IMATR2.IMATTT(1) = XMATR2 602 : SEGDES, xMATR2 603 : ENDIF 604 : * 605 : * Creation des rigidites calculees 606 : * 607 : NRIGE=7 608 : NRIGEL=1 609 : SEGINI, RI1 610 : RI1.MTYMAT = LETYPE 611 : RI1.IFORIG = IFOMOD 612 : RI1.IMGEO1 = 0 613 : RI1.IMGEO2 = 0 614 : RI1.COERIG = 1.D0 615 : RI1.IRIGEL(1,1) = MELEME 616 : RI1.IRIGEL(2,1) = 0 617 : RI1.IRIGEL(3,1) = DESCR 618 : RI1.IRIGEL(4,1) = xMATR1 619 : RI1.IRIGEL(5,1) = NIFOUR 620 : RI1.IRIGEL(6,1) = 0 621 : RI1.IRIGEL(7,1) = 2 622 : SEGDES, RI1 623 : IF (MODCOM) THEN 624 : SEGINI, RI2 = RI1 625 : RI2.IRIGEL(4,1) = xMATR2 626 : SEGDES, RI2 627 : ELSE 628 : RI2 = 0 629 : SEGSUP, MLCHP2 630 : ENDIF 631 : * 632 : iriout1 = ri1 633 : iriout2 = ri2 634 : 635 : 250 continue 636 : mrigid = ipri2 637 : segact mrigid 638 : nrigel = coerig(/1) 639 : if (nrigel.lt.1) goto 290 640 : typmod = ' ' 641 : 642 : nrigmat =100 643 : kgmat = 0 644 : segini prigmat 645 : 646 : KRIGEL = 0 647 : nrigel = irigel(/2) 648 : nrige = irigel(/1) 649 : segini ri1 650 : ri1.mtymat = mtymat 651 : nrige0 = nrigel 652 : 653 : kige = 0 654 : kige1 = 100 655 : nrigel = kige1 656 : segini ri2 657 : ri2.mtymat = mtymat 658 : 659 : DO ire = 1,nrige0 660 : meleme = irigel (1,ire) 661 : segact meleme 662 : if (itypel.ne.22) then 662 : call erreur(977) 663 : return 664 : endif 665 : nbelem = num(/2) 666 : nbele0 = nbelem 667 : descr = irigel(3,ire) 668 : segact descr 669 : nligrp0 = noelep(/1) 670 : nligrd0 = noeled(/1) 671 : nligrp = nligrp0 + nmapmo 672 : nligrd = nligrd0 + nmapmo 673 : 674 : nbnn = num(/1) 675 : nbsous = 0 676 : nbref = 0 677 : segini ipt2 678 : ipt2.itypel = itypel 679 : nbelem = 1 680 : nbnn = num(/1) + nmapmo 681 : segini ipt1 682 : ipt1.itypel = itypel 683 : ri1.coerig(ire) = coerig(ire) 684 : kele = 0 685 : 686 : xmatr1 = irigel(4,ire) 687 : segact xmatr1 688 : nelrig0 = xmatr1.re(/3) 689 : nelrig = nelrig0 + nmapmo 690 : segini xmatr2 691 : DO iele = 1,nbele0 692 : ie2 = min(iele,nelrig0) 693 : * xmatr1 = imatr1.imattt(ie2) 694 : * segact xmatr1 695 : nligrp = nligrp0 + nmapmo 696 : nligrd = nligrd0 + nmapmo 697 : nelrig=1 698 : segini des2,xmatri 699 : des2.lisinc(1) = 'LX' 700 : des2.lisdua(1) = 'FLX' 701 : des2.noelep(1) = 1 702 : des2.noeled(1) = 1 703 : * le premier point correspond aux multiplicateurs 703 : CALL CREPO1 (ZERO, ZERO, ZERO, IPTS) 704 : ipt1.num(1,1) = ipts 705 : kgrp = 1 706 : kirp = 1 707 : do ipmo = 1,nmapmo 708 : coepmo(ipmo) = 0. 709 : enddo 710 : do igrp = 2,nligrp0 711 : jno = noelep(igrp) 712 : motinc = lisinc(igrp) 713 : IP1 = num(jno,iele) 714 : * recherche association noeud physique - points support déformée 715 : do ilmat = 1,kgmat 716 : if(lrigmat(ilmat,1).eq.ip1) goto 315 717 : enddo 718 : 719 : kgmat = kgmat+1 720 : ilmat = kgmat 721 : if (kgmat.gt.nrigmat) then 722 : nrigmat = nrigmat + 100 723 : segadj prigmat 724 : endif 725 : kpb = 0 726 : jg = 100 727 : segini mlent3 728 : lrigmat(kgmat,1) = ip1 729 : do ikmo = 1, nmapmo 730 : ichp1 = defpmo(ikmo) 730 : call ecrcha('NOMU') 730 : call ecrcha('MAIL') 730 : call ecrobj('CHPOINT ',ichp1) 730 : call extrai 730 : call lirobj('MAILLAGE',iuu1,1,iretou) 730 : call ecrobj('MAILLAGE',iuu1) 730 : call ecrobj('POINT ',IP1) 730 : call DANS 730 : call lirlog(l3,1,iretou) 731 : if(l3) then 732 : kpb = kpb + 1 733 : if (kpb.gt.jg) then 734 : jg = jg + 100 735 : segadj mlent3 736 : endif 737 : mlent3.lect(kpb) = ikmo 738 : endif 739 : enddo 740 : jg = kpb 741 : segadj mlent3 742 : if (kpb.gt.0) then 743 : lrigmat(ilmat,2) = mlent3 744 : else 745 : lrigmat(ilmat,2) = 0 746 : segsup mlent3 747 : endif 748 : 749 : 315 continue 750 : ilr3 = lrigmat(ilmat,2) 751 : if (ilr3.eq.0) goto 253 752 : mlent3 = ilr3 753 : segact mlent3 754 : * selection selon nom composante 755 : mlmat = 0 756 : do lmo = 1,9 757 : if (motinc.eq.lcod(lmo)) mlmat = lmo+2 758 : enddo 759 : if (mlmat.eq.0) then 760 : * write(6,*) 'coefs pour cette composante non trouves' 761 : goto 253 762 : endif 763 : if (lrigmat(ilmat,mlmat).ne.0) then 764 : pcompo = lrigmat(ilmat,mlmat) 765 : segact pcompo 766 : nipmod = valmod(/1) 767 : do ilg = 1,nipmod 768 : lkmo = mlent3.lect(ilg) 769 : coepmo(lkmo) = (valmod(ilg)* xmatr1.re(1,igrp,ie2)) 770 : & + coepmo(lkmo) 771 : enddo 772 : else 773 : jg = mlent3.lect(/1) 774 : nipmod = jg 775 : segini pcompo 776 : mcol = motinc 777 : do ilg = 1,nipmod 778 : lkmo = mlent3.lect(ilg) 779 : ichp1 = defpmo(lkmo) 779 : CALL EXTRA9(ICHP1,ip1,motinc,KERRE,XFLOT) 780 : coepmo(lkmo) = (xflot * xmatr1.re(1,igrp,ie2)) 781 : & + coepmo(lkmo) 782 : valmod(ilg) = xflot 783 : enddo 784 : lrigmat(ilmat,mlmat) = pcompo 785 : endif 786 : 787 : 253 continue 788 : enddo 789 : 790 : xmaut1 = 0.d0 791 : do kpmo = 1,nmapmo 792 : xmaut1 = max(xmaut1,dabs(coepmo(kpmo))) 793 : enddo 794 : 795 : * synthèse 796 : do igrp = 2,nligrp0 797 : jno = noelep(igrp) 798 : motinc = lisinc(igrp) 799 : IP1 = num(jno,iele) 800 : lr2 = .false. 801 : do jgmat = 1,kgmat 802 : if(lrigmat(jgmat,1).eq.ip1) goto 325 803 : enddo 804 : c write(6,*) 'bizarre, point dans l element non repertorie' 804 : call erreur(5) 805 : return 806 : 325 continue 807 : mlmat = 0 808 : do lmo = 1,9 809 : if (motinc.eq.lcod(lmo)) mlmat = lmo+2 810 : enddo 811 : if (mlmat.eq.0) lr2 = .true. 812 : if (lrigmat(jgmat,mlmat).eq.0) lr2 = .true. 813 : if(lr2) then 814 : jirp = 0 815 : do iirp = 1,kgrp 816 : if (ipt1.num(iirp,1).eq.ip1) jirp = iirp 817 : enddo 818 : c recopie 819 : kgrp = kgrp + 1 820 : if (jirp.ne.0) then 821 : des2.noelep(kgrp) = des2.noelep(jirp) 822 : des2.noeled(kgrp) = des2.noeled(jirp) 823 : else 824 : kirp = kirp + 1 825 : ipt1.num(kirp,1) = ip1 826 : des2.noelep(kgrp) = kirp 827 : des2.noeled(kgrp) = kirp 828 : endif 829 : des2.lisinc(kgrp) = lisinc(igrp) 830 : des2.lisdua(kgrp) = lisdua(igrp) 831 : re(1,kgrp,1) = xmatr1.re(1,igrp,ie2) 832 : re(kgrp,1,1) = re(1,kgrp,1) 833 : endif 834 : * 835 : enddo 836 : 837 : do kpmo = 1,nmapmo 838 : if (dabs(coepmo(kpmo)).gt.xlopre*xmaut1) then 839 : kirp = kirp + 1 840 : kgrp = kgrp + 1 841 : ipt1.num(kirp,1) = lmapmo(kpmo) 842 : des2.noelep(kgrp) = kirp 843 : des2.noeled(kgrp) = kirp 844 : des2.lisinc(kgrp) = compmo(kpmo) 845 : if (compmo(kpmo).eq.'ALFA') des2.lisdua(kgrp) = 'FALF' 846 : if (compmo(kpmo).eq.'BETA') des2.lisdua(kgrp) = 'FBET' 847 : re(1,kgrp,1) = coepmo(kpmo) 848 : re(kgrp,1,1) = re(1,kgrp,1) 849 : endif 850 : enddo 851 : * 852 : lirl = .false. 853 : if (kirp.ne.num(/1)) then 854 : lirl = .true. 855 : else 856 : do io = 1,kirp 857 : if (num(io,iele).ne.ipt1.num(io,1)) lirl=.true. 858 : enddo 859 : endif 860 : c creation d'un irigel 861 : if (lirl) then 862 : kige = kige + 1 863 : if (kige.gt.kige1) then 864 : nrigel = kige1 + 100 865 : segadj ri2 866 : kige1 = nrigel 867 : endif 868 : nbelem = 1 869 : nbnn = kirp 870 : segini ipt3 871 : ipt3.itypel = itypel 872 : do io =1,nbnn 873 : ipt3.num(io,1) = ipt1.num(io,1) 874 : enddo 875 : nligrp = kgrp 876 : nligrd = kgrp 877 : nelrig=1 878 : segadj xmatri,des2 879 : * segini imatr3 880 : * imatr3.imattt(1) = xmatri 881 : segdes ipt3,des2,xmatri 882 : RI2.IRIGEL(1,kige) = IPT3 883 : RI2.IRIGEL(3,kige) = DES2 884 : RI2.IRIGEL(4,kige) = xmatri 885 : RI2.IRIGEL(2,kige) = 0 886 : RI2.IRIGEL(5,kige) = irigel(5,ire) 887 : RI2.IRIGEL(6,kige) = irigel(6,ire) 888 : ri2.coerig(kige) = coerig(ire) 889 : else 890 : * relation non modifiee pour cet element 891 : kele = kele + 1 892 : do ig = 1,nligrp0 893 : ipt2.num(ig,kele) = ipt1.num(ig,1) 894 : enddo 895 : * imatr2.imattt(kele) = xmatr1 896 : * kich : a tester 897 : do ju = 1,kgrp 898 : xmatr2.re(1,ju,kele) = re(1,ju,1) 899 : xmatr2.re(ju,1,kele) = re(ju,1,1) 900 : enddo 901 : segsup xmatri,des2 902 : endif 903 : ENDDO 904 : 905 : nbelem = kele 906 : nelrig = kele 907 : nligrd=xmatr2.re(/1) 908 : nligrp=xmatr2.re(/2) 909 : if (nbelem.gt.0) then 910 : segadj ipt2 911 : segadj xmatr2 912 : krigel = krigel + 1 913 : RI1.IRIGEL(1,krigel) = IPT2 914 : RI1.IRIGEL(3,krigel) = irigel(3,ire) 915 : RI1.IRIGEL(4,krigel) = xmatr2 916 : RI1.IRIGEL(2,krigel) = 0 917 : RI1.IRIGEL(5,krigel) = irigel(5,ire) 918 : RI1.IRIGEL(6,krigel) = irigel(6,ire) 919 : segdes ipt2,xmatr2 920 : else 921 : segsup ipt2 922 : endif 923 : segsup ipt1 924 : ENDDO 925 : 926 : iriout = 0 927 : nrigel = krigel 928 : segadj ri1 929 : nrigel = kige 930 : segadj ri2 931 : segdes mrigid,ri1,ri2 932 : if (kige.eq.0) segsup ri2 933 : if (krigel.eq.0) segsup ri1 934 : if (kige.gt.0.and.krigel.gt.0) then 935 : c write(6,*) 'fus', ri1,ri2,kige,krigel 935 : call fusrig(ri1,ri2,iriout) 936 : segsup ri1, ri2 937 : return 938 : endif 939 : if (kige.gt.0) iriout = ri2 940 : if (krigel.gt.0) iriout = ri1 940 : if (iriout.eq.0) call erreur(-5) 941 : c write(6,*) 'iriout', iriout 942 : 943 : 290 continue 944 : if (iriout.ne.0) iriout3 = iriout 945 : if (iriout1.ne.0) iriout3 = iriout1 946 : if (iriout.ne.0.and.iriout1.ne.0) then 946 : call fusrig(iriout, iriout1,iriout3) 947 : ri1 = iriout 948 : ri2 = iriout1 949 : segsup ri1,ri2 950 : endif 951 : 951 : call ecrobj('RIGIDITE',iriout3) 951 : if (modcom) call ecrobj('RIGIDITE',iriout2) 952 : 953 : goto 999 954 : 955 : 199 continue 956 : segsup descr,meleme,mlchp1,mlchp2 956 : call erreur(5) 957 : return 958 : 959 : 999 continue 960 : 961 : do ich = 1, ichaml(/1) 962 : mchaml = ichaml(ich) 963 : segact mchaml 964 : do iel =1,ielval(/1) 965 : melval = ielval(iel) 966 : segdes melval 967 : enddo 968 : segdes mchaml 969 : enddo 970 : do k1 = 1,kmodel(/1) 971 : imodel = kmodel(k1) 972 : segdes imodel 973 : enddo 974 : segsup mmode1, mchelm 975 : if (plcf.ne.0) segsup plcf 976 : return 977 : END 978 : 979 :
© Cast3M 2003 - All rights reserved.
Disclaimer