Télécharger la plaquette Cast3M


fr en

Télécharger mocon1.eso

Retour à la liste

Numérotation des lignes : 

mocon1

  1. C MOCON1    SOURCE    PV        21/07/08    21:15:01     11068          
  2. C  Mise en forme du maillage de contact:
  3. C conversion du deuxieme maillage en points + mult lagrange contact + mult lagrange frottement
  4.  
  5.       SUBROUTINE MOCON1(MELEME,laprop,itct)
  6. *
  7.  
  8.       IMPLICIT INTEGER(I-N)
  9.       IMPLICIT REAL*8 (A-H,O-Z)
  10.  
  11.  
  12. -INC PPARAM
  13. -INC CCOPTIO
  14. -INC SMCOORD
  15. -INC SMELEME
  16. -INC CCGEOME
  17.       segment icpr(nbpts)                 
  18.  
  19. *  itct = 0: frocable on fait 2 mults    par noeud cable
  20.  
  21.  
  22.  
  23.       segact meleme
  24.       if (lisous(/1).ne.0) call erreur(16)
  25.       if (ierr.ne.0) return
  26. **    write(6,*) ' laprop dans mocon1 ',laprop
  27.       nbsous=0
  28.       nbref=0
  29.       nbelem=num(/2)
  30.       nbnn=num(/1)
  31. *  nombre de lx par noeud
  32.        nbl=1
  33.        if(laprop.eq.3) then
  34.          if(idim.eq.2) nbl=2
  35.          if(idim.eq.3) nbl=3
  36.        endif
  37. *
  38. *  formulation forte
  39.       if (itct.eq.2) goto 1000
  40. *  remplissage icpr pour avoir le nombre de noeuds
  41.        segini icpr
  42.        icp=0
  43.        do j=1,nbelem
  44.        do i=1,nbnn
  45.          ip=num(i,j)
  46.          if(icpr(ip).eq.0) then
  47.            icp=icp+1
  48.            icpr(ip)=icp
  49.          endif
  50.        enddo
  51.        enddo
  52.        nbnn=nbl+1
  53.        nbelem=icp
  54.        nbsous=0
  55.        nbref=0
  56. *  on transforme le maillage en point. Dans le maillage cree,
  57. *  en 1 le multiplicateur de lagrange associe au contact
  58. *  en 2 le (ou les dans le cas faible) points physique
  59. *  en nbnn-1 et nbnn les multiplicateur associe au frottement
  60.        segini ipt1
  61.        ipt1.itypel=22
  62.        do i=1,nbpts
  63.        if(icpr(i).ne.0) then
  64.         ip=icpr(i)
  65.         ipt1.num(2,ip)=i
  66.         nbpts=nbpts+1
  67.         ipt1.num(1,ip)=nbpts
  68.         if(nbl.ge.2) then
  69.            nbpts=nbpts+1
  70.            ipt1.num(3,ip)=nbpts
  71.         endif
  72.         if(nbl.ge.3) then
  73.            nbpts=nbpts+1
  74.            ipt1.num(4,ip)=nbpts
  75.         endif
  76.        endif
  77.        enddo
  78. *  dedoublement (ou trimplement)  des elements dans le cas frocable
  79.        if (itct.eq.0) then
  80.          nbelem=nbelem*idim
  81. **         write(6,*) ' nouveau nbelem ',nbelem
  82.          segadj ipt1
  83.          do j=nbelem/idim,1,-1
  84.           ipt1.num(1,(j-1)*idim+1)=ipt1.num(1,j)
  85.           ipt1.num(2,(j-1)*idim+1)=ipt1.num(2,j)
  86.           if(nbl.ge.2) ipt1.num(3,(j-1)*idim+1)=ipt1.num(3,j)
  87.           if(nbl.ge.3) ipt1.num(4,(j-1)*idim+1)=ipt1.num(4,j)
  88.           nbpts=nbpts+1
  89.           ipt1.num(1,(j-1)*idim+2)=nbpts
  90.           ipt1.num(2,(j-1)*idim+2)=ipt1.num(2,j)
  91.           if(nbl.ge.2) then
  92.            nbpts=nbpts+1
  93.            ipt1.num(3,(j-1)*idim+2)=nbpts
  94.           endif
  95.           if(nbl.ge.3) then
  96.            nbpts=nbpts+1
  97.            ipt1.num(4,(j-1)*idim+2)=nbpts
  98.           endif
  99.  
  100.           if(idim.eq.3) then
  101.           nbpts=nbpts+1
  102.           ipt1.num(1,(j-1)*idim+3)=nbpts
  103.           ipt1.num(2,(j-1)*idim+3)=ipt1.num(2,j)
  104.           if(nbl.ge.2) then
  105.            nbpts=nbpts+1
  106.            ipt1.num(3,(j-1)*idim+3)=nbpts
  107.           endif
  108.           if(nbl.ge.3) then
  109.            nbpts=nbpts+1
  110.            ipt1.num(4,(j-1)*idim+3)=nbpts
  111.           endif
  112.           endif
  113.          enddo
  114.        endif
  115.  
  116. *  maintenant remplir les coordonnees
  117.        segadj mcoord
  118.        il1=0
  119.        il2=0
  120.        do j=1,ipt1.num(/2)
  121.         ip=ipt1.num(2,j)
  122.         il=ipt1.num(1,j)
  123.         if(nbl.ge.2) il1=ipt1.num(3,j)
  124.         if(nbl.ge.3) il2=ipt1.num(4,j)
  125.         do id=1,idim+1
  126.          xc=xcoor((ip-1)*(idim+1)+id)
  127.          xcoor((il-1)*(idim+1)+id)=xc
  128.          if(il1.ne.0)
  129.      >       xcoor((il1-1)*(idim+1)+id)=xc
  130.          if(il2.ne.0)
  131.      >       xcoor((il2-1)*(idim+1)+id)=xc
  132.         enddo
  133.        enddo
  134.        segsup icpr
  135.        segdes meleme
  136.        meleme=ipt1
  137.        return
  138. *  formulation faible
  139.  1000 continue
  140. *
  141. *  on garde le deuxieme maillage en rajoutant les multiplicateurs de lagrange aux elements
  142.        nbnno=nbnn
  143.        nbnn=nbnn+nbl
  144. **     write(6,*) 'nbnno nbnn ',nbnno,nbnn
  145.        segini ipt1
  146.        ipt1.itypel=22
  147.        nbptso=nbpts
  148.        nbpts=nbpts+nbelem*nbl
  149.        segadj mcoord
  150.        do j=1,nbelem
  151.         xc=0.d0
  152.         yc=0.d0
  153.         zc=0.d0
  154.         dc=0.d0
  155.         do i=1,nbnno
  156.          ipt1.num(i+1,j)=num(i,j)
  157.          ip=num(i,j)
  158.          xc=xc+xcoor((idim+1)*(ip-1)+1)
  159.          yc=yc+xcoor((idim+1)*(ip-1)+2)
  160.          zc=zc+xcoor((idim+1)*(ip-1)+3)
  161.          dc=dc+xcoor((idim+1)*(ip-1)+4)
  162.         enddo
  163.         xc=xc/nbnno
  164.         yc=yc/nbnno
  165.         zc=zc/nbnno
  166.         dc=dc/nbnno
  167.         nbptso=nbptso+1
  168.         ipt1.num(1,j)=nbptso
  169.          xcoor((idim+1)*(nbptso-1)+1)=xc
  170.          xcoor((idim+1)*(nbptso-1)+2)=yc
  171.          xcoor((idim+1)*(nbptso-1)+3)=zc
  172.          if(idim.eq.3) xcoor((idim+1)*(nbptso-1)+4)=dc
  173.         if(nbl.ge.2) then
  174.         nbptso=nbptso+1
  175.         ipt1.num(nbnno+2,j)=nbptso
  176.          xcoor((idim+1)*(nbptso-1)+1)=xc
  177.          xcoor((idim+1)*(nbptso-1)+2)=yc
  178.          xcoor((idim+1)*(nbptso-1)+3)=zc
  179.          if(idim.eq.3) xcoor((idim+1)*(nbptso-1)+4)=dc
  180.         endif
  181.         if(nbl.ge.3) then
  182.         nbptso=nbptso+1
  183.         ipt1.num(nbnno+3,j)=nbptso
  184.          xcoor((idim+1)*(nbptso-1)+1)=xc
  185.          xcoor((idim+1)*(nbptso-1)+2)=yc
  186.          xcoor((idim+1)*(nbptso-1)+3)=zc
  187.          if(idim.eq.3) xcoor((idim+1)*(nbptso-1)+4)=dc
  188.         endif
  189.        enddo
  190.        segdes meleme
  191.        meleme=ipt1
  192.        return
  193.  
  194.  
  195.  
  196.  
  197.       end
  198.  
  199.  
  200.  
  201.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales