Télécharger la plaquette Cast3M


fr en

Télécharger numopt.eso

Retour à la liste

Numérotation des lignes : 

numopt

  1. C NUMOPT    SOURCE    OF166741  24/10/21    21:15:20     12042          
  2. C   RACINE DE LA NUMEROTATION      POUR LA SORTIE SUR FAC
  3. C
  4.       SUBROUTINE NUMOPT(MELEME,ICPR,NODES)
  5.       IMPLICIT INTEGER(I-N)
  6.  
  7. -INC PPARAM
  8. -INC CCOPTIO
  9. -INC CCPRECO
  10.  
  11. -INC SMELEME
  12. -INC SMCOORD
  13.  
  14.       CHARACTER*4 MVAL
  15.       DATA MVAL/'NOOP'/
  16.       SEGMENT MEMJT1(NKON)
  17.       SEGMENT MEMJT2(NKON)
  18.       SEGMENT IPOME(NODES+1)
  19.       SEGMENT JNT(NODES)
  20.       SEGMENT JMEM(NODES)
  21.       SEGMENT ICPR(nbpts)
  22.       SEGMENT ICPREN(nbpts)
  23.       SEGMENT IDCP(nbpts)
  24.       SEGMENT JOINT(NODES)
  25.       SEGMENT NEWJT(NODES)
  26.  
  27.       icpren=0
  28.       segact meleme
  29. *  cas trivial ou le maillage est fait de points isoles
  30.       CALL LIRMOT(MVAL,1,IRET,0)
  31.        itypto=0
  32.        ipt1=meleme
  33.        do k=1,max(1,lisous(/1))
  34.         if (lisous(/1).ne.0) ipt1=lisous(k)
  35.         segact ipt1
  36.         itypto = itypto+ipt1.itypel-1
  37.        enddo
  38.  
  39.       if (itypto.eq.0.or.iret.eq.1.or.nucrou.eq.-1) then
  40.        segact icpr*MOD
  41.        ia=0
  42.        ipt1=meleme
  43.        do k=1,max(1,lisous(/1))
  44.         if (lisous(/1).ne.0) ipt1=lisous(k)
  45.         segact ipt1
  46.        do j=  ipt1.num(/2),1,-1
  47.        do i=1,ipt1.num(/1)
  48.         if (icpr(ipt1.num(i,j)).eq.0) then
  49.          ia=ia+1
  50.          icpr(ipt1.num(i,j))=ia
  51.         endif
  52.        enddo
  53.        enddo
  54.        enddo
  55.        nodes=ia
  56.        return
  57.       endif
  58.  
  59.       nucro=0
  60. *  si plus de 4000 noeuds dans le meleme numérotation NS
  61.       segact icpr*mod
  62.       do i=1,icpr(/1)
  63.        icpr(i)=0
  64.       enddo
  65.        ia=0
  66.        ib=0
  67.        ipd=0
  68.        ipt1=meleme
  69.        ifrot=0
  70.        do 60 isous=1,max(1,lisous(/1))
  71.        if (lisous(/1).ne.0) ipt1=lisous(isous)
  72.        segact ipt1
  73.        if (ipt1.itypel.eq.-49) then
  74. *         write (6,*) ' frottement '
  75.           ifrot=1
  76.        endif
  77.        NBNN1=ipt1.num(/1)
  78.        NBEL1=ipt1.num(/2)
  79.        do  i=1,NBEL1
  80.         do  j=1,NBNN1
  81.          iel1=ipt1.num(j,i)
  82.          if (icpr(iel1).eq.0) then
  83.           ia=ia+1
  84.           icpr(iel1)=ia
  85.          endif
  86.         enddo
  87.        enddo
  88.  
  89.        if (NBNN1.gt.ib .AND. NBEL1.GT.0.and.abs(ipt1.itypel).ne.49) then
  90.          ib=NBNN1
  91.          ipd=icpr(ipt1.num(1,1))
  92.        endif
  93.   60  continue
  94.       if (ia.lt.ib*2 .and.ia.gt.100) then
  95. **     write (6,*) ' superelement detecte ',ia,ib
  96.        ia=1
  97.       else
  98.        ipd=0
  99.       endif
  100.       do i=1,icpr(/1)
  101.        icpr(i)=0
  102.       enddo
  103. *  choix de la numerotation en fonction du nombre de noeuds
  104. **    write(6,*) ' numopt ia ',ia
  105.       if (ia.gt.4000) nucro=2
  106.       if (ia.le.4000) nucro=0
  107. *  Reverse CM en iteratif
  108.       IF(NUCROU.EQ.1) nucro=0
  109.       IF(NUCROU.EQ.2) nucro=2
  110. *  en cas de frottement forcément nested dissection pour placer correctement les mult de frot
  111.       if (ifrot.eq.1) nucro=2
  112. *  ajuster le nombre de zero consecutif
  113.       if (nucro.eq.0) izrosf=16
  114.       if (nucro.eq.2) izrosf=48
  115.       if(ipd.ne.0) izrosf=48
  116.       CALL LIRMOT(MVAL,1,IRET,0)
  117.       if (iret.eq.0) then
  118. *  a t'on deja fait l'operation ???
  119.       segini icpren,idcp
  120.       incren=0
  121.       do isous=2,lisous(/1)
  122.        ipt3=lisous(isous)
  123.        segact ipt3
  124. * SG 2011/10/07 on a échangé l'ordre des boucles pour le remettre dans
  125. *               l'ordre "usuel"
  126.        do j=1,ipt3.num(/2)
  127.        do i=1,ipt3.num(/1)
  128.         ip=ipt3.num(i,j)
  129.         if (icpren(ip).eq.0) then
  130.           incren=incren+1
  131.           icpren(ip)=incren
  132.         endif
  133.        enddo
  134.        enddo
  135.       enddo
  136.       do i=1,icpren(/1)
  137.        if (icpren(i).ne.0) idcp(icpren(i))=i
  138.       enddo
  139.  
  140.       do 150 ip=1,npreco
  141.       if (prenum(ip).ne.0) then
  142.        ipt1=prenum(ip)
  143.        segact ipt1
  144. *      write (6,*) ' prenum ',ip,lisous(/1),ipt1.lisous(/1)
  145.        if (lisous(/1)+1.eq.ipt1.lisous(/1)) then
  146. *SG Attention ! On a osé stocker le type de renumérotation qui avait été
  147. *SG   calculée dans le itypel
  148.         ipt2=ipt1.lisous(ipt1.lisous(/1))
  149.         segact ipt2
  150.         nucroa=ipt2.itypel
  151.         segdes ipt2
  152.         if (nucro.ne.nucroa) goto 11
  153. *  verif maillage identique
  154. *      write (6,*) ' maillage 1 ',lisous(/1)
  155. *      write (6,*) ' maillage 2 ',ipt1.lisous(/1)
  156.         do isous=2,lisous(/1)
  157.         ipt3=lisous(isous)
  158.         ipt4=ipt1.lisous(isous)
  159.         segact ipt3,ipt4
  160. *       write (6,*) ipt3.num(/1),ipt4.num(/1),ipt3.num(/2),ipt4.num(/2)
  161. *       write (6,*) ' itypel ',ipt3.itypel,ipt4.itypel
  162.         if (ipt3.itypel.ne.ipt4.itypel) goto 10
  163.         if (ipt3.num(/1).ne.ipt4.num(/1)) goto 10
  164.         if (ipt3.num(/2).ne.ipt4.num(/2)) goto 10
  165. * SG 2011/10/07 on a échangé l'ordre des boucles pour le remettre dans
  166. *               l'ordre "usuel"
  167.         do j=1,ipt3.num(/2)
  168.         do i=1,ipt3.num(/1)
  169. *        write (6,*) i,j,ipt3.num(i,j),
  170. *    >     icpren(ipt3.num(i,j)),ipt4.num(i,j)
  171.          if (icpren(ipt3.num(i,j)).ne.ipt4.num(i,j)) goto 10
  172.         enddo
  173.         enddo
  174.         segdes ipt3,ipt4
  175.         enddo
  176.         if (iimpi.eq.2022) write (6,*)
  177.      >     ' preconditionnement de la numerotation ',ip,
  178.      >     (prenum(i),i=1,30)
  179.         ipt2=ipt1.lisous(ipt1.lisous(/1))
  180.         segact ipt2
  181.         segact icpr*MOD
  182.         nodes=ipt2.num(/2)
  183. *       write (6,*) ' icpr(/1),nodes ',icpr(/1),nodes
  184.         do i=1,nodes
  185.            if (idcp(ipt2.num(1,i)).ne.0) icpr(idcp(ipt2.num(1,i)))=i
  186.         enddo
  187.         ia=nodes
  188.         segdes ipt1,ipt2
  189. *     write (6,*) ' apres precond dans numopt '
  190. *     write (6,*) (icpr(i),i=1,icpr(/1))
  191.         return
  192.  
  193.   10   continue
  194.        segdes ipt3,ipt4
  195.   11   continue
  196.        endif
  197.       endif
  198.  150  continue
  199.       if (nucro.eq.1) then
  200.        call numop1(MELEME,ICPR,NODES)
  201.        goto 300
  202.       endif
  203.       if (nucro.eq.2) then
  204.        call numop2(MELEME,ICPR,NODES)
  205.        segact icpr*mod
  206.        goto 300
  207.       endif
  208.       endif
  209.       if (iret.eq.1) call refus
  210.       IENORM=2000000000
  211.       NODES=nbpts
  212.       SEGACT ICPR*MOD
  213.       SEGACT MELEME
  214.       DO  I=1,ICPR(/1)
  215.       ICPR(I)=0
  216.       ENDDO
  217.       IPT1=MELEME
  218.       IKOU=0
  219.       DO 202 IO=1,MAX(1,LISOUS(/1))
  220.       IF (LISOUS(/1).NE.0) THEN
  221.        IPT1=LISOUS(IO)
  222.        SEGACT IPT1
  223.       ENDIF
  224.       DO 203 J=1,IPT1.NUM(/2)
  225. *  inversion car on reinerse a la fin
  226.       DO 204 I=1,IPT1.NUM(/1)
  227.       IJ=IPT1.NUM(I,J)
  228.       IF (ICPR(IJ).NE.0) GOTO 204
  229.        IKOU=IKOU+1
  230.       ICPR(IJ)=IKOU
  231.  204  CONTINUE
  232.  203  CONTINUE
  233.  202  CONTINUE
  234.       NODES=IKOU
  235.       SEGINI JNT,JMEM,NEWJT,JOINT,IPOME
  236.       DO I=1,NODES
  237.       JMEM(I)=0
  238.       enddo
  239.       IPT1=MELEME
  240.       DO 3 IO=1,MAX(1,LISOUS(/1))
  241.       IF (LISOUS(/1).NE.0) IPT1=LISOUS(IO)
  242.       DO  J=1,IPT1.NUM(/2)
  243.       DO  I=1,IPT1.NUM(/1)
  244.       JMEM(ICPR(IPT1.NUM(I,J)))=JMEM(ICPR(IPT1.NUM(I,J)))+1
  245.       enddo
  246.       enddo
  247.    3  CONTINUE
  248.       IPOME(1)=0
  249.       DO 6 I=1,NODES
  250.       IPOME(I+1)=IPOME (I) + JMEM(I)
  251.    6  CONTINUE
  252.       NKON=IPOME(NODES+1)
  253.       SEGINI MEMJT1,memjt2
  254.       DO  J=1,NODES
  255.       JMEM(J)=0
  256.       enddo
  257.       IPT1=MELEME
  258.       DO 101 IO=1,MAX(1,LISOUS(/1))
  259.       IF (LISOUS(/1).NE.0) IPT1=LISOUS(IO)
  260.       DO  I=1,IPT1.NUM(/2)
  261.       DO  J=1,IPT1.NUM(/1)
  262.       IND=ICPR(IPT1.NUM(J,I))
  263.       JMEM(IND)=JMEM(IND)+1
  264.       MEMJT1(IPOME(IND)+JMEM(IND))=I
  265.       MEMJT2(IPOME(IND)+JMEM(IND))=IO
  266.       enddo
  267.       enddo
  268.  101  CONTINUE
  269.       IDIFF=1
  270.       CALL OPTNUM(MELEME,MEMJT1,memjt2,JMEM,JNT,NEWJT,JOINT,IDIFF,IPOME,
  271.      # ICPR,NODES,NOOPTI,ipd,nucro)
  272. *  ne pas inverse si superelement pour conserver sa numerotation
  273.       if (ipd.eq.0) then
  274.       DO 110 I=1,NODES
  275.       JNT(I)=NODES-JNT(I)+1
  276.  110  CONTINUE
  277.       endif
  278. C  PERMUTER LES COORDONNEES ET CORRIGER IPAD
  279.       DO 111 I=1,icpr(/1)
  280.       IF (ICPR(I).EQ.0) GOTO 111
  281.       ICPR(I)=JNT(ICPR(I))
  282.  111  CONTINUE
  283.       IF(NOOPTI .EQ.0)  GO TO 116
  284.       IA=0
  285.       DO 115 I=1,icpr(/1)
  286.       IF(ICPR(I).EQ.0) GO TO 115
  287.       IA=IA+1
  288.       ICPR(I)=IA
  289.   115 CONTINUE
  290.   116 CONTINUE
  291.       SEGSUP MEMJT1,MEMJT2,JMEM,NEWJT,JOINT,IPOME
  292.       SEGSUP JNT
  293. * Mise en commentaire du retuen par SG car comme cela on préconditionne
  294. * aussi  le Cuthill-McKee
  295. *      RETURN
  296.  300  continue
  297. *  sauver le resultat au cas ou on veuille recommencer
  298.       segini,ipt1=meleme
  299.       nbsous=ipt1.lisous(/1)+1
  300.       nbref=ipt1.lisref(/1)
  301.       nbnn=ipt1.num(/1)
  302.       nbelem=ipt1.num(/2)
  303.       segadj ipt1
  304.       nbnn=1
  305.       nbsous=0
  306.       nbref=0
  307.       nbelem=nodes
  308.       segini ipt2
  309. *SG Attention ! On ose stocker le type de renumérotation qui a été
  310. *SG   calculée dans le itypel
  311.       ipt2.itypel=nucro
  312.       ia=0
  313.       do 120 i=1,icpr(/1)
  314.        if (icpr(i).gt.nodes.or.icpr(i).eq.0) goto 120
  315.        ia=ia+1
  316.        ipt2.num(1,icpr(i))=icpren(i)
  317.  120  continue
  318. *     call ecmail(ipt2)
  319.       ipt1.lisous(ipt1.lisous(/1))=ipt2
  320.       segdes ipt2
  321.       do 125 is=1,ipt1.lisous(/1)-1
  322.        ipt2=ipt1.lisous(is)
  323.        segini,ipt3=ipt2
  324.        segdes ipt2
  325. * SG 2011/10/07 on a échangé l'ordre des boucles pour le remettre dans
  326. *               l'ordre "usuel"
  327.        do il=1,ipt3.num(/2)
  328.        do ip=1,ipt3.num(/1)
  329.         ipt3.num(ip,il)=icpren(ipt3.num(ip,il))
  330.        enddo
  331.        enddo
  332.        segdes ipt3
  333.        ipt1.lisous(is)=ipt3
  334.  125  continue
  335.       if (ipt1.lisous(/1).eq.1) goto 131
  336.       segdes ipt1
  337.       do 130 ip=npreco-1,1,-1
  338.       prenum(ip+1)=prenum(ip)
  339.  130  continue
  340.       prenum(1)=ipt1
  341.  131  continue
  342.       if (icpren.ne.0) segsup icpren,idcp
  343.       segact icpr*mod
  344. *     write (6,*) ' fin de numopt '
  345. *     write (6,*) (icpr(i),i=1,icpr(/1))
  346. *
  347.       RETURN
  348.       END
  349.  
  350.  
  351.  

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