Esope source files

   1 : C GRAD1X    SOURCE    BP208322  12/12/20    21:15:19     7632           
   2 : C   
   3 :       subroutine GRAD1X (IMODEL,IVADEP,LRE,IVAGRA,NGRA,
   4 :      &                   MINT1,MINT2,IIPDPG,IOK)
   5 : c
   6 : c     CALCUL DU GRADIENT DANS LE CAS D'ELEMENTS XFEM (MFR=63)
   7 : c
   8 : c     largement inspiré de epsix.eso + bgrmas.eso
   9 : C     
  10 : C
  11 : C*********************************************************
  12 : C       PARTIE DECLARATIVE
  13 : C*********************************************************
  14 : 
  15 : C
  16 :       IMPLICIT REAL*8(A-H,O-Z)
  17 : C
  18 :       PARAMETER  (NDDLMAX=30,NBNIMAX=10)
  19 : C
  20 :       PARAMETER  (NBENRMAX=5)
  21 :       DIMENSION  MLRE(NBENRMAX+1)
  22 : C
  23 : -INC CCOPTIO
  24 : -INC SMCOORD
  25 : -INC SMMODEL
  26 : -INC SMCHAML
  27 : -INC SMINTE
  28 : -INC SMELEME
  29 : -INC SMLREEL
  30 : C
  31 :       POINTEUR    MCHEX1.MCHELM,MINT1.MINTE,MINT2.MINTE
  32 : C
  33 :       SEGMENT WRK1
  34 :        REAL*8 XE(3,NBBB)
  35 :        REAL*8 XDDL(LRE),GRADI(NGRA)
  36 :       ENDSEGMENT
  37 : c
  38 :       SEGMENT WRK2
  39 :        REAL*8 SHPWRK(6,NBNO),BGR(NGRA,LRE)
  40 :        REAL*8 LV7WRK(NBENRMA2,2,6,NBBB)
  41 :       ENDSEGMENT
  42 : C
  43 :       SEGMENT MPTVAL
  44 :         INTEGER      IPOS(NS),NSOF(NS)
  45 :         INTEGER      IVAL(NCOSOU)
  46 :         CHARACTER*16 TYVAL(NCOSOU)
  47 :       ENDSEGMENT
  48 : C
  49 :       SEGMENT MRACC
  50 :          INTEGER TLREEL(NBENRMA2,NBI)
  51 :       ENDSEGMENT
  52 : C
  53 : C      write (*,*) '############################'
  54 : C      write (*,*) '#####  DEBUT DE GRAD1X  #####'
  55 : C      write (*,*) '############################'
  56 : C
  57 : C
  58 : C*********************************************************
  59 : c
  60 : C       RECUPERATION + ACTIVATIONS + VALEURS UTILES
  61 : c
  62 : C*********************************************************
  63 : c
  64 :       IOK = 0
  65 : C++++ Recup de la geometrie deja activée par grad1 ++++++
  66 :       MELEME= IMAMOD
  67 : C     nbre de noeuds par element, nbre d elements
  68 :       NBNN1 = NUM(/1)
  69 :       NBEL1 = NUM(/2)
  70 : 
  71 : 
  72 : C++++ RECUP DES INFOS EF ++++++++++++++++++++++++++++++++
  73 :       MELE = NEFMOD
  74 :       NGAU1 = MINT1.POIGAU(/1)
  75 : c C     sous decoupage et points de Gauss de l element geometrique de base
  76 : c       IF(MINT2.NE.0) SEGACT,MINT2
  77 : c       NGAU2 = MINT2.POIGAU(/1)
  78 : 
  79 : 
  80 : C++++ Recup des infos d enrichissement +++++++++++++++++++
  81 : c     recup du MCHEX1 d'enrichissement
  82 :       MCHAM1 = 0
  83 :       NOBMO1 = IVAMOD(/1)
  84 :       if (NOBMO1.ne.0) then
  85 :        do iobmo1=1,NOBMO1
  86 :         if ((TYMODE(iobmo1)).eq.'MCHAML') then
  87 :          MCHEX1 = IVAMOD(iobmo1)
  88 :          segact,MCHEX1
  89 :           if ((MCHEX1.TITCHE).eq.'ENRICHIS') then
  90 :             MCHAM1 = MCHEX1.ICHAML(1)
  91 :             segact,MCHAM1
  92 :             segdes,MCHEX1
  93 :             goto 1000
  94 :           endif
  95 :          segdes,MCHEX1
  96 :         endif
  97 :        enddo
  98 :        write(*,*) 'Le modele est vide (absence d enrichissement)'
  99 :       else
 100 :       write(*,*) 'Il n y a pas de MCHEML enrichissement dans le Modele'
 101 :       endif
 102 : 
 103 :  1000 continue
 104 : c     niveau d enrichissement(s) du modele (ddls std u exclus)
 105 : c     NBENR1= 0 si std, 1 si H seul,  2 si H+F1,  3 si H+F1+F2, etc...
 106 :       if (MCHAM1.ne.0) then
 107 :         NBENR1= MCHAM1.IELVAL(/1)
 108 :       else
 109 :         NBENR1 = 0
 110 :       endif
 111 : C      write(*,*) 'niveau d enrichissement(s) du modele',NBENR1
 112 : c
 113 : C
 114 : C
 115 : C*********************************************************
 116 : C       INITIALISATIONS...
 117 : C*********************************************************
 118 : c
 119 : c     preparation des tables avec:
 120 : 
 121 :       if (NBENR1.ne.0) then
 122 :       do ienr=1,NBENR1
 123 : c        -le nombre d'inconnues de chaque sous-zone
 124 : c         determinee depuis le nombre de fonction de forme
 125 : c ienr=  1: U+H(1+1=2), 2: U+H+F1(2+4=6), 3: U+H+F1+F2(6+4=10)
 126 :          nbniJ = 2 + ((ienr-1)*4)
 127 :          MLRE(1+ienr) = IDIM*NBNN1*nbniJ
 128 :       enddo
 129 :       endif
 130 : C     Tables + longues car 1er indice correspond au fontion de forme std
 131 :       MLRE(1)    = IDIM*NBNN1*1
 132 : 
 133 :       if (NBENR1.lt.(NBENRMAX+1)) then 
 134 :         do ienr=(NBENR1+1),(NBENRMAX)
 135 :            MLRE(1+ienr) = 0
 136 :         enddo
 137 :       endif
 138 : c
 139 : c      ...DU SEGMENT WRK1
 140 :       NBENRMA2 = NBENRMAX
 141 :       NBBB = NBNN1
 142 :       SEGINI,WRK1
 143 : 
 144 : c      ...DU SEGMENT WRK2
 145 : c      NBNO = NBNI
 146 :       NBNO = LRE/IDIM
 147 :       SEGINI,WRK2
 148 : C
 149 : c      ...DU SEGMENT MRACC
 150 :       NBENRMA2 = NBENRMAX
 151 :       NBI = NBNN1
 152 :       segini,MRACC
 153 : C
 154 : C     du nombre d erreur sur le noms de composantes
 155 :       NBERR1=0
 156 : 
 157 : 
 158 : 
 159 : C*********************************************************
 160 : C
 161 : C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>  BOUCLE SUR LES ELEMENTS
 162 : C
 163 :       DO 2000 J=1,NBEL1
 164 : C      write(*,*) '========= EF',J,' /',NBEL1,'========='
 165 : C
 166 : C
 167 : C*********************************************************
 168 : C     POUR CHAQUE ELEMENT, ...
 169 : C*********************************************************
 170 : C
 171 : C++++ ON RECUPERE LES COORDONNEES DES NOEUDS DE L ELEMENT
 171 :       CALL DOXE(XCOOR,IDIM,NBNN1,NUM,J,XE)
 172 : C
 173 : c
 174 : C++++ NBENRJ = niveau d enrichissement de l element ++++
 175 : C     =0 si EF std   =1 si U+H   =2 si U+H+F1   =3 si U+H+F1+F2
 176 :       NBENRJ=0
 177 :       if(NBENR1.ne.0) then
 178 :       do IENR=1,NBENR1
 179 :         MELVA1 = MCHAM1.IELVAL(IENR)
 180 :         segact,MELVA1
 181 :         do I=1,NBNN1
 182 :          mlree1 = MELVA1.IELCHE(I,J)
 183 : C        on en profite pour remplir MRACC table de raccourcis pour cet element
 184 :          TLREEL(IENR,I) = mlree1
 185 :          if (mlree1.ne.0)  then
 186 :             NBENRJ = max(NBENRJ,IENR)
 187 : C           et on active la listreel
 188 :             segact,mlree1
 189 :          endif
 190 :         enddo
 191 :         segdes,MELVA1
 192 :       enddo
 193 :       endif
 194 :       if (NBENRMA2.gt.NBENR1) then
 195 :         do IENR2=(NBENR1+1),NBENRMA2
 196 :         do I=1,NBNN1
 197 :            TLREEL(IENR2,I) = 0
 198 :         enddo
 199 :         enddo
 200 :       endif
 201 : C
 202 : c
 203 : c++++ quelques indices pour dimensionner au plus juste
 204 : c     nbre total de ddl de l'élément considéré
 205 : c       NLIGRD = MLRE(1+NBENRJ)
 206 : c       NLIGRD = MLRE(1+NBENRJ)
 207 :       NDDL = MLRE(1+NBENRJ)
 208 : c     nbre fonction de forme=((Ni_std+Ni_enrichi)*nbnoeud)=Ninconnues/idim
 209 :       NBNI = (MLRE(1+NBENRJ)) / IDIM
 210 : *     nbre d inconnue / noeud
 211 :       NCOMP = NDDL/NBNN1
 212 : C      write(6,*) 'EF',J,' NBENRJ=',NBENRJ,
 213 : C     &'donc',NDDL,' ddls et ',NBNI,' fonctions de forme'
 214 : c
 215 : c
 216 : C
 217 : C++++ ON CALCULE XDDL ++++
 218 :       MPTVAL = IVADEP
 219 :       NCOSOU = IVAL(/1)
 220 : C
 221 :       INITOT = 0
 222 : C-->> BOUCLE SUR LES niveaux d'ENRICHISSEMENTS (0:U 1:A 2:BCDE1 3:BCDE2)
 223 :       DO IENR=0,NBENRJ
 224 : *nbre de fonction(s) de ce niveau d'enrichissement (=1 si std ou H, =4 pour F1,2,...)
 225 :         if (IENR .le. 1) then
 226 :          NBNIENR = 1
 227 :         else
 228 :          NBNIENR = 4
 229 :         endif
 230 : C---->> BOUCLE SUR LES fonctions de forme de ce type d'enrichissement
 231 :         DO INI=1,NBNIENR
 232 :           INITOT = INITOT + 1
 233 : C------>> BOUCLE SUR LA DIMENSION
 234 :           DO 2220 KDIM=1,IDIM
 235 :             ICOMP = (INITOT-1)*IDIM + KDIM
 236 :             
 237 : c         --cas ou on n'a pas trouvé assez de composantes--
 238 :             if(ICOMP.GT.NCOSOU) GOTO 2221
 239 :             
 240 :             MELVAL = IVAL(ICOMP)
 241 : c         --cas ou on a pas trouvé cette composante dans cette zone du 
 242 : c           chpoint solution devenu mchaml  --
 243 :             if(MELVAL.eq.0)then
 244 : c             Avait on besoin de cette composante?
 245 : c             oui, si c'est une composante obligatoire
 246 :               if(IENR.eq.0) goto 2991
 247 : c             oui, si l'un des noeuds est enrichi
 248 :               do I=1,NBNN1
 249 :                 if(TLREEL(IENR,I).ne.0) goto 2991
 250 :               enddo
 251 : c             non, si c'est facultatif et qu'on n'est pas enrichi -> on saute
 252 :               goto 2220
 253 : c           ->AVERTISSEMENT puis on saute 
 254 :  2991         NBERR1=NBERR1+1
 255 :               if(IIMPI.lt.1) goto 2220
 256 : c               write(IOIMP,*) 'PB OPERATEUR GRAD :'
 257 :               write(IOIMP,991) ICOMP,IENR,INI,KDIM
 258 :   991         format(2X,'ABSENCE DANS LE CHPOINT DEPLACEMENT DE LA ',I3,
 259 :      $                  'ieme COMPOSANTE (enrichissement',I3,
 260 :      $                   ', fonction',I3,', direction ',I3,
 261 :      $                   ') NECESSAIRE POUR L UN DES NOEUDS SUIVANTS :')
 262 :               write(IOIMP,*)'  noeuds :',(NUM(iou,J),iou=1,NBNN1)
 263 :               goto 2220
 264 :             endif
 265 : 
 266 : C---------->> BOUCLE SUR LES NOEUDS
 267 :             DO I=1,NBNN1
 268 :               IQ = ((INITOT-1)*NBNN1*IDIM) + ((I-1)*IDIM) + KDIM
 269 :               XDDL(IQ) = VELCHE(I,J)
 270 :             ENDDO
 271 : 
 272 :  2220     CONTINUE
 273 :         ENDDO
 274 :       ENDDO
 275 :       
 276 : c   --cas normal (toutes les composantes souhaitees etaient presentes)--
 277 :       GOTO 2223      
 278 : 
 279 : c   --cas ou on n'a pas trouvé assez de composantes--
 280 :  2221 CONTINUE
 281 :       if (IIMPI.ge.1) then
 282 :         WRITE(IOIMP,2222) J,NCOSOU,NCOMP
 283 :  2222   FORMAT('OPERATEUR GRAD : ELEMENT ',I6,       
 284 :      &  ' LE CHAMP DE DEPLACEMENT CONTIENT ',I3,' COMPOSANTES',
 285 :      &  ' PAR NOEUD AU LIEU DE ',I3,' ATTENDUES')
 286 :       endif
 287 :       NDDL=NCOSOU*NBNN1
 288 :       NBENRJ=IENR
 289 :         
 290 :  2223 CONTINUE
 291 : c   --cas ou on a une ou des erreurs--
 292 :       IF (NBERR1.gt.0.and.J.eq.NBEL1) THEN
 293 :         write(IOIMP,*) 'OPERATEUR GRAD : ABSENCE DANS LE CHPOINT ',
 294 :      &   'DEPLACEMENT DE CERTAINES INCONNUES ATTENDUES PAR LE MODELE'
 295 :       ENDIF
 296 : 
 297 : c
 298 : c
 299 : c rem: il serait probablement interessant au niveau du tems cpu 
 300 : c  d'utiliser moins de pts de Gauss lorsque l element est élastique
 301 : c  On pourrait par ex. utiliser MINT2 = infell(12) qui contient 
 302 : c  le segment d'integration de l'EF std (QUA4 par ex.)
 303 : *      if((NBENRJ.eq.0).and.(MINT2.ne.0)) then
 304 : *         MINTE = MINT2
 305 : *         NBPGAU= NGAU2
 306 : *      else
 307 :          MINTE = MINT1
 308 :          NBPGAU= NGAU1
 309 : *      endif
 310 : c
 311 : C
 312 : C*********************************************************
 313 : C
 314 : C>>>>>>>>>>  BOUCLE SUR LES POINTS DE GAUSS
 315 : C
 316 :       DO 2500 KGAU=1,NBPGAU
 317 : C
 318 : C*********************************************************
 319 : C     Initialisation à 0
 320 : C*********************************************************
 321 : 
 322 : c ZERO ne serait-il pas facultatif?
 322 :       CALL ZERO(SHPWRK,6,NBNO)
 323 :       XGAU = 0.D0
 324 :       YGAU = 0.D0
 325 :       ZGAU = 0.D0
 326 : C
 327 :       i6zz = 3
 328 :       IF(IDIM.EQ.3) i6zz = 4
 329 : C
 330 :       do ienr7=1,NBENRJ
 331 :        do inod7=1,NBNN1
 332 :         do i6=1,i6zz
 333 :            LV7WRK(ienr7,1,i6,inod7) = 0.D0
 334 :            LV7WRK(ienr7,2,i6,inod7) = 0.D0
 335 :         enddo
 336 :        enddo
 337 :       enddo
 338 : 
 339 : 
 340 : c*********************************************************
 341 : c     Calcul des fonction de forme std dans repere local
 342 : c*********************************************************
 343 : 
 344 : ccccc BOUCLE SUR LES NOEUDS ccccccccccccccccccccccccccccc
 345 : c     (et donc sur les Ni std)
 346 :       DO 2510 I=1,NBNN1
 347 : 
 348 : C++++ Calcul des Ni std  
 349 : c     (rappel: 1:Ni  2:Ni,qsi  3:Ni,eta  avec i=1,4)
 350 :       SHPWRK(1,I) = SHPTOT(1,I,KGAU)
 351 :       SHPWRK(2,I) = SHPTOT(2,I,KGAU)
 352 :       SHPWRK(3,I) = SHPTOT(3,I,KGAU)
 353 :       IF(IDIM.EQ.3) SHPWRK(4,I) = SHPTOT(4,I,KGAU)
 354 : 
 355 : C++++ Calcul des coordonnees dans le repere global
 356 :       XGAU = XGAU + (SHPWRK(1,I)*XE(1,I))
 357 :       YGAU = YGAU + (SHPWRK(1,I)*XE(2,I))
 358 :       IF(IDIM.EQ.3) ZGAU = ZGAU + (SHPWRK(1,I)*XE(3,I))
 359 : 
 360 :  2510 CONTINUE
 361 : ccccc fin de BOUCLE SUR LES NOEUDS ccccccccccccccccccccccc
 362 : c      if(J.eq.77) write(6,*) 'xgau_',KGAU,' =',XGAU,YGAU,ZGAU
 363 : 
 364 : 
 365 : 
 366 : c*********************************************************
 367 : c     Passage des fonctions de forme std dans repere global
 368 : c*********************************************************
 369 : 
 370 : C++++ CALCUL DES  Ni,x Ni,y (i=1,4) + CALCUL DE det(J)
 370 :       CALL JACOBI(XE,SHPWRK,IDIM,NBNN1,DJAC)
 371 : c       if(J.eq.1.and.KGAU.eq.1) 
 372 : c      & write(*,*) 'Ni(i=1,4)=',(SHPWRK(1,iou),iou=1,NBNN1)
 373 : 
 374 : c*********************************************************
 375 : c     Si on est pas du tout enrichi on peut sauter une grosse partie
 376 : c*********************************************************
 377 :       if(NBENRJ.eq.0) goto 2999
 378 : 
 379 : c*********************************************************
 380 : c     Calcul des level set + leurs derivees dans repere global
 381 : c*********************************************************
 382 : 
 383 : ccccc BOUCLE SUR LES enrichissements ccccccccccccccccccc
 384 :       do 2520 ienr=1,NBENRJ
 385 : 
 386 : c       MELVA1=MCHAM1.IELVAL(IENR)
 387 : c       segact,MELVA1
 388 : 
 389 : ccccc  BOUCLE SUR LES NOEUDS ccccccccccccccccccccccccccc
 390 :        DO 2521 I=1,NBNN1
 391 : 
 392 : C++++ Le I eme noeud est-il ienr-enrichi?
 393 :         mlree1= TLREEL(IENR,I)
 394 : c        if(J.eq.77) write(6,*) 'ienr,I',ienr,I,' mlree1=',mlree1
 395 :         if(mlree1.eq.0)  goto 2521
 396 : 
 397 : 
 398 : c       Calcul du repere local de fissure(=PSI,PHI)
 399 : c       (rappel: 1,1:psi  1,2:phi  2,1 psi,x  ...etc...)
 400 :         do 2522 inode=1,NBNN1
 401 : c         pour le H-enrichissement, on n a pas gardé PSI (inutile)
 402 :           if (ienr.ne.1) then
 403 : c          valeur de PSI au inode^ieme noeud
 404 :            xpsi1 = mlree1.PROG(inode)
 405 : c          qu on multiplie par la valeur de Ni^std au pt de G considéré
 406 :            LV7WRK(ienr,1,1,I)= LV7WRK(ienr,1,1,I) 
 407 :      &      + (SHPWRK(1,inode)*xpsi1)
 408 :            LV7WRK(ienr,1,2,I)= LV7WRK(ienr,1,2,I) 
 409 :      &      + (SHPWRK(2,inode)*xpsi1)
 410 :            LV7WRK(ienr,1,3,I)= LV7WRK(ienr,1,3,I) 
 411 :      &      + (SHPWRK(3,inode)*xpsi1)
 412 :       IF(IDIM.EQ.3) LV7WRK(ienr,1,4,I)= LV7WRK(ienr,1,4,I) 
 413 :      &      + (SHPWRK(4,inode)*xpsi1)
 414 : c          valeur de PHI au inode^ieme noeud
 415 :            xphi1 = mlree1.PROG(NBNN1+inode)
 416 :           else
 417 :            xphi1 = mlree1.PROG(inode)
 418 :           endif
 419 :           LV7WRK(ienr,2,1,I)= LV7WRK(ienr,2,1,I) 
 420 :      &      + (SHPWRK(1,inode)*xphi1)
 421 :           LV7WRK(ienr,2,2,I)= LV7WRK(ienr,2,2,I) 
 422 :      &      + (SHPWRK(2,inode)*xphi1)
 423 :           LV7WRK(ienr,2,3,I)= LV7WRK(ienr,2,3,I) 
 424 :      &      + (SHPWRK(3,inode)*xphi1)
 425 :       IF(IDIM.EQ.3) LV7WRK(ienr,2,4,I)= LV7WRK(ienr,2,4,I) 
 426 :      &      + (SHPWRK(4,inode)*xphi1)
 427 :  2522   continue
 428 : 
 429 :  2521  continue
 430 : ccccc  fin de BOUCLE SUR LES NOEUDS ccccccccccccccccccccccc
 431 : 
 432 : 
 433 :  2520 CONTINUE
 434 : ccccc fin de BOUCLE SUR LES enrichissements cccccccccccccccc
 435 : 
 436 : c     on a construit
 437 : C     LV7WRK(ienr, PSI/PHI, valeur/deriveeparqsi/pareta, iNOEUD)
 438 : 
 439 : 
 440 : c*********************************************************
 441 : c     Ajout des fonctions de forme d enrichissement
 442 : c     + leurs derivees dans repere global
 443 : c*********************************************************
 443 :       CALL SHAPX(XGAU,YGAU,ZGAU,MELE,LV7WRK,NBENRMAX,NBENRJ,
 444 :      &TLREEL,J,SHPWRK,IRET)
 445 : c      if(J.eq.77) write(6,*) 'retour de shapx'
 446 :      
 447 : 
 448 : 
 449 : c     retour a la partie commune aux EF enrichi et non enrichi
 450 :  2999 continue
 451 : 
 452 : C*********************************************************
 453 : C       CALCUL DE BGR
 454 : C*********************************************************
 455 : c ZERO ne serait-il pas facultatif?
 455 :         call ZERO(BGR,9,NDDL)
 456 :         KB=1
 457 : C       boucle sur tous les Ni
 458 :         DO 3001 II=1,NBNI
 459 : 
 460 : c         partie commune 2D 3D
 461 :           BGR(1,KB)   = SHPWRK(2,II)
 462 :           BGR(2,KB)   = SHPWRK(3,II)
 463 :           BGR(4,KB+1) = SHPWRK(2,II)
 464 :           BGR(5,KB+1) = SHPWRK(3,II)
 465 : 
 466 : c         cas 3D
 467 :           IF (IDIM.EQ.3) THEN
 468 :             BGR(3,KB)   = SHPWRK(4,II)
 469 :             BGR(6,KB+1) = SHPWRK(4,II)
 470 :             BGR(7,KB+2) = SHPWRK(2,II)
 471 :             BGR(8,KB+2) = SHPWRK(3,II)
 472 :             BGR(9,KB+2) = SHPWRK(4,II)
 473 :           ENDIF
 474 : 
 475 : c         cas 2D PLAN GENE
 476 :           IF (IFOUR.EQ.-3) THEN
 477 :             IREF=(IIPDPG-1)*(IDIM+1)
 478 :             BGR(9,KB)   = 1.D0
 479 :             BGR(9,KB+1) = XCOOR(IREF+1)-XGAU
 480 :             BGR(9,KB+2) = YGAU-XCOOR(IREF+2)
 481 :           ENDIF
 482 : 
 483 :           KB = KB + IDIM
 484 : 
 485 :  3001   CONTINUE
 486 : C
 487 : c
 488 : 
 489 : C*********************************************************
 490 : C       CALCUL DE  BGR.q
 491 : C*********************************************************
 492 : 
 493 : C  APPEL A LA PROCEDURE DE CALCUL
 493 :         CALL BGRDEP(BGR,NGRA,XDDL,NDDL,GRADI)
 494 : c
 495 : 
 496 : C*********************************************************
 497 : C       ECRITURE DES DIFFERENTES COMPOSANTES DU GRADIENT
 498 : C*********************************************************
 499 :         MPTVAL=IVAGRA
 500 :         DO i=1,NGRA
 501 :           MELVAL=IVAL(i)
 502 :           IGMN=MIN(KGAU,VELCHE(/1))
 503 :           IBMN=MIN(J,VELCHE(/2))
 504 :           VELCHE(IGMN,IBMN)=GRADI(i)
 505 :         ENDDO
 506 : c
 507 : c
 508 :  2500 CONTINUE
 509 : C FIN DE BOUCLE SUR LES POINTS DE GAUSS <<<<<<<<<<<<<<
 510 : C
 511 : c
 512 : c     quelques desactivations
 513 :       do IENR=1,NBENR1
 514 :         do I=1,NBNN1
 515 :          mlree1 = TLREEL(IENR,I)
 516 :          if(mlree1.ne.0)  segdes,mlree1
 517 :         enddo
 518 :       enddo
 519 : c
 520 : c
 521 : c
 522 :  2000 CONTINUE
 523 : C FIN DE BOUCLE SUR LES ELEMENTS <<<<<<<<<<<<<<<<<<<<<
 524 : c
 525 : c
 526 : 
 527 : C*********************************************************
 528 : C       SUPPRESSION ET DESACTIVATION DE SEGMENTS
 529 : C*********************************************************
 530 : C
 531 :       SEGSUP,WRK1,WRK2
 532 :       SEGSUP,MRACC
 533 : 
 534 : c      segdes,MINT1 fait dans grad1
 535 : C       if(MINT2.ne.0)  segdes,MINT2
 536 : c     on met mint2 a zero pour eviter pb dans grad1
 537 :       MINT2 = 0
 538 : c
 539 : c     pas de retour en erreur pour linstant
 540 :       IOK=1
 541 : c      write(6,*) 'iok=',IOK
 542 : C
 543 :       RETURN
 544 :       END
 545 :  
 546 :  
 547 :  
 548 :  
 549 :