Sources Esope | Cast3M

duali2
C DUALI2    SOURCE    GF238795  18/02/01    21:15:21     9724           
C  DUALISE LE RESULTAT DE SURF POUR LE MAILLAGE PAR POLYGONE
C
      SUBROUTINE DUALI2(FER,XPRO,XPROJ1,IPT2,NUMELG,NDEB,NUMNP)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
      LOGICAL  PORDO
      INTEGER INTD
      SEGMENT /FER/(NFI(ITT),MAI(IPP),ITOUR)
      SEGMENT XPRO
        REAL*8 XPROJ(3,1)
      ENDSEGMENT
      POINTEUR XPROJ1.XPRO
      SEGMENT ILIST(NBNN)
      SEGMENT INB(NUMNP)
-INC SMELEME
      POINTEUR POLY.MELEME, POLY1.MELEME
*
      INTD=0
        DO 84, NUCOT = 1, ITOUR
*
          IDEB = MAI(NUCOT)
          IFIN = MAI(NUCOT+1)-1
*
          DO 84, IP2 = IDEB, IFIN
*
  84      CONTINUE
 
      IAUX=XPRO
      XPRO=XPROJ1
      XPROJ1=IAUX
      SEGINI INB
*  ON CREE UN NOEUD AU CENTRE DE GRAVITE DE CHAQUE TRIANGLE
      DO 15 I=1,NUMELG
*
        XPROJ(1,NDEB+I-1)=0.
        XPROJ(2,NDEB+I-1)=0.
        XPROJ(3,NDEB+I-1)=0.
        DO 10 J=1,3
          IP=IPT2.NUM(J,I)
          INB(IP)=INB(IP)+1
          XPROJ(1,NDEB+I-1)=XPROJ(1,NDEB+I-1)+XPROJ1.XPROJ(1,IP)
          XPROJ(2,NDEB+I-1)=XPROJ(2,NDEB+I-1)+XPROJ1.XPROJ(2,IP)
          XPROJ(3,NDEB+I-1)=XPROJ(3,NDEB+I-1)+XPROJ1.XPROJ(3,IP)
 10     CONTINUE
        XPROJ(1,NDEB+I-1)=XPROJ(1,NDEB+I-1)/3
        XPROJ(2,NDEB+I-1)=XPROJ(2,NDEB+I-1)/3
        XPROJ(3,NDEB+I-1)=XPROJ(3,NDEB+I-1)/3
 15   CONTINUE
*  ON CONSTRUIT LES ELEMENTS
      NBNN=0
      DO 20 IP=1,NUMNP
       NBNN=MAX(INB(IP),NBNN)
       INB(IP)=0
  20  CONTINUE
*
      SEGINI ILIST
      NBELEM=NUMNP
      NBSOUS=0
      NBREF=0
      SEGINI MELEME
      ITYPEL=32
      DO 35 I=1,NUMELG
       DO 30 J=1,3
        IP=IPT2.NUM(J,I)
        INB(IP)=INB(IP)+1
        NUM(INB(IP),IP)=I
  30   CONTINUE
  35  CONTINUE
*
      NUMNP  = NUMELG + NDEB - 1
      NUMELG = NBELEM
*
*  MAINTENANT IL FAUT REPASSER LES ELEMENTS POUR METTRE LES NOEUDS
*  DANS LE BON SENS ET S'OCCUPER DES BORDS
*
      DO 100 INT=1,NBELEM
*
*       Ordonnancement
*
        NUSP = 0
        PORDO = .FALSE.
*
*       TANT QUE LE POLYGONE N'EST PAS ENTIEREMENT ORDONNEE
*
  50    CONTINUE
*
*         Boucle sur tous les triangles voisins
*
          DO 70 I=1,INB(INT)
*
            ICT = NUM(I,INT)
*
*           Boucle sur les sommets du triangle associé
*
            DO 60 K=1,3
              IF (IPT2.NUM(K,ICT).EQ.INT) THEN
*
*               C'est le centre du polygone
*
                INT1 = IPT2.NUM(MOD(K,3)+1,ICT)
                INT2 = IPT2.NUM(MOD(K+1,3)+1,ICT)
*
                IF (NUSP.EQ.0) THEN
*
*                 Pas encore de sommets mémorisés
*
                  IF (INT.LT.NDEB) THEN
*
*                   Le centre du polygone est sur le coté
*
                    INP3 = NUSOM(INT1, INT, FER, NDEB)
                    INP4 = NUSOM(INT2, INT, FER, NDEB)
*
                    IF (INP3.NE.0) THEN
*
*                     Premier sommet du polygone
*
                      ILIST(1) = INP3
                      ILIST(2) = ICT + NDEB - 1
                      INTF = INT2
                      NUSP = 2
*
                      IF (INP4.NE.0) THEN
*
*                       Le polygone est triangulaire
*
                        ILIST(3) = INP4
                        PORDO = .TRUE.
                        NUSP = 3
*
                      ENDIF
*
                    ELSEIF (INP4.NE.0) THEN
*
*                     Premier sommet du polygone
*
                      ILIST(1) = INP4
                      ILIST(2) = ICT + NDEB - 1
                      INTF = INT1
                      NUSP = 2
*
                    ENDIF
                  ELSE
*
*                   Le centre du polygone est au milieu de la surface
*
                    ILIST(1) = ICT + NDEB - 1
                    INTD = INT1
                    INTF = INT2
                    NUSP = 1
*
                  ENDIF
*
                ELSE
*
*                 Des noeuds sont deja memorisés
*
                  IF (INT1.EQ.INTF.OR.INT2.EQ.INTF) THEN
*
                    NUSP = NUSP+1
                    ILIST (NUSP) = ICT + NDEB - 1
*
                    IF (INT1.EQ.INTF) THEN
                      INTF = INT2
                    ELSE IF (INT2.EQ.INTF) THEN
                      INTF = INT1
                    ENDIF
*
                    IF (INTF.EQ.INTD) THEN
*
*                     Polygone fermé
*
                      PORDO = .TRUE.
*
                    ENDIF
*
                    INP3 = NUSOM(INTF, INT, FER, NDEB)
*
                    IF (INP3.NE.0) THEN
*
*                     Le deux sommets sont voisins sur la frontiere
*                     => on ferme le polygone
*
                      NUSP = NUSP+1
                      ILIST (NUSP) = INP3
                      PORDO = .TRUE.
*
                    ENDIF
*
                  ENDIF
*
                ENDIF
*
              ENDIF
*
  60        CONTINUE
*
  70      CONTINUE
*
        IF (.NOT.PORDO) GOTO 50
*
*       Stockage du maillage dans un segment MELEME
*
        IF (INT.EQ.1) THEN
*
*         Initialisation du pointeur chapeau du maillage
*
          NBNN   = 0
          NBELEM = 0
          NBREF  = 0
          NBSOUS = 1
          SEGINI POLY1
*
        ELSE
*
*         Recherche si un polygone a NUSP cotés existe deja dans MELEME
*
          NBELEM = 0
*
          DO 80 I=1, POLY1.LISOUS(/1)
*
            POLY = POLY1.LISOUS(I)
*
            IF (POLY.NUM(/1).EQ.NUSP) THEN
*
              NBELEM = POLY.NUM(/2)+1
              NBNN   = NUSP
              NBSOUS = 0
              NBREF  = 0
*
              SEGADJ POLY
              GOTO 81
*
            ENDIF
*
  80      CONTINUE
  81      CONTINUE
*
          IF (NBELEM.EQ.0) THEN
*
            NBNN = 0
            NBELEM = 0
            NBREF  = 0
            NBSOUS = POLY1.LISOUS(/1)+1
            SEGADJ POLY1
*
          ENDIF
*
        ENDIF
*
        IF (NBELEM.EQ.0) THEN
*
*         Creation de l'element a NUSP cote
*
          NBELEM = 1
          NBNN   = NUSP
          NBSOUS = 0
          NBREF  = 0
*
          SEGINI POLY
*
          NBSOUS = POLY1.LISOUS(/1)
          POLY1.LISOUS(NBSOUS) = POLY
          POLY.ITYPEL = 32
*
        ENDIF
*
*       Recopie des données dans le MELEME
*
        DO 90 I = 1, NUSP
*
          POLY.NUM(I, NBELEM) = ILIST(I)
*
  90    CONTINUE
*
 100  CONTINUE
*
*     Recopie du nouveau MELEME dans l'ancien
*
      IPT3 = IPT2
      IPT2 = POLY1
*
      SEGSUP IPT3
*
      END