Télécharger la plaquette Cast3M


fr en

Télécharger noepr2.eso

Retour à la liste

Numérotation des lignes : 

noepr2

  1. C NOEPR2    SOURCE    PV090527  24/01/18    21:15:01     10699          
  2. C****************************************************************************
  3. C****************************************************************************
  4. C*************NOEPERI ..NOEuds PERIpheriques*********************************
  5. C****************************************************************************
  6. C****************************************************************************
  7.  
  8. C  NOEPERI part de PIVOT,lui associe l=1, associe l=l+1 a ses voisins directs,
  9. C  repart des voisins directs pour associer un l a leur voisins....
  10. C  LONG=max(l).
  11. C  NRELONG(I)=l ,NOELON contient les noeuds tels que l=LONG.
  12. *  noepr2 rend en prime la distance a la frontiere, sa taille et le
  13. *  desequilibre des domaines
  14. *  si nbtot n'est pas nul, noepr2 s'arrete des qu'il a trouve la
  15. *  frontiere (il a remplis nbtot/2) pts. nrelong est alors remis a
  16. *  zero
  17. *  Octobre 2014: croissance simultanee de deux zones basees sur les pivots impairs ou pairs
  18. *
  19.         SUBROUTINE NOEPR2(IADJ,JADJC,PIVOT,NRELONG,NOELON,isens,dimlon,
  20.      >  NODES,IPOS,NBTOT,lfront,lfron,londim,fcout,lpiv,iccon,icouch)
  21.       IMPLICIT INTEGER(I-N)
  22.       IMPLICIT REAL*8 (A-H,O-Z)
  23. -INC CCREEL
  24.         integer iadj(*),jadjc(*)
  25.         integer pivot(lpiv)
  26.  
  27.         INTEGER LONG,DIMLON,Y,X,dimini,diminj
  28.         integer nrelong(*),noelon(*)
  29.  
  30.         integer ipos(*)
  31.         integer londim(*)
  32.         INTEGER NODES
  33.  
  34.         INTEGER diml1,DIML2,diml3,diml1i,diml2f,diml1f,diml2i
  35.         LOGICAL bool1,bool2,pretest
  36.  
  37.         pretest=.false.
  38.         MDOMN=-1
  39.         do i=1,lpiv
  40.           if (pivot(i).ne.0) then
  41.            MDOMN=IPOS(PIVOT(i)+NODES)
  42.            goto 1
  43.           endif
  44.         enddo
  45.    1    continue
  46.         mode=1
  47.         if (nbtot.eq.0) mode=2
  48.         if (nbtot.eq.-1) mode=3
  49.         if (nbtot .gt. nodes)  call erreur(5)
  50.         nodz=nodes
  51.         if (nbtot.gt.0) nodz=nbtot
  52.         DIMLON=0
  53.         LONG=3
  54. *       write (6,*) ' noepr2 lpiv pivot ',lpiv,(pivot(i),i=1,lpiv)
  55.         nbz1=0
  56.         nbz2=0
  57. *  insertion pivot impair zone 1
  58.         DIML1=0
  59.         do 10 i=1,lpiv,2
  60. **       write (6,*) ' i pivot(i) ',i,pivot(i)
  61.          if (pivot(i).le.0) goto 10
  62.          if (ipos(pivot(i)+nodes).ne.mdomn) goto 10
  63. *        if (i.eq.1) write (6,*) 'noepr2 mode i pivot nrelong',
  64. *    >    i,pivot(i),nrelong(pivot(i))
  65.          if (nrelong(pivot(i)).ne.0) goto 10
  66.          diml1=diml1+1
  67.          noelon(diml1)=pivot(i)
  68.          nrelong(pivot(i))=1
  69.          nbz1=nbz1+1
  70.   10   continue
  71. *  insertion pivot pair zone 2
  72.         DIML2=nodz+1
  73.         do 11 i=2,lpiv,2
  74.          if (pivot(i).le.0) goto 11
  75. *        write (6,*) 'noepr2 mode i pivot nrelong',
  76. *    >    i,pivot(i),nrelong(pivot(i))
  77.          if (ipos(pivot(i)+nodes).ne.mdomn) goto 11
  78.          if (nrelong(pivot(i)).ne.0) goto 11
  79.          diml2=diml2-1
  80.          noelon(diml2)=pivot(i)
  81.          nrelong(pivot(i))=3
  82.          nbz2=nbz2+1
  83.   11   continue
  84.        if (diml1.eq.0.and.diml2.eq.nodz+1) then
  85.          fcout=xgrand
  86.          lfront=1
  87.          return
  88.        endif
  89. *
  90. *  croissance des deux zones
  91. *
  92.        icr1=2
  93.        icr2=3
  94.        bool1=.true.
  95.        bool2=.true.
  96.        diml1i=1
  97.        diml2f=nodz
  98.        icouch=0
  99.   20   continue
  100. ** evaluation frontiere
  101.        if (bool1.and.(icr1*nbz1.le.icr2*nbz2.or..not.bool2)) then
  102.        if (pretest) then
  103.        do 22 i=diml1i,diml1
  104.             nz1=0
  105.             nz2=0
  106.             x=noelon(i)
  107.             if (nrelong(x).eq.2) goto 22
  108.             DO 29 J=IADJ(X),IADJ(X+1)-1
  109.               Y=JADJC(J)
  110.               nrly=nrelong(y)
  111. **            if (mdomn.ne.ipos(y+nodes))goto 29
  112.               if (nrly.eq.1) nz1=nz1+1
  113.               if (nrly.eq.3) nz2=nz2+1
  114.   29        continue
  115.               if (nz2.gt.nz1) then
  116.                 nbz1=nbz1-1
  117.                 nrelong(x)=2
  118.               endif
  119. *             write (6,*) 'frontiere 1 => 2',nz1,nz2
  120.   22   continue
  121.        endif
  122.         bool1=.false.
  123.         diml1f=diml1
  124. *       write (6,*) ' diml1i diml1f ',diml1i,diml1f
  125.         do 30 i=diml1i,diml1f
  126.             x=noelon(i)
  127.             if(nrelong(x).eq.2) goto 30
  128.             DO 40 J=IADJ(X),IADJ(X+1)-1
  129.                Y=JADJC(J)
  130.                if (mdomn.ne.ipos(y+nodes))goto 40
  131.                  IF(NRELONG(Y).EQ.0) then
  132. *  insertion zone 1
  133.                    diml1=diml1+1
  134.                    if (.not.bool1) icouch=icouch+1
  135.                    noelon(diml1)=y
  136.                    nrelong(y)=1
  137.                    nbz1=nbz1+1
  138.                    BOOL1=.true.
  139. *     if (mode.ne.2) write (6,*) ' insertion dans 1 de ',diml1,y
  140.                    goto 40
  141.                  elseif (NRELONG(Y).eq.1) then
  142. *  zone 1
  143. *     if (mode.ne.2) write (6,*) ' deja en zone 1      ',y
  144.                    goto 40
  145.                  elseif (nrelong(y).eq.2) then
  146. *  frontiere
  147. *     if (mode.ne.2) write (6,*) ' deja en frontiere 1 ',y
  148.                    goto 40
  149.                  elseif (nrelong(y).eq.3) then
  150. *  zone 2 <<> frontiere
  151. *     if (mode.ne.2) write (6,*) ' passage 1 sur frontiere de ',y
  152.                    nrelong(y)=2
  153.                    nbz2=nbz2-1
  154.                    pretest=.true.
  155.                    goto 40
  156.                  endif
  157.   40        continue
  158.   30    continue
  159.         diml1i=diml1f+1
  160.        endif
  161.        if (bool2.and.(icr1*nbz2.le.icr2*nbz1.or..not.bool1)) then
  162.        if (pretest) then
  163.        do 23 i=diml2f,diml2,-1
  164.             nz1=0
  165.             nz2=0
  166.             x=noelon(i)
  167.             if (nrelong(x).eq.2) goto 23
  168.             DO 28 J=IADJ(X+1)-1,IADJ(X),-1
  169.               Y=JADJC(J)
  170.               nrly=nrelong(y)
  171. **            if (mdomn.ne.ipos(y+nodes))goto 28
  172.               if (nrly.eq.1) nz1=nz1+1
  173.               if (nrly.eq.3) nz2=nz2+1
  174.   28        continue
  175.               if (nz1.gt.nz2) then
  176.                 nbz2=nbz2-1
  177.                 nrelong(x)=2
  178.               endif
  179. *             write (6,*) 'frontiere 3 => 2',nz1,nz2
  180.   23   continue
  181.        endif
  182.         bool2=.false.
  183.         diml2i=diml2
  184. *       write (6,*) ' diml2f diml2i ',diml2f,diml2i
  185.         do 50 i=diml2f,diml2i,-1
  186.             x=noelon(i)
  187.             if (nrelong(x).eq.2) goto 50
  188.             DO 60 J=IADJ(X+1)-1,IADJ(X),-1
  189.                Y=JADJC(J)
  190.                if (mdomn.ne.ipos(y+nodes))goto 60
  191.                  nrly=nrelong(y)
  192.                  IF(NRLY.EQ.0) then
  193. *  insertion zone 3
  194.                    diml2=diml2-1
  195.                    if (.not.bool2) icouch=icouch+1
  196.                    noelon(diml2)=y
  197.                    nrelong(y)=3
  198.                    nbz2=nbz2+1
  199.                    BOOL2=.true.
  200. *     if (mode.ne.2) write (6,*) ' insertion dans 2 de ',diml2,y
  201.                    goto 60
  202.                  elseif (NRLY.eq.3) then
  203. *  zone 3
  204. *     if (mode.ne.2) write (6,*) ' deja en zone 2      ',y
  205.                    goto 60
  206.                  elseif (nrly.eq.2) then
  207. *  frontiere
  208. *     if (mode.ne.2) write (6,*) ' deja en frontiere 2 ',y
  209.                    goto 60
  210.                  elseif (nrly.eq.1) then
  211. *  zone 1 ==> frontiere
  212. *     if (mode.ne.2) write (6,*) ' passage 2 sur frontiere de ',y
  213.                    nrelong(y)=2
  214.                    nbz1=nbz1-1
  215.                    pretest=.true.
  216.                    goto 60
  217.                  endif
  218.   60        continue
  219.   50    continue
  220.         diml2f=diml2i-1
  221.        endif
  222.        if (bool1.or.bool2) goto 20
  223.        nbtotn=diml1+nodz+1-diml2
  224. *  dans le cas ou la zone n'est pas connexe on va completer par
  225. *  la partie non connexe
  226.         if (mode.ne.2.and.nbtot.gt.nbtotn) then
  227. *       write (6,*) ' ajout autres composantes connexes ',diml1,diml2,
  228. *    >  nbtotn,nbtot,nodz
  229. *  si pas trop de noeuds
  230. *  on commence par examiner le voisinage de la frontiere car c'est moins cher
  231. ****    if (mode.ne.3) goto 21
  232.         diml1i=diml1
  233.         do 200 i=1,nbtot
  234.          x=noelon(i)
  235.          if (x.eq.0) goto 200
  236.          if (nrelong(x).ne.2) goto 200
  237.             DO 210 J=IADJ(X),IADJ(X+1)-1
  238.                Y=JADJC(J)
  239.                if (mdomn.ne.ipos(y+nodes))goto 210
  240.                  IF(NRELONG(Y).EQ.0) then
  241.                     DIML1=DIML1+1
  242.                     NOELON(DIML1)=Y
  243.                     NRELONG(Y)=1
  244.                     nbtotn=nbtotn+1
  245.                     if (nbtot.le.nbtotn) goto 21
  246.                endif
  247.  210        continue
  248.  200    continue
  249.           if (diml1.eq.diml1i) goto 220
  250. *  Si on a rajoute des noeuds, on fait aussi leurs voisinage
  251.  230    continue
  252.           nodi=diml1i+1
  253.           diml1i=diml1
  254.           nodf=diml1
  255.         do 240 i=nodi,nodf
  256.          x=noelon(i)
  257.             DO 250 J=IADJ(X),IADJ(X+1)-1
  258.                Y=JADJC(J)
  259.                if (mdomn.ne.ipos(y+nodes))goto 250
  260.                  IF(NRELONG(Y).EQ.0) then
  261.                     DIML1=DIML1+1
  262.                     NOELON(DIML1)=Y
  263.                     NRELONG(Y)=1
  264.                     nbtotn=nbtotn+1
  265.                     if (nbtot.le.nbtotn) goto 21
  266.                endif
  267.  250        continue
  268.  240    continue
  269.         if (diml1.ne.diml1i) goto 230
  270. 220     continue
  271.         if (mode.ne.3) goto 21
  272.  
  273. *       write (6,*) ' ajout 2autres composantes connexes ',diml1,diml2,
  274. *    >  nbtotn,nbtot,nodz
  275.   21    continue
  276.         endif
  277. * Si on n'a toujours pas notre compte, on balaye tout
  278.         if ((iccon.eq.1.and.mode.eq.2).or. 
  279.      >       (mode.eq.1.and.nbtot.gt.nbtotn)) then
  280.         do Y=1,NODES
  281.           IF((NRELONG(Y).EQ.0).AND.(MDOMN.EQ.IPOS(Y+NODES))) THEN
  282.                     DIML1=DIML1+1
  283.                     NOELON(DIML1)=Y
  284.                     NRELONG(Y)=1
  285.                     nbtotn=nbtotn+1
  286.                     if (mode.eq.1.and.nbtot.le.nbtotn) goto 215
  287.           ENDIF
  288.         enddo
  289.   215   continue
  290.         endif
  291.        nbtot=nbtotn
  292. *
  293. *  remise au carre des deux zones et de la frontiere
  294. *
  295.        nbtot=diml1+nodz+1-diml2
  296. *      write (6,*) ' nodes diml1 diml2 ',nodes,diml1,diml2,nbtot
  297. *  classer a la fin par permutation dans zone 1 la frontiere
  298.        diml1f=diml1
  299.        do 100 i=1,diml1-1
  300.          x=noelon(i)
  301.            if (nrelong(x).eq.2) then
  302.              do 105 j=diml1f,i+1,-1
  303.                y=noelon(j)
  304.                if (nrelong(y).eq.1) goto 106
  305.  105         continue
  306.              goto 107
  307.  106        continue
  308.             diml1f=j-1
  309.             noelon(i)=y
  310.             noelon(j)=x
  311.            endif
  312.  100   continue
  313.  107   continue
  314. *  classer au debut par permutation dans zone 2 la frontiere
  315.        diml2i=diml2
  316.        do 110 i=nodz,diml2+1,-1
  317.          x=noelon(i)
  318.            if (nrelong(x).eq.2) then
  319.              do 115 j=diml2i,i-1
  320.                y=noelon(j)
  321.                if (nrelong(y).eq.3) goto 116
  322.  115         continue
  323.              goto 117
  324.  116       continue
  325.            diml2i=j+1
  326.            noelon(i)=y
  327.            noelon(j)=x
  328.            endif
  329.  110   continue
  330.  117   continue
  331. * suprimer le trou entre les deux
  332.            if (nodz.ne.nbtot) then
  333.             do 120 i=diml2,nodz
  334.              noelon(i-nodz+nbtot)=noelon(i)
  335.  120        continue
  336.            endif
  337.  
  338. * decompte des longueurs
  339.            n1bz1=0
  340.            n1bz2=0
  341.            nbfr=0
  342.            do 130 i=1,nbtot
  343.             x=noelon(i)
  344.             nrlx=nrelong(x)
  345.             if (nrlx.eq.1) n1bz1=n1bz1+1
  346.             if (nrlx.eq.2) nbfr=nbfr+1
  347.             if (nrlx.eq.3) n1bz2=n1bz2+1
  348.  130       continue
  349. **         if (nbfr.ne.nbfra) write (6,*) ' nbfr nbfra ',nbfr,nbfra
  350. *           write (6,*) 'nbtot n1bz1 nbfr n1bz2',nbtot,n1bz1,nbfr,n1bz2
  351. *
  352. *  finalisation des infos
  353. *
  354.        londim(1)=n1bz1
  355.        londim(2)=n1bz1+nbfr
  356.        londim(3)=n1bz1+nbfr+n1bz2
  357.        long=3
  358.        dimlon=londim(3)
  359.         ideseq=abs(n1bz1-n1bz2)
  360.         xbtot=nbtot
  361.         xbfr=nbfr
  362.         xdeseq=ideseq
  363. **      fcout=(xbfr        )**4  +  (xdeseq )**3
  364.  
  365.         fcout= 2*xbtot*xbfr+xbfr*xbfr*xbfr/3.0 + xdeseq**2 
  366.      >     -icouch**3.0
  367.         fcout=max(0.d0,fcout)
  368. **      write (6,*) 'xbtot',nbtot,'xbfr',nbfr,'xdeseq',ideseq,'fcout',
  369. **   >   fcout,'icouch',icouch
  370.  
  371.         if (nbtot/(xdeseq+1).le.2)
  372.      >     fcout=fcout+xbtot*xbtot*xbtot*xbtot
  373.         if (nbfr.eq.0.and.lpiv.ne.1.and.pivot(1).ne.0.and.
  374.      >        pivot(2).ne.0) then
  375. **       write (6,*) ' frontiere vide ',n1bz1,n1bz2,nbtot,mdomn
  376.          if (mode.ne.2.and.n1bz1*n1bz2.ne.0) fcout=1
  377.         endif
  378.  
  379. *       write (6,*) ' noepr2 ',n1bz1,n1bz2,nbfr,nbtot,
  380. *    >        noelon(1),noelon(nbtot),fcout
  381.         isens=2
  382.         if (n1bz1.lt.n1bz2) isens=1
  383. *  mise a zero eventuelle de nrelong
  384.         if (mode.ne.2) then
  385.          do i=1,dimlon
  386.           x=noelon(i)
  387.           if (x.ne.0) nrelong(x)=0
  388.          enddo
  389.         endif
  390.         lfron = long
  391.         lfront=2
  392.         if (n1bz2.eq.0.or.n1bz1.eq.0) lfront=1
  393.         if (mode.eq.1) nbtot=nodz
  394.         return
  395.         end
  396.  
  397.  
  398.  
  399.  
  400.  
  401.  
  402.  
  403.  
  404.  
  405.  
  406.  
  407.  
  408.  
  409.  
  410.  
  411.  
  412.  
  413.  
  414.  
  415.  
  416.  
  417.  
  418.  
  419.  
  420.  
  421.  
  422.  
  423.  
  424.  
  425.  
  426.  
  427.  
  428.  
  429.  
  430.  
  431.  
  432.  
  433.  
  434.  

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