Télécharger la plaquette Cast3M


fr en

Télécharger frig3c.eso

Retour à la liste

Numérotation des lignes : 

frig3c

  1. C FRIG3C    SOURCE    CB215821  23/01/25    21:15:15     11573          
  2.       SUBROUTINE FRIG3C (maifro,IPRIGI,IPCHJE,IPRIG2)
  3.  
  4.       IMPLICIT INTEGER(I-N)
  5.       IMPLICIT REAL*8(A-H,O-Z)
  6.  
  7. *  Ce sous-programme calcule la raideur de frottement en 3D.
  8. *  il a besoin pour cela du maillage de frottement et de la raideur
  9. *  de contact (ou la raideur totale si c'est plus simple)
  10.  
  11. -INC PPARAM
  12. -INC CCOPTIO
  13. -INC CCREEL
  14. -INC CCGEOME
  15. -INC SMCHPOI
  16. -INC SMELEME
  17. -INC SMRIGID
  18. -INC SMCOORD
  19.  
  20. *  icpr lx du contact   ==> lx du frottement
  21.       segment icpr1(nbpts)
  22.       segment icpr2(nbpts)
  23. *  xjeu champs de jeux initiaux
  24.       segment xjeu(nbpts)
  25. *
  26. *
  27. *  creation et remplissage de icpr
  28. *
  29.       segini icpr1,icpr2
  30.       nbp=0
  31.       meleme=maifro
  32.       segact meleme
  33.       ipt1=meleme
  34.       do is=1,max(1,lisous(/1))
  35.        if (lisous(/1).ne.0) ipt1=lisous(is)
  36.        segact ipt1
  37.        if (ipt1.itypel.ne.22) then
  38.          write (6,*) ' ipt1.itypel ',ipt1.itypel
  39.          call erreur(16)
  40.        endif
  41.        if (ierr.ne.0) return
  42.        do iel=1,ipt1.num(/2)
  43.          il=ipt1.num(1,iel)
  44.          if (icpr1(il).eq.0) then
  45.           nbp=nbp+1
  46.           icpr1(il)=ipt1.num(ipt1.num(/1)-1,iel)
  47.           icpr2(il)=ipt1.num(ipt1.num(/1),iel)
  48.          endif
  49.          if(icpr1(il).ne.ipt1.num(ipt1.num(/1)-1,iel)) call erreur(5)
  50.        enddo
  51.       enddo
  52.  
  53. * remplissage du champ de jeux (demi-frottement si jeu non nul)
  54.  
  55.       segini xjeu
  56.       mchpoi = IPCHJE
  57.       segact mchpoi
  58.       iOK=0
  59.       do 15 isoupo = 1, ipchp(/1)
  60.         msoupo = ipchp(isoupo)
  61.         segact msoupo
  62.         DO 16 i=1,nocomp(/2)
  63.           IF (NOCOMP(i).NE.'FLX ') GOTO 16
  64.           mpoval=ipoval
  65.           segact mpoval
  66.           ipt8=igeoc
  67.           segact ipt8
  68.           DO 17 j=1,vpocha(/1)
  69.             xjeu(ipt8.num(1,j))=vpocha(j,i)
  70.   17      CONTINUE
  71.           iOK=1
  72.   16    CONTINUE
  73.   15  continue
  74.       IF (iOK.NE.1) THEN 
  75.         MOTERR(1:4)='FLX '
  76.         MOTERR(5:8)='DEPI'
  77.         CALL ERREUR(77) 
  78.       ENDIF
  79.       IF (ierr.ne.0) return
  80. *
  81. *  boucle sur les raideurs de contact pour les transformer en frottement
  82. *
  83.       mrigid=iprigi
  84.       segact mrigid
  85.       segini,ri1=mrigid
  86.          er1 = 1.d0
  87.          er2 = 2.71828182845904523536
  88.          er3 = xpi
  89.          sr=sqrt(er1**2+er2**2+er3**2)
  90.          er1=er1/sr
  91.          er2=er2/sr
  92.          er3=er3/sr
  93.       do 10 ir=1,irigel(/2)
  94.         ri1.irigel(1,ir)=0
  95.         ri1.irigel(4,ir)=0
  96.         meleme=irigel(1,ir)
  97.         segact meleme
  98.         ipt1=0
  99.         if (itypel.ne.22) goto 10
  100.         if (irigel(6,ir).eq.0)  goto 10
  101.         nbsous=0
  102.         nbref=0
  103.         nbnn=num(/1)
  104.         nbelem=num(/2)*2
  105.         segini,ipt1
  106.         ipt1.itypel=22
  107.         xmatri=irigel(4,ir)
  108.         segact xmatri
  109.         nligrd=re(/1)
  110.         nligrp=re(/2)
  111.         nelrig=re(/3)*2
  112.         segini,xmatr1
  113.         do iel=1,num(/2)
  114.          il=num(1,iel)
  115. *  coefficient multiplicateur suivant le jeu par rapport a la taille de l'element
  116. *  taille de l'element au carre
  117.          iel1=2*iel-1
  118.          if1=icpr1(il)
  119.          ipt1.num(1,iel1)=if1
  120.          ipt1.icolor(iel1)=icolor(iel)
  121.          do in=2,ipt1.num(/1)
  122.           ipt1.num(in,iel1)=num(in,iel)
  123.          enddo
  124.          sre=sqrt(re(1,nligrp-2,iel)**2+re(1,nligrp-1,iel)**2+
  125.      >            re(1,nligrp,iel)**2)
  126.          do ic=2,re(/2),3
  127.           srev=sqrt(re(1,ic,iel)**2+re(1,ic+1,iel)**2+re(1,ic+2,iel)**2)
  128.           xmatr1.re(1,ic,iel1)=(er2*re(1,ic+2,iel)-er3*re(1,ic+1,iel))
  129.           xmatr1.re(1,ic+1,iel1)=(er3*re(1,ic,iel)-er1*re(1,ic+2,iel))
  130.           xmatr1.re(1,ic+2,iel1)=(er1*re(1,ic+1,iel)-er2*re(1,ic,iel))
  131.           srep=sqrt(xmatr1.re(1,ic,iel1)**2+xmatr1.re(1,ic+1,iel1)**2+
  132.      >              xmatr1.re(1,ic+2,iel1)**2)
  133.           if (srep.ne.srep) write(6,*) ' frig3c ',re(1,ic,iel),
  134.      >         re(1,ic+1,iel),re(1,ic+2,iel)
  135. *  iel1 orthogonal a iel et a 1 e pi
  136.           if (srep/srev.gt.1d-3) then
  137.            xmatr1.re(1,ic,iel1)=xmatr1.re(1,ic,iel1)*srev/srep
  138.            xmatr1.re(1,ic+1,iel1)=xmatr1.re(1,ic+1,iel1)*srev/srep
  139.            xmatr1.re(1,ic+2,iel1)=xmatr1.re(1,ic+2,iel1)*srev/srep
  140.           else
  141.            write(6,*) ' frig3c second choix ',srep,srev
  142.           srev=sqrt(re(1,ic,iel)**2+re(1,ic+1,iel)**2+re(1,ic+2,iel)**2)
  143.            srep=sqrt(re(1,ic+2,iel)**2+re(1,ic,iel)**2)
  144.              if (srep.lt.xpetit) srep=1.d0
  145.              xmatr1.re(1,ic,iel1)=-re(1,ic+2,iel)*srev/srep
  146.              xmatr1.re(1,ic+1,iel1)=0.d0
  147.              xmatr1.re(1,ic+2,iel1)=re(1,ic,iel)*srev/srep
  148.           endif
  149. *        write(6,*) ' re ',re(1,ic,iel),re(1,ic+1,iel),re(1,ic+2,iel)
  150. *        write(6,*)        xmatr1.re(1,ic,iel1),xmatr1.re(1,ic+1,iel1),
  151. *    >          xmatr1.re(1,ic+2,iel1)
  152.          enddo
  153. **             xmatr1.re(1,ic,iel1)=-re(1,ic+1,iel)*srev/srep
  154. **             xmatr1.re(1,ic+1,iel1)=re(1,ic,iel)*srev/srep
  155. **             xmatr1.re(1,ic+2,iel1)=0.d0
  156. **          enddo
  157.          do ic=2,re(/1)
  158.           xmatr1.re(ic,1,iel1)=xmatr1.re(1,ic,iel1)
  159.          enddo
  160.          iel2=2*iel
  161.          if2=icpr2(il)
  162.          ipt1.num(1,iel2)=if2
  163.          ipt1.icolor(iel2)=icolor(iel)
  164.          do in=2,ipt1.num(/1)
  165.           ipt1.num(in,iel2)=num(in,iel)
  166.          enddo
  167.          sre=sqrt(re(1,nligrp-2,iel)**2+re(1,nligrp-1,iel)**2+
  168.      >            re(1,nligrp,iel)**2)
  169.          do ic=2,re(/2),3
  170.           if (sre.lt.xpetit) sre=1.d0
  171. *  iel2 orthogonal a iel et iel1
  172.            xmatr1.re(1,ic,iel2)=
  173.      >                 (re(1,nligrp-1,iel)*xmatr1.re(1,ic+2,iel1)-
  174.      >                  re(1,nligrp  ,iel)*xmatr1.re(1,ic+1,iel1))/sre
  175.            xmatr1.re(1,ic+1,iel2)=
  176.      >                 (re(1,nligrp  ,iel)*xmatr1.re(1,ic  ,iel1)-
  177.      >                  re(1,nligrp-2,iel)*xmatr1.re(1,ic+2,iel1))/sre
  178.            xmatr1.re(1,ic+2,iel2)=
  179.      >                 (re(1,nligrp-2,iel)*xmatr1.re(1,ic+1,iel1)-
  180.      >                  re(1,nligrp-1,iel)*xmatr1.re(1,ic,  iel1))/sre
  181.          enddo
  182.          do ic=2,re(/1)
  183.           xmatr1.re(ic,1,iel2)=xmatr1.re(1,ic,iel2)
  184.          enddo
  185.         enddo
  186.         segdes xmatri
  187.         ri1.irigel(1,ir)=ipt1
  188.         ri1.irigel(4,ir)=xmatr1
  189.         ri1.irigel(6,ir)=2
  190.  
  191.   10  continue
  192.       segdes mrigid
  193. *
  194. *  boucle de compaction du resultat
  195. *
  196.       mrigid=ri1
  197.       irr=0
  198.       do 100 ir=1,irigel(/2)
  199.        meleme=irigel(1,ir)
  200.        xmatri=irigel(4,ir)
  201.        if (meleme.eq.0) goto 100
  202.        ill=0
  203.        do iel=1,num(/2)
  204.         if (num(1,iel).ne.0) then
  205.          ill=ill+1
  206.          if (ill.ne.0) then
  207.           do in=1,num(/1)
  208.            num(in,ill)=num(in,iel)
  209.           enddo
  210.           icolor(ill)=icolor(iel)
  211.           do ic=1,re(/1)
  212.            re(1,ic,ill)=re(1,ic,iel)
  213.            re(ic,1,ill)=re(ic,1,iel)
  214.           enddo
  215.          endif
  216.         endif
  217.        enddo
  218.        if (ill.eq.0) goto 100
  219.        if (ill.ne.num(/2)) then
  220.         nbsous=0
  221.         nbref=0
  222.         nbnn=num(/1)
  223.         nbelem=ill
  224.         segadj meleme
  225.        endif
  226. **     write (6,*) ' meleme sortie dans frig2c '
  227. **     call ecmail(meleme,0)
  228.  
  229.  
  230.        irr=irr+1
  231.        if (irr.ne.ir) then
  232.         do ir1=1,irigel(/1)
  233.          irigel(ir1,irr)=irigel(ir1,ir)
  234.         enddo
  235.         coerig(irr)=coerig(ir)
  236.        endif
  237.  100  continue
  238.       nrigel=irr
  239.       if (irigel(/2).ne.irr) segadj mrigid
  240.       iprig2=mrigid
  241. **    call prrigi(mrigid,1)
  242.       segsup icpr1,icpr2,xjeu
  243.       return
  244.       end
  245.  
  246.  
  247.  
  248.  
  249.  
  250.  
  251.  
  252.  
  253.  
  254.  
  255.  
  256.  
  257.  
  258.  
  259.  
  260.  
  261.  
  262.  
  263.  

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