Download excit1.eso

Back to the list

   1 : C EXCIT1    SOURCE    PV        13/09/26    21:16:37     7827           
   2 : 
   3 :       SUBROUTINE EXCIT1(MRIGID,MCHPO1,MCHPO2,MCHPO3,MLENTI,indic,RI2,
   4 :      &                  ipt8,ITYP)
   5 : 
   6 :       IMPLICIT INTEGER(I-N)
   7 :       IMPLICIT REAL*8(A-H,O-Z)
   8 : 
   9 : *  modif : declarer inactifs les blocages redondants s'ils sont dans 
  10 : *  ipt8. Si celui-ci est non vide, on enleve les blocages inclus 
  11 : *  dedans et on s'en tient la.
  12 : 
  13 : -INC CCOPTIO
  14 : 
  15 : -INC SMRIGID
  16 : -INC SMCHPOI
  17 : -INC SMELEME
  18 : -INC SMCOORD
  19 : -INC SMLENTI
  20 : 
  21 :       SEGMENT SNOMIN 
  22 :         CHARACTER*4 NOMIN(0) 
  23 :       ENDSEGMENT 
  24 :       SEGMENT ISOU(IRIGEL(/2)) 
  25 :       SEGMENT ICPR(XCOOR(/1)/idimp1)
  26 :       SEGMENT JCPR(XCOOR(/1)/idimp1)
  27 :       SEGMENT KCPR(XCOOR(/1)/idimp1)
  28 :       SEGMENT YCPR(XCOOR(/1)/idimp1)
  29 :       SEGMENT ITRAV
  30 :         real*8 DVAL(NIN,NPO),DIN(NPO) 
  31 :         integer IEXI(NIN,NPO)
  32 :       endsegment
  33 :       SEGMENT  IPASS(ITEMP) 
  34 :       segment trav1
  35 :         integer iav(nbo)
  36 :         real*8 viol(nbo)
  37 :       endsegment
  38 : C
  39 : C  ***  ICPR REFERENCIE DVAL 
  40 : C  ***  NOMIN CONTIENT LES INCONNUES REFERENCES PAR LES RELATIONS 
  41 : C  ***  DVAL CONTIENT LE RESULTAT DES DEPLACEMENTS ET DES LAMBDA 
  42 : C  ***  DIN CONTIENT LES VALEURS INITIALES DU SECOND MEMBRE 
  43 : C
  44 :       logical indic
  45 : C
  46 :       idimp1 = idim+1
  47 :       segact mcoord
  48 : C 
  49 :       xcrot  = 0.
  50 :       XMCRIT = 0.
  51 :       YMCRIT = 0.
  52 :       ymult  = 1.
  53 : ***   if (ityp.eq.3) ymult=1.D3
  54 :       if (ityp.eq.2) xcrot=0.250
  55 : 
  56 :       icpr=0
  57 :       jcpr=0
  58 :       ycpr=0
  59 :       kcpr=0
  60 :       itrav=0
  61 : 
  62 : * on a un maillage de relation(s) a enlever
  63 :       if (ipt8.ne.0) then
  64 :         segact,ipt8
  65 :         if (ipt8.itypel.ne.1) then
  66 :           segdes,ipt8
  66 :           call erreur(16)
  67 :           return
  68 :         endif
  69 :         if (ipt8.num(/2).ne.0) then
  70 :           segini kcpr
  71 :           do i=1,ipt8.num(/2)
  72 :             kcpr(ipt8.num(1,i))=1
  73 :           enddo
  74 :         endif
  75 :         segdes,ipt8
  76 :       endif
  77 : *
  78 : *  mchpo3 est le champ des LX sorti de excfro jcpr donne existence
  79 : *   et ycpr donne valeur
  80 : *     write (ioimp,*) ' mchpo3 dans excit1 ',mchpo3
  81 :       if (mchpo3.ne.0) then
  82 :         segini jcpr,ycpr
  83 :         segact mchpo3
  84 :         do 100 isoupo=1,mchpo3.ipchp(/1)
  85 :           msoupo=mchpo3.ipchp(isoupo)
  86 :           segact msoupo
  87 :           do 110 ic=1,nocomp(/2)
  88 :             if (nocomp(ic).eq.'LX  ') goto 120
  89 :  110      continue
  90 :           goto 140
  91 :  120      continue
  92 :           ipt2=igeoc
  93 :           mpoval=ipoval
  94 :           segact ipt2,mpoval
  95 :           do 130 iel=1,ipt2.num(/2)
  96 :             impt=ipt2.num(1,iel)
  97 :             jcpr(impt)=1
  98 :             ycpr(impt)=ycpr(impt)+vpocha(iel,ic)
  99 :  130      continue
 100 :           segdes,ipt2,mpoval
 101 :  140      continue
 102 :           segdes,msoupo
 103 :  100    continue
 104 :         segdes mchpo3
 105 :       endif
 106 : C
 107 : C** INITIALISATION DE NOMIN ET ICPR icpr repere les points appartenant 
 108 : C** a la matrice de blocage mini maxi ou ? 
 109 : C
 110 :       SEGINI,SNOMIN,ICPR
 111 :       NOMIN(**)='LX  ' 
 112 :       NPO=0 
 113 :       NBO=0 
 114 : *
 115 :       SEGACT MRIGID 
 116 :       DO 1 I=1,IRIGEL(/2) 
 117 :         ITY=IRIGEL(6,I)
 118 :         IF (ITY.EQ.0) GO TO 1
 119 : *  cas du frottement : petite verification
 120 :         if (ity.eq.2 .and. mchpo3.eq.0) then
 120 :           call erreur(721)
 121 :           return
 122 :         endif
 123 :         MELEME=IRIGEL(1,I) 
 124 :         SEGACT MELEME 
 125 :         NBO=NBO+NUM(/2) 
 126 :         IF(IIMPI.EQ.528) WRITE(IOIMP,3765) NBO 
 127 :  3765   FORMAT('  NBO  ',I5) 
 128 :         DO 2 J=1,NUM(/2) 
 129 :           DO 2 K=1,NUM(/1) 
 130 :             impt=NUM(K,J)
 131 :             IF (ICPR(impt).EQ.0) THEN
 132 :               NPO=NPO+1 
 133 :               ICPR(impt)=NPO
 134 :             ENDIF
 135 :     2   CONTINUE 
 136 :         SEGDES MELEME 
 137 :         DESCR=IRIGEL(3,I) 
 138 :         SEGACT DESCR 
 139 :         DO 3 J=1,LISINC(/2) 
 140 :           DO 4 K=1,NOMIN(/2) 
 141 :             IF (NOMIN(K).EQ.LISINC(J)) GO TO 3 
 142 :     4     CONTINUE 
 143 :           NOMIN(**)=LISINC(J) 
 144 :     3   CONTINUE 
 145 :         SEGDES DESCR 
 146 :  1    CONTINUE 
 147 :       NIN=NOMIN(/2) 
 148 : *
 149 :       IF(IIMPI.EQ.528) THEN
 150 :         WRITE(ioimp,90) NPO
 151 :         WRITE(IOIMP,9876)( NOMIN(i),i=1,NIN) 
 152 :   90    FORMAT(' ON VIENT DE PASSER LA BOUCLE 1  NPO  ',I5)
 153 :  9876   FORMAT('  NOMIN ' ,10(A4,1X)) 
 154 :       ENDIF
 155 : C 
 156 : C ****  ON REMPLIT DVAL AVEC LES VALEURS DE MCHPO2 
 157 : C   mchpo2 est le champ de deplacement propose par RESOU
 158 : C   dval(j,icpr(i)) est la j eme composante du deplacement du point i 
 159 : C   et iexi(j,icpr(i))=1
 160 : C
 161 :       SEGINI ITRAV 
 162 : c
 163 :       MCHPOI=MCHPO2 
 164 :       SEGACT MCHPOI 
 165 :       DO 5 I=1,IPCHP(/1) 
 166 :         MSOUPO=IPCHP(I) 
 167 :         SEGACT MSOUPO 
 168 :         MELEME=IGEOC 
 169 :         SEGACT MELEME 
 170 :         ITEMP=NOCOMP(/2) 
 171 :         SEGINI IPASS 
 172 :         DO 6 K=1,ITEMP
 173 :           DO 7 J=1,NIN 
 174 :             IF (NOMIN(J).EQ.NOCOMP(K)) THEN
 175 :               IPASS(K)=J 
 176 :               GO TO 6 
 177 :             ENDIF
 178 :    7      CONTINUE 
 179 :    6    CONTINUE 
 180 :         IF (IIMPI.EQ.528) WRITE(IOIMP,1555) (IPASS(KHU),KHU=1,ITEMP) 
 181 :  1555 FORMAT('   IPASS   ' ,9I10) 
 182 :         MPOVAL=IPOVAL 
 183 :         SEGACT MPOVAL 
 184 :         III=0 
 185 :         DO 8 J=1,ITEMP
 186 :           K=IPASS(J) 
 187 :           IF (K.EQ.0) GO TO 8 
 188 :           DO 9 L=1,NUM(/2) 
 189 :             IP=ICPR(NUM(1,L)) 
 190 :             IF (IP.EQ.0) GO TO 9 
 191 : *     if (k.eq.1.and.vpocha(l,j).eq.0.d0 ) then
 192 : *      write (6,*) 'LX nul dans excit1'
 193 : *      goto 9
 194 : *     endif
 195 :             DVAL(K,IP)=dval(k,ip)+VPOCHA(L,J)
 196 :             IEXI(K,IP)=1 
 197 :             III=III+1 
 198 :    9      CONTINUE 
 199 :    8    CONTINUE 
 200 :         SEGSUP IPASS 
 201 :         SEGDES MPOVAL,MSOUPO,MELEME 
 202 :    5  CONTINUE 
 203 :       SEGDES MCHPOI 
 204 : *
 205 :       IF (IIMPI.EQ.528) then
 206 :         WRITE(IOIMP,1556)III,(DVAL(1,i),i=1,NPO) 
 207 :         WRITE(IOIMP,101) 
 208 :  1556   FORMAT(' III  DVAL ',I6,/,(1X,10E12.5)) 
 209 :   101   FORMAT(' ON VIENT DE PASSER LA BOUCLE 5') 
 210 :       endif
 211 : C 
 212 : C  **** ON REMPLIT DIN PAR LES VALEURS DE MCHPO1 POUR LES LAMBDAS 
 213 : C    mchpo1 est le vecteur initial FLX ( deplacement initial)
 214 : C    din (icpr(i)) est le deplacement (FLX) initial du point i 
 215 : C
 216 :       MCHPOI=MCHPO1 
 217 :       SEGACT MCHPOI 
 218 :       DO 10 i=1,IPCHP(/1) 
 219 :         MSOUPO=IPCHP(I) 
 220 :         SEGACT MSOUPO 
 221 :         MELEME=IGEOC 
 222 :         SEGACT MELEME 
 223 :         MPOVAL=IPOVAL 
 224 :         SEGACT MPOVAL 
 225 :         DO 11 J=1,NOCOMP(/2) 
 226 :           IF (NOCOMP(J).NE.'FLX ') GO TO 11 
 227 :           DO 12 K=1,NUM(/2) 
 228 :             IP=ICPR(NUM(1,K)) 
 229 :             IF (IP.EQ.0) GO TO 12 
 230 :             DIN(IP)=VPOCHA(K,J) 
 231 :    12     CONTINUE 
 232 :    11   CONTINUE 
 233 :         SEGDES MPOVAL,MELEME,MSOUPO 
 234 :    10 CONTINUE 
 235 :       SEGDES MCHPOI 
 236 : *
 237 :       IF (IIMPI.EQ.528) then
 238 :         WRITE(IOIMP,102) 
 239 :         WRITE(IOIMP,666) (DIN(i),i=1,NPO) 
 240 :   102   FORMAT(' ON VIENT DE PASSER LA BOUCLE 10') 
 241 :   666   FORMAT('   DIN ',/,(1X,10E12.5)) 
 242 :       ENDIF
 243 : C 
 244 : C **** ON BOUCLE SUR LES RIGIDITE ET ON TESTE: 
 245 : C **** SI LE MULTI EXISTE DANS DVAL SON SIGNE PAR RAPPORT A IRIGEL(6,I) 
 246 : C **** SINON ON TESTE LA RELATION ELLE MEME A L'AIDE DES COEFF 
 247 : C **** DE LA MATRICE ET DE LA VALEUR DU LAMBDA INI (DANS DIN) 
 248 : C **** ON CREE UN TABLEAU CONTENANT LE NUMERO DES MATRICES A GARDER ET 
 249 : C **** LE NUMERO DE LA SOUS RIGIDITE 
 250 :       IPA=0 
 251 :       JG=NBO 
 252 :       SEGINI MLENTI,trav1 
 253 :       DO 313 I=1,IRIGEL(/2)
 254 :         ITY=IRIGEL(6,I)
 255 :         IF (ITY.EQ.0) GO TO 313
 256 :         MELEME=IRIGEL(1,I) 
 257 :         SEGACT MELEME 
 258 :         DESCR=IRIGEL(3,I) 
 259 :         SEGACT DESCR 
 260 :         DO 314 J=1,LISINC(/2)
 261 :           IF(LISINC(J).EQ.'LX  ') GO TO 315
 262 :   314   CONTINUE
 262 :         CALL ERREUR(5) 
 263 :         RETURN 
 264 :   315   CONTINUE
 265 :         JJ=NOELEP(J) 
 266 :         ITEMP=LISINC(/2) 
 267 :         SEGINI IPASS 
 268 :         DO 316 J=1,ITEMP
 269 :           DO 317 K=1,NIN
 270 :             IF(LISINC(J).NE.NOMIN(K)) GO TO 317
 271 :             IPASS(J)=K 
 272 :             GO TO 316
 273 :   317     CONTINUE
 273 :           CALL ERREUR (5) 
 274 :           RETURN 
 275 :   316   CONTINUE
 276 : *  RECHERCHE DU MAX DES LAMDAS POUR LE CRITERE DE DECOLLEMENT D'APPUI 
 277 : *  RECHERCHE DU MAX DES DEPLACEMENT POUR LE CRITERE DE TRAVERSEE D'APPUI
 278 : * ymcrit est le deplacement maxi de tous les points en relation unilateral
 279 : * xmcrit est le maxi des multiplicateurs existant dans le chpoin de depla
 280 : * remarque : ce chpoint de depla contient les multiplicateurs de contact
 281 : * (pression) 
 282 :         DO 30 J=1,NUM(/2) 
 283 :           IP=ICPR(NUM(JJ,J))
 284 :           DO 31 K=1,ITEMP
 285 :             IF(LISINC(K).EQ.'LX  ') GO TO 31 
 286 :             IPPP=NUM(NOELEP(K),J) 
 287 :             IPP=ICPR(IPPP) 
 288 :             YMCRIT=MAX(YMCRIT,ABS(DVAL(IPASS(K),IPP)))
 289 :    31     CONTINUE 
 290 : *  on rajoute le jeu dans le critere
 291 :           YMCRIT=MAX(YMCRIT,ABS(DIN(IP)))
 292 :           IF (IEXI(1,IP).EQ.0) GOTO 30 
 293 :           XMCRIT=MAX(XMCRIT,ABS(DVAL(1,IP))) 
 294 :   30    CONTINUE 
 295 :         SEGSUP IPASS
 296 : * on rajoute dans ymcrit la dimension de l'element
 297 :         do 32 j=1,num(/2)
 298 :           xr=0.
 299 :           yr=0.
 300 :           zr=0.
 301 :           xp=0.
 302 :           yp=0.
 303 :           zp=0.
 304 :           DO 190 K=1,ITEMP
 305 :             IF (LISINC(K).EQ.'LX  ') GO TO 190
 306 :             IPPP=NUM(NOELEP(K),J)
 307 :             ipp = idimp1*(ippp-1)
 308 :             xr = xcoor(ipp+1)
 309 :             yr = xcoor(ipp+2)
 310 :             if (idim.ge.3) zr = xcoor(ipp+3)
 311 :             goto 191
 312 :  190      continue
 313 :  191      continue
 314 :           DO 192 K=1,LISINC(/2)
 315 :             IF (LISINC(K).EQ.'LX  ') GO TO 192
 316 :             IPPP=NUM(NOELEP(K),J)
 317 :             ipp = idimp1*(ippp-1)
 318 :             xp = xcoor(ipp+1)
 319 :             yp = xcoor(ipp+2)
 320 :             if (idim.ge.3) zp=xcoor(ipp+3)
 321 :             ymcrit=max(ymcrit,abs(xp-xr),abs(yp-yr),abs(zp-zr))
 322 :  192      continue
 323 :   32    continue
 324 :  313  CONTINUE
 325 : *
 326 :       if (mchpo3.ne.0) then
 327 : *        write (6,*) ' excit1 xmcrit avant ',xmcrit
 328 :         do ip=1,ycpr(/1)
 329 :           xmcrit=max(xmcrit,abs(ycpr(ip)))
 330 :         enddo
 331 : *        write (6,*) ' excit1 xmcrit apres ',xmcrit
 332 :       endif
 333 : 
 334 : * test pv
 335 : *     do 320 j2=dval(/2)
 336 : *     do 320 j1=2,dval(/1)
 337 : *       ymcrit=max(abs(dval(j1,j2)),ymcrit)
 338 : * 320  continue
 339 : 
 340 :       XMCRIT=1D-16*XMCRIT
 341 :       YMCRIT=1D-10*YMCRIT
 342 :       ymcrit=max(xmcrit,ymcrit)
 343 :       xmcrit=xmcrit*1d-4
 344 : * Critere en cas de frottement
 345 :       yfcrit = YMCRIT*ymult
 346 : *
 347 :       DO 13 I=1,IRIGEL(/2) 
 348 :         ITY=IRIGEL(6,I)
 349 :         IF (ITY.EQ.0) GO TO 13
 350 :         MELEME=IRIGEL(1,I) 
 351 :         SEGACT MELEME 
 352 :         DESCR=IRIGEL(3,I) 
 353 :         SEGACT DESCR 
 354 :         xMATRI=IRIGEL(4,I) 
 355 :         SEGACT xMATRI 
 356 :         ITEMP=LISINC(/2) 
 357 :         SEGINI IPASS 
 358 :         DO 14 J=1,ITEMP
 359 :           IF (LISINC(J).EQ.'LX  ') GO TO 15 
 360 :    14   CONTINUE 
 360 :         CALL ERREUR(5) 
 361 :         RETURN 
 362 :    15   CONTINUE 
 363 :         JJ=NOELEP(J) 
 364 :         DO 16 J=1,ITEMP
 365 :           DO 17 K=1,NIN 
 366 :             IF(LISINC(J).NE.NOMIN(K)) GO TO 17 
 367 :             IPASS(J)=K 
 368 :             GO TO 16 
 369 :    17     CONTINUE 
 369 :           CALL ERREUR (5) 
 370 :           RETURN 
 371 :    16   CONTINUE 
 372 :         DO 18 J=1,NUM(/2) 
 373 :           IPA=IPA+1 
 374 :           IPT=NUM(JJ,J)
 375 :           IP=ICPR(IPT)
 376 :           ityz=0
 377 : *  eliminer condition redondantes pour le pas d'apres
 378 :           if (kcpr.ne.0) then
 379 :             if (kcpr(num(noelep(1),j)).ne.0) then
 380 :               goto 18
 381 :             endif
 382 :           endif
 383 : C 
 384 : C  ** CAS OU LE BLOQUAGE N'ETAIT PAS SOLLICITE 
 385 : C 
 386 :           IF (IEXI(1,IP).EQ.0) THEN 
 387 :             SS=0.D0
 388 :             DO 19 K=1,ITEMP
 389 :               if (ipass(k).eq.0) goto 19
 390 :               IF (LISINC(K).EQ.'LX  ') GO TO 19 
 391 :               IPPP=NUM(NOELEP(K),J) 
 392 :               IPP=ICPR(IPPP) 
 393 :               if (ipp.eq.0) goto 19
 394 :               SS = DVAL(IPASS(K),IPP) * RE(1,K,j) + SS
 395 :    19       CONTINUE 
 396 :             IF (IIMPI.EQ.528) WRITE(IOIMP,529) IPPP,SS 
 397 : 529   FORMAT( ' LIBRE  ',I5,2X,E12.5) 
 398 :             r_z  = ABS(din(ip))*1.D-10
 399 :             r_p1 = DIN(IP) + r_z
 400 :             r_m1 = DIN(IP) - r_z
 401 : *  cas du frottement
 402 :             if (ity.eq.2) then
 403 :               ityz=jcpr(ipt)
 404 : *              write(ioimp,*) 'ipt jcpr ycpr ',ipt,jcpr(ipt),ycpr(ipt)
 405 :               if (ityz.ne.0) ityz=sign(1.1D0,ycpr(ipt))
 406 : *              if (ityz.eq.0) write (6,*) ' frottement -1 tyz ',ityz
 407 : *              write(ioimp,*) '1 dans excite ityz ',ityz
 408 : *  apparamment il faut etre plus laxiste pour le frottement
 409 : *  peut etre pas finalement
 410 :               if (ityz.eq.0) goto 20
 411 :               IF (ITYz.EQ. 1.AND.SS.LE.r_p1+yfcrit) GOTO 20
 412 :               IF (ITYz.EQ.-1.AND.SS.GE.r_m1-yfcrit) GOTO 20
 413 :             endif
 414 :             IF (ITY.EQ. 1.AND.SS.LE.r_p1+YMCRIT) GOTO 20
 415 :             IF (ITY.EQ.-1.AND.SS.GE.r_m1-YMCRIT) GOTO 20
 416 :             LECT(IPA)=1 
 417 :             viol(ipa)=ss-din(ip)-ymcrit
 418 :        if (iimpi.eq.1967)
 419 :      >     write (6,*) ' ss din ymcrit nouveau blocage '
 420 :      $   ,ss,din(ip),ymcrit,i,ipa,viol(ipa),ipt,ityz,ityp
 421 : *     on a un (1) nouveau blocage on arrete
 422 :    20       CONTINUE 
 423 : C 
 424 : C  ** CAS OU LE BLOQUAGE ETAIT SOLLICITE  PETIT PROBLEME DE TEST DE 
 425 : C     PRECISION SUR LA VALEUR DE LA REACTION 
 426 : C 
 427 :           ELSE 
 428 :             iav(ipa)=1
 429 :             SS=DVAL(1,IP) 
 430 :             IF(IIMPI.EQ.528) WRITE(IOIMP,530) NUM(JJ,J),SS 
 431 :   530 FORMAT('  BLOQUER ' ,I5,2X,E12.5) 
 432 : *  cas du frottement
 433 :             if (ity.eq.2) then
 434 :               ityz=jcpr(ipt)
 435 :               if (ityz.ne.0) ityz=sign(1.1D0,ycpr(ipt))
 436 : *              if (ityz.eq.0) write (6,*) ' frottement -2 ityz ',ityz
 437 : *  apparamment il faut etre plus laxiste pour le frottement
 438 : *  peut etre pas finalement
 439 :               if (ityz.eq.0) goto 21
 440 :               IF(ITYz.EQ.1.AND.SS.LT. XMCRIT) GO TO 21
 441 :               IF(ITYz.EQ.-1.AND.SS.GT.-XMCRIT) GO TO 21
 442 :             endif
 443 :             IF(ITY.EQ.1.AND.SS.LT. XMCRIT) GO TO 21
 444 :             IF(ITY.EQ.-1.AND.SS.GT.-XMCRIT) GO TO 21
 445 :             LECT(IPA)=1 
 446 :             goto 18
 447 :    21       CONTINUE 
 448 :             viol(ipa)=ss+xmcrit
 449 : C      write (6,*) ' ss xmcrit ',ss,xmcrit,ityz,i,ipa,
 450 : C     >  num(3,j),num(4,j),num(5,j)
 451 :       if (iimpi.eq.1967)
 452 :      > write (6,*) ' appui disparait ',ipa,ss,xmcrit,ipt,ityz,ityp
 453 :           ENDIF 
 454 :    18   CONTINUE 
 455 :         SEGSUP IPASS 
 456 :         SEGDES MELEME,DESCR,xMATRI 
 457 :    13 CONTINUE 
 458 :       IF(IIMPI.EQ.528) WRITE(IOIMP,*) 'ON VIENT DE PASSER LA BOUCLE 13'
 459 : C 
 460 : C  **** DANS LECT ON DIT SI LA JEEME RIGI ELEMTAIRE EST A CONSERVER 
 461 : C 
 462 :       NRIGEL=IRIGEL(/2) 
 463 :       SEGINI RI2
 464 : C 
 465 : C  **** CAS OU IL N'Y A RIEN A GARDER  ON CREE UNE RIGIDITE VIDE 
 466 : C  **** POUR POUVOIR CONTINUER D'ITERER SI IL Y A LIEU 
 467 : C 
 468 :       IF (NRIGEL.EQ.0) THEN
 469 :         IF (IIMPI.EQ.528) WRITE(IOIMP,*) ' IL N''Y A RIEN A CREER'
 470 :         GO TO 50
 471 :       ENDIF 
 472 : *
 473 : *** ne garder que ce qui viole de plus de xcrot du max
 474 : *
 475 : * recherche du max
 476 :       violmf=0.
 477 :       violmd=0.
 478 :       imf=0
 479 :       imd=0
 480 :       do 40 i=1,nbo
 481 :         violm=abs(viol(i))
 482 :         if (iav(i).eq.0.and.lect(i).ne.0) then
 483 :           if (violm.gt.violmd) then
 484 :             violmd=violm
 485 :             imd=i
 486 :           endif
 487 :         else if (iav(i).eq.1.and.lect(i).ne.1) then
 488 :           if (violm.gt.violmf) then
 489 :             violmf=violm
 490 :             imf=i
 491 :           endif
 492 :         endif
 493 :   40  continue
 494 : *      write (6,*) ' ityp xcrot violmd violmf ',ityp,xcrot,
 495 : *     >            violmd,violmf,imd,imf
 496 : ****  if (ityp.eq.3) then
 497 : ****    do 42 i=1,nbo
 498 : ****      if (iav(i).eq.0.and.imd.ne.i) lect(i)=iav(i)
 499 : ****      if (iav(i).eq.1.and.(imd.ne.0.or.imf.ne.i)) lect(i)=iav(i)
 500 : *42     continue
 501 : ****  else
 502 : **      rvd = xcrot*violmd
 503 : **      rvf = xcrot*violmf
 504 :         rvd = 0.250*violmd
 505 :         rvf = 0.707*violmf
 506 :         if (ityp.eq.3) then
 507 :          rvd=rvd*1d-14
 508 :          rvf=rvf*1d-14
 509 :         endif
 510 :         do 41 i=1,nbo
 511 :           violm=abs(viol(i))
 512 :           if (imf.ne.0.and.iav(i).eq.0) lect(i)=iav(i)
 513 : **        if (imd.ne.0.and.iav(i).eq.1) lect(i)=iav(i)
 514 :           if (iav(i).eq.0.and.violm.ge.rvd) goto 41
 515 : **?          if (iav(i).eq.1.and.violm.ge.rvf) goto 41
 516 :           if (iav(i).eq.1) goto 41
 517 :           lect(i)=iav(i)
 518 :   41    continue
 519 : ****  endif
 520 : C 
 521 : C  ** IL FAUT CREER UNE RIGIDITE 
 522 : C 
 523 :       IPA= 0 
 524 :       IRI=0 
 525 :       DO 25 I =1,IRIGEL(/2) 
 526 :         IF(IRIGEL(6,I).EQ.0) GO TO 25
 527 :         MELEME=IRIGEL(1,I) 
 528 :         SEGACT MELEME 
 529 :         NELRIG=0 
 530 :         DO 27 J=1,NUM(/2) 
 531 :           IF(LECT(IPA+J).EQ.1) NELRIG=NELRIG+1 
 532 :    27   CONTINUE 
 533 :         if (nelrig.eq.0) goto 26
 534 :         IRI=IRI+1
 535 :         xMATRI=IRIGEL(4,I) 
 536 :         SEGACT xMATRI 
 537 :         nligrp=re(/2)
 538 :         nligrd=re(/1) 
 539 :         SEGINI xMATR1 
 540 :         RI2.IRIGEL(4,IRI)=xMATR1 
 541 :         DESCR=IRIGEL(3,I) 
 542 :         SEGINI,DES1=DESCR 
 543 :         SEGDES DES1 
 544 :         RI2.IRIGEL(3,IRI)=DES1 
 545 :         RI2.IRIGEL(5,IRI)=IRIGEL(5,I) 
 546 :         RI2.IRIGEL(6,IRI)=IRIGEL(6,I) 
 547 :         RI2.COERIG(IRI)=COERIG(I) 
 548 :         NBNN=NUM(/1) 
 549 :         NBELEM=NELRIG 
 550 :         NBSOUS=0 
 551 :         NBREF=0 
 552 :         SEGINI IPT1 
 553 :         IPT1.ITYPEL=ITYPEL 
 554 :         RI2.IRIGEL(1,IRI)=IPT1 
 555 :         I2=0 
 556 :         DO 28 J=1,NUM(/2) 
 557 :           IF(LECT(IPA+J).EQ.0) GO TO 28 
 558 :           I2=I2+1 
 559 :           DO 29 K=1,NUM(/1) 
 560 :             IPT1.NUM(K,I2)=NUM(K,J) 
 561 :  29       CONTINUE 
 562 :           do io=1,nligrp
 563 :             do iu=1,nligrd
 564 :               xmatr1.re(iu,io,i2)=re(iu,io,j)
 565 :             enddo
 566 :           enddo 
 567 :  28     CONTINUE 
 568 :         SEGDES xMATR1,xMATRI,IPT1 
 569 :  26     CONTINUE 
 570 :         IPA=IPA+NUM(/2) 
 571 :         SEGDES MELEME 
 572 :    25 CONTINUE 
 573 : *
 574 :       if (iri.ne.ri2.irigel(/2)) then
 575 :          nrigel=iri
 576 :          segadj ri2
 577 :       endif
 578 : 
 579 :    50 CONTINUE 
 580 :       SEGDES,RI2,MRIGID 
 581 : 
 582 : *  indice de retour
 583 :       indic = .true.
 584 :       do 55 i = 1, nbo
 585 :        if (iav(i).ne.lect(i)) indic = .false.
 586 :   55  continue
 587 :       SEGDES MLENTI
 588 : 
 589 :       SEGSUP SNOMIN,ICPR,ITRAV,trav1
 590 :       if (mchpo3.ne.0) segsup jcpr,ycpr
 591 :       if (kcpr.ne.0) segsup kcpr 
 592 : 
 593 :       RETURN 
 594 :       END 
 595 : 
 596 :  
 597 :  
 598 :  

© Cast3M 2003 - All rights reserved.
Disclaimer