C IMPOS2    SOURCE    FANDEUR   12/03/15    21:26:27     7314           

*  impo bloc en 2D

      SUBROUTINE IMPOS2

      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)

-INC CCOPTIO

-INC SMCOORD
-INC SMELEME
-INC SMRIGID
-INC SMCHPOI

      segment mfopa1
        integer lpoin(nfopa1,3)
        real*8  xpoin(nfopa1,2)
      endsegment

      character*4 modepl(4),moforc(4)

      data modepl /'UX  ','UY  ','UR  ','UZ  '/
      data moforc /'FX  ','FY  ','FR  ','FZ  '/
*
      idimp1 = IDIM + 1
*
*  Lecture du maillage support des conditions de contact
*
      call lirobj('MAILLAGE',ipt1,1,iretou)
      if (ierr.ne.0) return
*
*  Activation du maillage et petites verifications
*
      segact ipt1
      nbno1 = ipt1.num(/1)
      nbel1 = ipt1.num(/2)
      if (ipt1.lisous(/1).ne.0) call erreur(25)
      if (ipt1.itypel.ne.22) call erreur(16)
      if (nbno1.ne.4 .and. nbno1.ne.5) call erreur(16)
      if (ierr.ne.0) goto 900
*
*  Indice des ddls concernes en fonction du mode 2D
*
      imo = 1
      if (ifour.eq.0) imo = 3
*
*   Increment sur le nombre d'elements de contact retenus et utilises 
*   dans le chpoint (mpoval et igeoc) et rigidite pour agrandir les 
*   longueurs des segments adequats
      incnel = 2000
*
*  Creation de la raideur des conditions de contact
*  Remplissage de l'entete commun
*
      nrigel = 1
      segini mrigid
      ichole = 0
      imgeo1 = 0
      imgeo2 = 0
      isupeq = 0
      iforig = ifomod
      coerig(1)=1.d0
      mtymat='RIGIDITE'
*
*  Creation du chpoint de depi
*
      nat=1
      nsoupo=1
      segini mchpoi
      mtypoi='DEPIMP'
      mochde='engendré par impose'
      ifopoi=ifomod
      jattri(1)=2
*
      nc=1
      segini msoupo
      ipchp(1)=msoupo
      nocomp(1)='FLX '
      nocons(1)='                '
*
      nbnn  =1
      nbelem=incnel
      nbref =0
      nbsous=0
      segini ipt7
      ipt7.itypel=1
      igeoc=ipt7
*
      n=incnel
      segini mpoval
      ipoval=mpoval
*
*  Nombre de noeuds dans le chpoint (en totalite) : ipt7 et mpoval
      nelch = 0
*
*=======================================================================
* Formulation "faible" du contact : relation segment-segment (type 0)
*=======================================================================
      if (nbno1.ne.5) goto 500
*
*  Creation du meleme associe a la relation
*
      nbnn  =5
      nbelem=incnel
      nbsous=0
      nbref =0
      segini meleme
      itypel=22
      irigel(1,nrigel) = meleme
*
*  Creation du descriptif commun a toutes les raideurs
*
      nligrp=9
      nligrd=nligrp
      segini,descr
      lisinc(1)='LX  '
      lisdua(1)='FLX '
      noelep(1)=1
      noeled(1)=1
      do 100 i = 2, nligrp, 2
        lisinc(i  )=modepl(imo)
        lisinc(i+1)=modepl(imo+1)
        lisdua(i  )=moforc(imo)
        lisdua(i+1)=moforc(imo+1)
        noelep(i  )=(i+2)/2
        noelep(i+1)=noelep(i)
        noeled(i  )=noelep(i)
        noeled(i+1)=noelep(i)
 100  continue
      segdes,descr
      irigel(3,nrigel) = descr
*
*     creation du segment xmatri
*
      nelrig=incnel
      segini xmatri
      irigel(4,nrigel) = xmatri
*
*     ce qu'on cree est unilateral
*
      irigel(6,nrigel) = 1
*
*     ce qu'on cree est symetrique
*
      irigel(7,nrigel) = 0
*
*  Nombre d'elements crees dans meleme=irigel(nrigel,1), ipt7 et mpoval
      nelri0 = 0
*
*  boucle sur le maillage support
*
      do 110 iel = 1, nbel1
*
        ip1 = ipt1.num(2,iel)
        ip2 = ipt1.num(3,iel)
        ip3 = ipt1.num(4,iel)
        ip4 = ipt1.num(5,iel)
*        write(ioimp,*) iel,ip1,ip2,ip3,ip4,ipt1.num(1,iel)
*  verification (provisoire) que pas de noeuds doubles
        if (ip1.eq.ip3) goto 110
        if (ip1.eq.ip4) goto 110
        if (ip2.eq.ip3) goto 110
        if (ip2.eq.ip4) goto 110
*
        ipv = (ip1-1)*idimp1
        xp1 = xcoor(ipv+1)
        yp1 = xcoor(ipv+2)
        ipv = (ip2-1)*idimp1
        xp2 = xcoor(ipv+1)
        yp2 = xcoor(ipv+2)
        xp12 = xp2 - xp1
        yp12 = yp2 - yp1
*
        ipv = (ip3-1)*idimp1
        xp3 = xcoor(ipv+1)
        yp3 = xcoor(ipv+2)
        ipv = (ip4-1)*idimp1
        xp4 = xcoor(ipv+1)
        yp4 = xcoor(ipv+2)
        xp34 = xp4 - xp3
        yp34 = yp4 - yp3
*
*   orientations respectives correctes des 2 segments :
*   "normales de sens opposes" = produit scalaire negatif ou nul
        scal = xp12*xp34 + yp12*yp34
*        write(ioimp,*) xp12,yp12,xp34,yp34,scal
        if (scal.gt.0.) goto 110
*
*   critere d'acceptation de l'élément :
*   angles des diagonales
        xl13 = sqrt((xp3-xp1)**2+(yp3-yp1)**2)
        xl14 = sqrt((xp4-xp1)**2+(yp4-yp1)**2)
        xl23 = sqrt((xp3-xp2)**2+(yp3-yp2)**2)
        xl24 = sqrt((xp4-xp2)**2+(yp4-yp2)**2)
        sca312 = (xp3-xp1)*(xp3-xp2)+(yp3-yp1)*(yp3-yp2)
        sca412 = (xp4-xp1)*(xp4-xp2)+(yp4-yp1)*(yp4-yp2)
        sca134 = (xp1-xp3)*(xp1-xp4)+(yp1-yp3)*(yp1-yp4)
        sca234 = (xp2-xp3)*(xp2-xp4)+(yp2-yp3)*(yp2-yp4)
*        write(ioimp,*) sca312/(xl13*xl23),sca412/(xl14*xl24),
*    &                  sca134/(xl13*xl14),sca234/(xl23*xl24)
        if (sca312/(xl13*xl23).gt.0.7 .and.
     >      sca412/(xl14*xl24).gt.0.7 .and.
     >      sca134/(xl13*xl14).gt.0.7 .and.
     >      sca234/(xl23*xl24).gt.0.7) goto 110
*
* Quelle est la relation ???
*
*   direction du contact unitaire
*
        xr = xp12 - xp34
        yr = yp12 - yp34
        sr = sqrt(xr*xr + yr*yr)
        xr = xr/sr
        yr = yr/sr
*        write(ioimp,*) 'direction contact',xr,yr,yr,-xr
*
*   projection des points sur la direction du contact
*
        xl1 = xp1*xr + yp1*yr
        xl2 = xp2*xr + yp2*yr
        xl3 = xp3*xr + yp3*yr
        xl4 = xp4*xr + yp4*yr

*        write(ioimp,*) 'projection pts sur contact',xl1,xl2,xl3,xl4
*
*   calcul de l'intersection des projections
*
        xm1 = min(xl1,xl2)
        xm2 = max(xl1,xl2)
        xm3 = min(xl3,xl4)
        xm4 = max(xl3,xl4)
*        write(ioimp,*) ' xmi',xm1,xm2,xm3,xm4
*
*  critere de precision sur l'intersection
*
        xcr = max(xm2-xm1,xm4-xm3)*1.d-4
*
*  taille de l'intersection
*
        xi = max(xm1,xm3)
        xj = min(xm2,xm4)
        xl = xj - xi
*        write(ioimp,*) ' intersection',xi,xj,xl,xcr
*
        if (xl.le.xcr) goto 110
*
*   distance des points a leur projection
*
        d1 = (xp1-xl1*xr)*yr-(yp1-xl1*yr)*xr
        d2 = (xp2-xl2*xr)*yr-(yp2-xl2*yr)*xr
        d3 = (xp3-xl3*xr)*yr-(yp3-xl3*yr)*xr
        d4 = (xp4-xl4*xr)*yr-(yp4-xl4*yr)*xr
*
*  coordonnées paramétriques de l'intersection sur les segments 1-2 et 3-4
*
        xi2 = (xi-xl1) / (xl2-xl1)
        xi1 = (xl2-xi) / (xl2-xl1)

        xj2 = (xj-xl1) / (xl2-xl1)
        xj1 = (xl2-xj) / (xl2-xl1)

        xi4 = (xi-xl3) / (xl3-xl4)
        xi3 = (xl4-xi) / (xl3-xl4)

        xj4 = (xj-xl3) / (xl3-xl4)
        xj3 = (xl4-xj) / (xl3-xl4)
*
*        write(ioimp,*) ' xi1 xi2 xi3 xi4 xj1 xj2 xj3 xj4 xl '
*        write(ioimp,*)  xi1,xi2,xi3,xi4,xj1,xj2,xj3,xj4,xl
*        write(ioimp,*) ' d1 d2 d3 d4 ',d1,d2,d3,d4
*
* surface actuelle
*
        sc = ((xi1+xj1)*d1+(xi2+xj2)*d2+(xi3+xj3)*d3+(xi4+xj4)*d4)*0.5
*        write(ioimp,*) 'Surface actuelle :',sc
*
* on a un element ou imposer la relation a ajouter
*
        nelri0 = nelri0 + 1
        nelch  = nelch +1
*
*   on ajuste les differents segments
*
        if (nelri0.gt.nelrig) then
          nelrig = nelrig + incnel
          segadj,xmatri
          nbelem = nbelem + incnel
          nbnn = 5
          segadj,meleme
          nbnn = 1
          segadj,ipt7
          n = n + incnel
          segadj,mpoval
        endif
*
*   on range dans le meleme
*
        num(1,nelri0) = ipt1.num(1,iel)
        num(2,nelri0) = ip1
        num(3,nelri0) = ip2
        num(4,nelri0) = ip3
        num(5,nelri0) = ip4
*
*   on remplit le xmatri
*
*  lambda
        re(1,1,nelri0)= 0.d0
*  ip1
        re(2,1,nelri0)=  yr * (xi1+xj1) * 0.5
        re(3,1,nelri0)= -xr * (xi1+xj1) * 0.5
*  ip2
        re(4,1,nelri0)=  yr * (xi2+xj2) * 0.5
        re(5,1,nelri0)= -xr * (xi2+xj2) * 0.5
*  ip3
        re(6,1,nelri0)=  yr * (xi3+xj3) * 0.5
        re(7,1,nelri0)= -xr * (xi3+xj3) * 0.5
*  ip4
        re(8,1,nelri0)=  yr * (xi4+xj4) * 0.5
        re(9,1,nelri0)= -xr * (xi4+xj4) * 0.5
*
*  on transpose
        do 120 ic = 2, nligrp
          re(1,ic,nelri0)=re(ic,1,nelri0)
 120    continue
*
*  Le reste est nul
*
*  remplissage du champoint de depi et du maillage support
*
        ipt7.num(1,nelch)=ipt1.num(1,iel)
        vpocha(nelch,1) = -sc
*
 110   continue

*     write (ioimp,*) ' nb relation type 0 ',nbelem,n,nelri0

* Ajustement au plus juste puis desactivation des segments lies
* a la rigidite du type 0
      if (nelri0.ne.nelrig) then
        nelrig = nelri0
        segadj,xmatri
        nbelem = nelri0
        nbnn = 5
        segadj,meleme
      endif
      segdes,meleme,xmatri
*
* S'il n'y a pas d'elements en contact, alors pas de relation unilaterale
*     if (nelri0.eq.0) irigel(6,nrigel)=0
*
      GOTO 1000

*=======================================================================
* Formulation "forte" (standard) du contact :
*=======================================================================
 500  continue
*
*  Nombre de "cas pathologiques" de type noeud-noeud (type 1)
      nelpath = 0
      nfopa1  = incnel
      segini,mfopa1
*
* Relation type 2 : noeud-segment
*---------------------------------
*     creation du meleme associe a la relation
*
      nbnn   = 4
      nbelem = incnel
      nbsous = 0
      nbref  = 0
      segini meleme
      itypel=22
      irigel(1,nrigel)=meleme
*
*     creation du descriptif commun a toutes les raideurs
*
      nligrp=7
      nligrd=nligrp
      segini descr
      lisinc(1)='LX  '
      lisdua(1)='FLX '
      noelep(1)=1
      noeled(1)=1
      do 510 i=2,nligrp,2
        lisinc(i  )=modepl(imo)
        lisinc(i+1)=modepl(imo+1)
        lisdua(i  )=moforc(imo)
        lisdua(i+1)=moforc(imo+1)
        noelep(i  )=(i+2)/2
        noelep(i+1)=noelep(i)
        noeled(i  )=noelep(i)
        noeled(i+1)=noelep(i)
 510  continue
      segdes,descr
      irigel(3,nrigel)=descr
*
*     creation du segment xmatri
*
      nelrig=incnel
      segini xmatri
      irigel(4,nrigel)=xmatri
*
*     ce qu'on cree est unilateral
*
      irigel(6,nrigel)=1
*
*     ce qu'on cree est symetrique
*
      irigel(7,nrigel)=0
*
*  Nombre d'elements dans la rigidite de type 2
      nelri2=0
*
*  boucle sur le maillage support
*
      do 520 iel = 1,nbel1
        ip1 = ipt1.num(2,iel)
        ip2 = ipt1.num(3,iel)
        jp  = ipt1.num(4,iel)
*D       write(ioimp,*) iel,ip1,ip2,jp
*
*  verification que pas relation du point sur lui meme
        if (jp.eq.ip1) goto 520
        if (jp.eq.ip2) goto 520

        ipv = (ip1-1)*idimp1
        xp1 = xcoor(ipv+1)
        yp1 = xcoor(ipv+2)
        ipv = (ip2-1)*idimp1
        xp2 = xcoor(ipv+1)
        yp2 = xcoor(ipv+2)
        ipv = (jp-1)*idimp1
        xp  = xcoor(ipv+1)
        yp  = xcoor(ipv+2)

        x12 = xp2 - xp1
        y12 = yp2 - yp1
        x1p = xp - xp1
        y1p = yp - yp1
        x2p = xp - xp2
        y2p = yp - yp2

*  distance signee de p a la ligne 1-2
        dp12s = x1p * y12 - y1p * x12
*        write(ioimp,*) 'dist. signee',dp12s
*   verif que le point est du bon cote du segment
*       if (dp12s.lt.0.) goto 520

*  verification si autre point dans le cercle de selection
        d1p = x1p*x1p + y1p*y1p
        d2p = x2p*x2p + y2p*y2p
        scal = (x1p*x2p + y1p*y2p) / sqrt(d1p*d2p)
*        write(ioimp,*) 'cos(angle_1p2)',scal
        if (scal.gt.0.7) goto 520
*
* on a un point ou imposer la relation
* Quelle est la relation ???
*
*   direction de la relation
*
*   calcul du pt a mettre en relation avec ip : alpha ip2 + (1-alpha) ip1
*
        sr2 = x12*x12 + y12*y12
        alpha = (x1p*x12 + y1p*y12) / sr2
*
        sr = sqrt(sr2)
        xr = x12 / sr
        yr = y12 / sr
*  traitement des cas pathologiques (on est sur un point du segment)
        if (abs(alpha).le.1e-2) then
          nelpath = nelpath + 1
          if (nelpath.gt.nfopa1) then
            nfopa1 = nfopa1 + incnel
            segadj,mfopa1
          endif
          lpoin(nelpath,1)=ip1
          lpoin(nelpath,2)=jp
          lpoin(nelpath,3)=ipt1.num(1,iel)
          xpoin(nelpath,1)= yr
          xpoin(nelpath,2)=-xr
*          write(ioimp,*) ' cas pathologique 1 (2)',alpha,ip1,jp,yr,-xr
          goto 520
        else if (abs(1.-alpha).le.1e-2) then
          nelpath = nelpath + 1
          if (nelpath.gt.nfopa1) then
            nfopa1 = nfopa1 + incnel
            segadj,mfopa1
          endif
          lpoin(nelpath,1)=ip2
          lpoin(nelpath,2)=jp
          lpoin(nelpath,3)=ipt1.num(1,iel)
          xpoin(nelpath,1)= yr
          xpoin(nelpath,2)=-xr
*          write(ioimp,*) ' cas pathologique 1 (1)',alpha,ip1,jp,yr,-xr
          goto 520
        endif
*
* on a un element ou imposer la relation a ajouter
*
        nelri2 = nelri2 + 1
        nelch = nelch + 1
*
*   on ajuste les differents segments
*
        if (nelri2.gt.nelrig) then
          nelrig = nelrig + incnel
          segadj,xmatri
          nbelem = nbelem + incnel
          nbnn = 4
          segadj,meleme
          nbnn = 1
          segadj,ipt7
          n = n + incnel
          segadj,mpoval
        endif
*
*   on range dans le meleme
*
        num(1,nelri2)=ipt1.num(1,iel)
        num(2,nelri2)=ip1
        num(3,nelri2)=ip2
        num(4,nelri2)=jp
*
*   on remplit le xmatri
*
*  lambda
        re(1,1,nelri2)= 0.d0
*  ip1
        re(2,1,nelri2)=  yr * (1.-alpha)
        re(3,1,nelri2)= -xr * (1.-alpha)
*  ip2
        re(4,1,nelri2)=  yr * alpha
        re(5,1,nelri2)= -xr * alpha
*  jp
        re(6,1,nelri2)= -yr
        re(7,1,nelri2)=  xr
*  on transpose
        re(1,2,nelri2) = re(2,1,nelri2)
        re(1,3,nelri2) = re(3,1,nelri2)
        re(1,4,nelri2) = re(4,1,nelri2)
        re(1,5,nelri2) = re(5,1,nelri2)
        re(1,6,nelri2) = re(6,1,nelri2)
        re(1,7,nelri2) = re(7,1,nelri2)
*  le reste est nul
*
*  remplissage du champoint de depi
*
        ipt7.num(1,nelch) = ipt1.num(1,iel)
        vpocha(nelch,1) = dp12s / sr
*D        write(ioimp,*) ' jeu type 2 ',nelri2,nelch,vpocha(nelch,1),
*D     &                 ipt7.num(1,nelch),ip1,ip2,jp
*
 520  continue
*
*D      write (ioimp,*) ' nb relation type 2 ',nelri2
*
* Ajustement au plus juste puis desactivation des segments lies
* a la rigidite du type 2
      if (nelri2.ne.nelrig) then
        nelrig = nelri2
        segadj,xmatri
        nbelem = nelri2
        nbnn = 4
        segadj,meleme
      endif
      segdes,meleme,xmatri
*
*  si il n'y a rien on dit que pas unilateral pour pas passer dans unilater
*     if (nbelem.eq.0) irigel(6,nrigel) = 0
*
*  Relation du type 1 : noeud-noeud
*-----------------------------------
      if (nelpath.eq.0) goto 700
*  Il y en a au maximum nelpath !
*
*  Ajustement du chpoint et du maillage support associe
*
      nbnn   = 1
      nbelem = nelch + nelpath
      nbsous = 0
      nbref  = 0
      segadj,ipt7
      n = nelch + nelpath
      segadj,mpoval
*
*   creation de la raideur de conditions
*
      nrigel = irigel(/2)+1
      segadj mrigid
      coerig(nrigel)=1.d0
*
*   creation du meleme associe a la relation
*
      nbsous = 0
      nbref  = 0
      nbnn   = 3
      nbelem = nelpath
      segini meleme
      itypel=22
      irigel(1,nrigel)=meleme
*
*     creation des descriptif commun a toutes les raideurs
*
      nligrp=5
      nligrd=nligrp
      segini descr
      lisinc(1)='LX  '
      lisdua(1)='FLX '
      noelep(1)=1
      noeled(1)=1
      do 610 i=2, nligrp,2
        lisinc(i  )=modepl(imo)
        lisinc(i+1)=modepl(imo+1)
        lisdua(i  )=moforc(imo)
        lisdua(i+1)=moforc(imo+1)
        noelep(i  )=(i+2)/2
        noelep(i+1)=noelep(i)
        noeled(i  )=noelep(i)
        noeled(i+1)=noelep(i)
 610  continue
      segdes,descr
      irigel(3,nrigel)=descr
*
*     creation du segment imatri
*
      nelrig = nelpath
      segini xmatri
      irigel(4,nrigel)=xmatri
*
*     ce qu'on cree est unilateral
*
      irigel(6,nrigel)=1
*
*     ce qu'on cree est symetrique
*
      irigel(7,nrigel)=0
*
      nelri1 = 0
*
      do 620 iel=1,nelpath
        ip1 = lpoin(iel,1)
        if (ip1.eq.0) goto 620
*D       write (ioimp,*) 'point en relation point point ',ip1
        lpoin(iel,1)=0
        xp1=xcoor((ip1-1)*idimp1+1)
        yp1=xcoor((ip1-1)*idimp1+2)
        jp =lpoin(iel,2)
        xp =xcoor((jp -1)*idimp1+1)
        yp =xcoor((jp -1)*idimp1+2)
        il1=lpoin(iel,3)
        xn = xpoin(iel,1)
        yn = xpoin(iel,2)
        sn1 = sqrt(xn*xn+yn*yn)
*        write (6,*) ' xn yn avant 630 ',xn,yn
        do 630 jel=iel+1,nelpath
          if (lpoin(jel,1).ne.ip1.or.lpoin(jel,2).ne.jp) goto 630
*          write (ioimp,*) 'point en double '
          lpoin(jel,1)=0
          xn=xn+xpoin(jel,1)
          yn=yn+xpoin(jel,2)
 630    continue
        do 640 jel=iel+1,nelpath
          if (lpoin(jel,2).ne.ip1.or.lpoin(jel,1).ne.jp) goto 640
*          write (6,*) 'point en double '
          lpoin(jel,1)=0
          xn=xn-xpoin(jel,1)
          yn=yn-xpoin(jel,2)
 640    continue
*D       write(ioimp,*) ' xn yn sn ',xn,yn,sn
        sn=sqrt(xn*xn+yn*yn)
        if (sn.lt.sn1*1.e-3) then
          call erreur(942)
          return
        endif
        xn=xn/sn
        yn=yn/sn
*D       write (ioimp,*) ' xn yn sn apres 640 ',xn,yn,sn,ip1,jp
*
        nelri1 = nelri1 + 1
        nelch = nelch + 1
*
*   on range dans le meleme
*
        num(1,nelri1)=il1
        num(2,nelri1)=ip1
        num(3,nelri1)=jp
*
*   on remplit le xmatri
*
*  lambda
        re(1,1,nelri1)= 0.
*  ip1
        re(2,1,nelri1)=  xn
        re(3,1,nelri1)=  yn
*  jp
        re(4,1,nelri1)= -xn
        re(5,1,nelri1)= -yn
*  on transpose
        re(1,2,nelri1) = re(2,1,nelri1)
        re(1,3,nelri1) = re(3,1,nelri1)
        re(1,4,nelri1) = re(4,1,nelri1)
        re(1,5,nelri1) = re(5,1,nelri1)
*  le reste est nul
*
*  remplissage du champoint de depi
*
        ipt7.num(1,nelch) = il1
        scal=(xp-xp1)*xn+(yp-yp1)*yn
        vpocha(nelch,1) = scal
*        write(ioimp,*) ' jeu type 1 ',scal,ip1,jp,il1,il2
*
 620  continue
*
*      write (ioimp,*) ' nb relation type 1 ',nelri1
*
* Ajustement au plus juste puis desactivation des segments lies
* a la rigidite du type 1
      if (nelri1.ne.nelrig) then
        nelrig = nelri1
        segadj,xmatri
        nbelem = nelri1
        nbnn = 3
        segadj,meleme
      endif
      segdes,meleme,xmatri
*
*     if (nelri1.eq.0) irigel(6,nrigel)=0

 700  CONTINUE
*
      segsup,mfopa1
*      GOTO 1000
*
*----------------------------
*    on renvoie le resultat
*----------------------------
 1000 CONTINUE
*
* Ajustement au plus juste du chpoint de depi (jeu) : mpoval et ipt7
* puis desactivation du chpoint
      if (n.ne.nelch) then
        n = nelch
        segadj,mpoval
        nbnn   = 1
        nbelem = nelch
        nbsous = 0
        nbref  = 0
        segadj,ipt7
      endif
      segdes,mpoval,ipt7
      segdes,msoupo,mchpoi
* Desctivation de la matrice de raideur de contact
      segdes,mrigid
*
* Ecriture des objets resultats
      call ecrobj('RIGIDITE',mrigid)
      call ecrobj('CHPOINT',mchpoi)

 900  continue
      segdes,ipt1

      return
      end