Télécharger la plaquette Cast3M


fr en

Télécharger part7.eso

Retour à la liste

Numérotation des lignes : 

part7

  1. C PART7     SOURCE    GOUNAND   26/01/09    21:15:47     12442          
  2. ************************************************************************
  3. * NOM         : PART7
  4. * DESCRIPTION : Sous-programme dedie a la separation en composantes
  5. *               connexes d'un maillage + regles supplementaires de
  6. *               separation en differentes zones
  7. ************************************************************************
  8. * APPELE PAR : part.eso ; ccon.eso (obsolete)
  9. ************************************************************************
  10. * ENTREES :: MEL1  = pointeur sur le maillage a partitionner
  11. *            KLI   > 0 si option 'LIGN'
  12. *            KFA   > 0 si option 'FACE'
  13. *            KMA   > 0 si option 'MAIL'
  14. *            MEL2  = pointeur sur le maillage separateur (option 'MAIL')
  15. *            KAN   > 0 si option 'ANGL'
  16. *            ANG   = valeur seuil pour l'angle (option 'ANGL')
  17. *            ITQ   > 0 si mot-cle 'TELQ' present (option 'ANGL')
  18. *            KESCL > 0 si besoin des indices SOUSTYPE et CREATEUR
  19. * SORTIES :: ITAB  = pointeur vers la table de partitionnement
  20. ************************************************************************
  21.       SUBROUTINE PART7(MEL1,KLI,KFA,KMA,MEL2,KAN,ANG,ITQ,ITAB,KESCL)
  22.  
  23.       IMPLICIT INTEGER(I-N)
  24.       IMPLICIT REAL*8 (A-H,O-Z)
  25.  
  26.  
  27. -INC PPARAM
  28. -INC CCOPTIO
  29. -INC CCGEOME
  30. -INC SMCOORD
  31. -INC SMELEME
  32. -INC SMTABLE
  33. -INC SMCHAML
  34. -INC SMMODEL
  35.  
  36.       SEGMENT JMEM(NODES+1)
  37. C     JMEM CONTIENT LE NOMBRE D'ELEMENTS AUQUEL APPARTIENT CHAQUE NOEUD
  38. C     PUIS LA POSITION DU PREMIER ELEMENT DANS IMEMO ET LMEMO
  39.  
  40.       SEGMENT ICPR(nbpts)
  41. C     ICPR(I) DONNE LE NUMERO LOCAL (DANS LES TABLEAUX DE LA PRESENTE
  42. C     SUBROUTINE) DU I-EME NOEUD GLOBAL (DANS LA TABLE MCOORD)
  43.  
  44.       SEGMENT IMEMO(NBV),LMEMO(NBV)
  45. C     CONTIENT LA LISTE DES ELEMENTS APPARTENANT AU NOEUD 1, 2, 3...
  46. C     (IMEMO => NUMERO DE L'ELEMENT ET LMEMO => NUMERO DU LISOUS)
  47.  
  48.       SEGMENT LISIND(NBS)
  49. C     POINTE VERS LES SEGMENTS INDIC ASSOCIES A CHAQUE SOUS-MAILLAGE
  50.  
  51.       SEGMENT JMEM2(NODES2+1),ICPR2(nbpts),IMEMO2(NBV2),
  52.      &        LMEMO2(NBV2)
  53. C     IDEM, MAIS POUR LE MAILLAGE SEPARATEUR
  54.  
  55.       SEGMENT INDIC(NBEL)
  56. C     INDICATEUR DU NUMERO DE ZONE
  57.  
  58.       SEGMENT LISCO1(NELTOT),LISCO2(NELTOT)
  59. C     LISTES DES ELEMENTS VOISINS
  60.  
  61.       SEGMENT LISIN(NNOMAX)
  62. C     LISTE DES NOEUDS A L'INTERFACE ENTRE DEUX ELEMENTS VOISINS
  63.  
  64.       SEGMENT MIELVA
  65.          INTEGER IELVAX(N1)
  66.          INTEGER IELVAY(N1)
  67.          INTEGER IELVAZ(N1)
  68.       ENDSEGMENT
  69. C     POINTEURS VERS LES SEGMENTS MELVAL (OPTION 'ANGL')
  70.  
  71.       LOGICAL KDIM1,KDIM2,KDIM3
  72.       INTEGER NNOMAX
  73.  
  74.  
  75.  
  76.  
  77. *     +---------------------------------------------------------------+
  78. *     |                                                               |
  79. *     |                 I N I T I A L I S A T I O N S                 |
  80. *     |                                                               |
  81. *     +---------------------------------------------------------------+
  82.  
  83.  
  84.  
  85. *     VERIFICATION QUE LE MAILLAGE EST COMPATIBLE AVEC LES OPTIONS
  86. *     FOURNIES
  87. *     ************************************************************
  88.       NNOMAX=0
  89.       MELEME=MEL1
  90.       SEGACT,MELEME
  91.       IPT1=MELEME
  92.  
  93.       IDIM1=0
  94.       IDIM2=0
  95.       IDIM3=0
  96.       DO ISO=1,MAX(1,LISOUS(/1))
  97.          IF (LISOUS(/1).GT.1) THEN
  98.             IPT1=LISOUS(ISO)
  99.             SEGACT,IPT1
  100.          ENDIF
  101.  
  102.          ITY=IPT1.ITYPEL
  103.          NNOMAX=MAX(NNOMAX,NBNNE(ITY))
  104.  
  105. *         KDIM1=(LDLR(ITY).EQ.1.AND.ITY.NE.12.AND.ITY.NE.13)
  106. *         KDIM2=(ITY.EQ.KSURF(ITY))
  107. *         KDIM3=(LDLR(ITY).EQ.3.AND.ITY.NE.30.AND.ITY.NE.31)
  108.          KDIM1=(LDLR(ITY).EQ.1)
  109.          KDIM2=(LDLR(ITY).EQ.2)
  110.          KDIM3=(LDLR(ITY).EQ.3)
  111.  
  112. *        Element special type 'MULT' ou 'SUPE'
  113.          IF (.NOT.(KDIM1.OR.KDIM2.OR.KDIM3)) THEN
  114.             CALL ERREUR(16)
  115.             RETURN
  116.          ENDIF
  117.  
  118.          IF (KDIM1) IDIM1=IDIM1+1
  119.          IF (KDIM2) IDIM2=IDIM2+1
  120.          IF (KDIM3) IDIM3=IDIM3+1
  121.       ENDDO
  122.  
  123.       IF ((KFA.GT.0.AND.(IDIM1.GT.0.OR.IDIM3.GT.0)).OR.
  124.      &    (KLI.GT.0.AND.(IDIM2.GT.0.OR.IDIM3.GT.0)).OR.
  125.      &    (KAN.GT.0.AND.IDIM3.GT.0)) THEN
  126.          CALL ERREUR(16)
  127.          RETURN
  128.       ENDIF
  129.  
  130.  
  131. *     OPTION 'ANGL' => CREATION DES TABLEAUX DONNANT LE VECTEUR
  132. *     NORMAL/TANGENT A CHAQUE ELEMENT
  133. *     *********************************************************
  134.  
  135.       IF (KAN.GT.0) THEN
  136.  
  137. *        Transformation en un maillage lineaire
  138.          CALL ECROBJ('MAILLAGE',MELEME)
  139.          CALL CHANLI
  140.          IF (IERR.NE.0) RETURN
  141.          CALL LIROBJ('MAILLAGE',MELEME,1,IRET)
  142.          IF (IERR.NE.0) RETURN
  143.  
  144. *        Creation d'un objet MMODEL temporaire (le type est sans
  145. *        importance)
  146.          CALL ECRCHA('POUT')
  147.          CALL ECRCHA('COQ4')
  148.          CALL ECRCHA('COQ3')
  149.          CALL ECRCHA('ELASTIQUE')
  150.          CALL ECRCHA('MECANIQUE')
  151.          CALL ECROBJ('MAILLAGE',MELEME)
  152.          CALL MODELI
  153.          IF (IERR.NE.0) RETURN
  154.          CALL LIROBJ('MMODEL  ',IPMODL,1,IRET)
  155.          CALL ACTOBJ('MMODEL  ',IPMODL,1)
  156.          IF (IERR.NE.0) RETURN
  157.  
  158. *        Calcul des vecteurs directeurs
  159.          SEGACT,MCOORD
  160.          CALL JACONO(IPMODL,1,IPCHE5,IRET)
  161.          IF (IERR.NE.0) RETURN
  162.  
  163. *        On ne garde qu'une seule valeur par element (=> type GRAVITE)
  164.          CALL CHASUP(IPMODL,IPCHE5,IPCHE2,IRET,2)
  165.          IF (IERR.NE.0) RETURN
  166.          IF (IRET.NE.0) THEN
  167.             CALL ERREUR(IRET)
  168.             RETURN
  169.          ENDIF
  170.  
  171. *        On recupere le MELEME pour etre sur d'avoir le bon ordre de
  172. *        description des elements
  173.          CALL ECRCHA('MAIL')
  174.          CALL ECROBJ('MCHAML',IPCHE2)
  175.          CALL EXTRAI
  176.          IF (IERR.NE.0) RETURN
  177.          CALL LIROBJ('MAILLAGE',MELEME,1,IRET)
  178.          IF (IERR.NE.0) RETURN
  179.  
  180. *        Suppression de l'objet MMODEL
  181.          MMODEL=IPMODL
  182.          SEGACT,MMODEL
  183.          DO K=1,KMODEL(/1)
  184.             IMODEL=KMODEL(K)
  185.             SEGSUP,IMODEL
  186.          ENDDO
  187.          SEGSUP,MMODEL
  188.  
  189. *        Remplissage du segment MIELVA (pointeurs vers les MELVAL)
  190.          MCHELM=IPCHE2
  191.          SEGACT,MCHELM
  192.          N1=ICHAML(/1)
  193.          SEGINI,MIELVA
  194.          DO I=1,N1
  195.             MCHAML=ICHAML(I)
  196.             SEGACT,MCHAML
  197.  
  198.             IELVAX(I)=IELVAL(1)
  199.             MELVAL=IELVAX(I)
  200.             SEGACT,MELVAL
  201.  
  202.             IELVAY(I)=IELVAL(2)
  203.             MELVAL=IELVAY(I)
  204.             SEGACT,MELVAL
  205.  
  206.             IF (IDIM.EQ.3) THEN
  207.                IELVAZ(I)=IELVAL(3)
  208.                MELVAL=IELVAZ(I)
  209.                SEGACT,MELVAL
  210.             ENDIF
  211.  
  212.             SEGSUP,MCHAML
  213.          ENDDO
  214.          SEGSUP,MCHELM
  215.  
  216.       ENDIF
  217.  
  218.  
  219. *     CORRESPONDANCE ENTRE LES NUMEROTATIONS LOCALE/GLOBALE
  220. *     *****************************************************
  221.  
  222.       SEGINI,ICPR
  223.  
  224.       SEGACT,MELEME
  225.       NBSOUS=LISOUS(/1)
  226.       NBS=MAX(1,NBSOUS)
  227.       IPT1=MELEME
  228.  
  229. *     Boucle sur les eventuels sous-maillages
  230.       IKOU=0
  231.       DO 100 IO=1,NBS
  232.          IF (NBSOUS.GT.0) THEN
  233.             IPT1=LISOUS(IO)
  234.             SEGACT,IPT1
  235.          ENDIF
  236.  
  237. *        Remplissage du tableau de correspondance ICPR
  238.          DO 150 J=1,IPT1.NUM(/2)
  239.             DO 151 I=1,IPT1.NUM(/1)
  240.                IJ=IPT1.NUM(I,J)
  241.                IF (ICPR(IJ).EQ.0) THEN
  242.                   IKOU=IKOU+1
  243.                   ICPR(IJ)=IKOU
  244.                ENDIF
  245.  151        CONTINUE
  246.  150     CONTINUE
  247.  
  248.  100  CONTINUE
  249.  
  250. *     Nombre de noeuds distincts dans le maillage
  251.       NODES=IKOU
  252.  
  253. *     MAILLAGE vide => on sort des maintenant
  254.       IF (NODES.EQ.0) THEN
  255.          M=0
  256.          SEGINI,MTABLE
  257.          ITAB=MTABLE
  258.          MLOTAB=0
  259.          GOTO 9999
  260.       ENDIF
  261.  
  262. *
  263. *     IDENTIFICATION DES ELEMENTS OU APPARAISSENT TOUS LES NOEUDS
  264. *          =>  IMEMO = NUMERO DE ELEMENT
  265. *          =>  LMEMO = NUMERO DU SOUS-MAILLGE
  266. *     JMEM(I)+1 IDENTIFIE LA POSITION DANS IMEMO/LMEMO DU PREMIER
  267. *     ELEMENT ASSOCIE AU NOEUD I
  268. *     ***************************************************************
  269.  
  270.       SEGINI,JMEM
  271.  
  272. *     On compte combien de fois chaque noeud apparait dans le maillage
  273.       DO 200 IO=1,NBS
  274.          IF (NBSOUS.GT.0) IPT1=LISOUS(IO)
  275.          DO 250 J=1,IPT1.NUM(/2)
  276.             DO 251 I=1,IPT1.NUM(/1)
  277.                IJ=ICPR(IPT1.NUM(I,J))
  278.                JMEM(IJ)=JMEM(IJ)+1
  279.  251        CONTINUE
  280.  250     CONTINUE
  281.  200  CONTINUE
  282. *
  283. *     On en deduit par cumul la position de depart dans IMEMO/LMEMO
  284.       DO 290 I=1+1,NODES+1
  285.          JMEM(I)=JMEM(I)+JMEM(I-1)
  286. 290   CONTINUE
  287.       NBV=JMEM(NODES)
  288. *
  289. *     Remplissage de IMEMO et LMEMO
  290.       SEGINI,IMEMO,LMEMO
  291.       DO 300 IO=1,NBS
  292.          IF (NBSOUS.GT.0) IPT1=LISOUS(IO)
  293.          DO 350 J=1,IPT1.NUM(/2)
  294.             DO 351 I=1,IPT1.NUM(/1)
  295.                IJ=ICPR(IPT1.NUM(I,J))
  296.                IMEMO(JMEM(IJ))=J
  297.                LMEMO(JMEM(IJ))=IO
  298.                JMEM(IJ)=JMEM(IJ)-1
  299.  351        CONTINUE
  300.  350     CONTINUE
  301.  300  CONTINUE
  302.  
  303.  
  304. *     OPTION 'MAIL' => ON REMPLIT DE LA MEME MANIERE ICPR2, JMEM2,
  305. *     IMEMO2 ET LMEMO2
  306. *     ************************************************************
  307.  
  308.       IF (KMA.GT.0) THEN
  309.  
  310. *        Tableau ICPR2
  311.          SEGINI,ICPR2
  312.          IPT2=MEL2
  313.          SEGACT,IPT2
  314.          NBSOU2=IPT2.LISOUS(/1)
  315.          NBS2=MAX(1,NBSOU2)
  316.          IPT5=IPT2
  317.          IKOU=0
  318.          DO 400 IO=1,NBS2
  319.             IF (NBSOU2.GT.0) THEN
  320.                IPT5=IPT2.LISOUS(IO)
  321.                SEGACT,IPT5
  322.             ENDIF
  323.             DO 401 J=1,IPT5.NUM(/2)
  324.                DO 402 I=1,IPT5.NUM(/1)
  325.                   IJ=IPT5.NUM(I,J)
  326.                   IF (ICPR2(IJ).EQ.0) THEN
  327.                      IKOU=IKOU+1
  328.                      ICPR2(IJ)=IKOU
  329.                   ENDIF
  330.  402           CONTINUE
  331.  401        CONTINUE
  332.  400     CONTINUE
  333.          NODES2=IKOU
  334.  
  335. *        MAILLAGE vide => l'option 'MAIL' est desactivee
  336.          IF (NODES2.EQ.0) THEN
  337.             DO 410 IO=1,NBS2
  338.                IF (NBSOU2.GT.0) IPT5=IPT2.LISOUS(IO)
  339.                SEGDES,IPT5
  340. 410         CONTINUE
  341.             SEGSUP,ICPR2
  342.             KMA=0
  343.             GOTO 499
  344.          ENDIF
  345.  
  346. *        Tableaux JMEM2, IMEMO2 et LMEMO2
  347.          SEGINI,JMEM2
  348.          DO 420 IO=1,NBS2
  349.             IF (NBSOU2.GT.0) IPT5=IPT2.LISOUS(IO)
  350.             DO 421 J=1,IPT5.NUM(/2)
  351.                DO 422 I=1,IPT5.NUM(/1)
  352.                   IJ=ICPR2(IPT5.NUM(I,J))
  353.                   JMEM2(IJ)=JMEM2(IJ)+1
  354.  422           CONTINUE
  355.  421        CONTINUE
  356.  420     CONTINUE
  357.          DO 430 I=1+1,NODES2+1
  358.             JMEM2(I)=JMEM2(I)+JMEM2(I-1)
  359.  430     CONTINUE
  360.          NBV2=JMEM2(NODES2)
  361.          SEGINI,IMEMO2,LMEMO2
  362.          DO 440 IO=1,NBS2
  363.             IF (NBSOU2.GT.0) IPT5=IPT2.LISOUS(IO)
  364.             DO 441 J=1,IPT5.NUM(/2)
  365.                DO 442 I=1,IPT5.NUM(/1)
  366.                   IJ=ICPR2(IPT5.NUM(I,J))
  367.                   IMEMO2(JMEM2(IJ))=J
  368.                   LMEMO2(JMEM2(IJ))=IO
  369.                   JMEM2(IJ)=JMEM2(IJ)-1
  370.  442           CONTINUE
  371.  441        CONTINUE
  372.  440     CONTINUE
  373.  
  374. *        On cree aussi le segment LISIN qui servira plus bas
  375.          SEGINI,LISIN
  376.  
  377.       ENDIF
  378. 499   CONTINUE
  379.  
  380.  
  381. *     CREATION D'UN SEGMENT INDIC POUR CHAQUE LISOUS
  382. *     **********************************************
  383.  
  384.       SEGINI,LISIND
  385.       NELTOT=0
  386.       DO 500 IO=1,NBS
  387.          IF (NBSOUS.NE.0) IPT1=LISOUS(IO)
  388.          NBEL=IPT1.NUM(/2)
  389.          NELTOT=NELTOT+NBEL
  390.          SEGINI,INDIC
  391.          LISIND(IO)=INDIC
  392.  500  CONTINUE
  393.       SEGINI LISCO1,LISCO2
  394.  
  395.  
  396.  
  397. *     +---------------------------------------------------------------+
  398. *     |                                                               |
  399. *     |          C O N S T R U C T I O N   D E S   Z O N E S          |
  400. *     |                                                               |
  401. *     +---------------------------------------------------------------+
  402.  
  404.       IOC=1
  405.       IELC=0
  406.  
  407. *     LABEL 1000 : PARCOURS DES SEGMENTS INDIC A LA RECHERCHE D'UN
  408. *                  ELEMENT ENCORE NON ATTRIBUE
  409. *     ************************************************************
  410. 1000  CONTINUE
  411.  
  413.       IELC=IELC+1
  414.       DO 1010 IO=IOC,NBS
  415.          IF (NBSOUS.NE.0) IPT1=LISOUS(IO)
  416.          INDIC=LISIND(IO)
  417.          DO 1020 IEL=IELC,IPT1.NUM(/2)
  418.             IF (INDIC(IEL).EQ.0) GOTO 1030
  419. 1020     CONTINUE
  420.          IELC=1
  421. 1010  CONTINUE
  422.  
  423. C     TOUS LES ELEMENTS ONT ETE CLASSES => ON A FINI
  424.       GOTO 1500
  425.  
  426. *     ON A TROUVE UN ELEMENT DE DEPART D'UNE NOUVELLE ZONE
  427. *     => ON VA ETENDRE AUX ELEMENTS VOISINS
  428. *     ****************************************************
  429. 1030  CONTINUE
  430.       IOC=IO
  431.       IELC=IEL
  432.  
  433. *     On attribue une zone a l'element trouve
  435.  
  436. *     ILRMP = Nombre d'elements ajoutes a LISCO1/LISCO2
  437. *     ILEXT = Nombre d'elements parcourus dans LISCO1/LISCO2
  438.       ILRMP=1
  439.       ILEXT=1
  440.  
  441. *     On reinitialise LISCO1/LISCO2 avec seulement cet element
  442.       LISCO1(ILRMP)=IO
  443.       LISCO2(ILRMP)=IEL
  444.  
  445. *     BOUCLE DE REMPLISSAGE DE LISCO1/LISCO2, DE VOISIN EN VOISIN
  446. *     ***********************************************************
  447.  
  448. *     Label 1120 => element suivant dans les listes LISCO1/LISCO2
  449. 1120  CONTINUE
  450.       IF (ILEXT.GT.ILRMP) GOTO 1130
  451.  
  452.       ION=LISCO1(ILEXT)
  453.       IEL=LISCO2(ILEXT)
  454.       IF (NBSOUS.NE.0) IPT1=LISOUS(ION)
  455.       IF (KAN.GT.0) THEN
  456. *        Vecteur directeur de l'element 1 (option 'ANGL')
  457.          MELVAL=IELVAX(ION)
  458.          IELMIN=MIN(IEL,MELVAL.VELCHE(/2))
  459.          X1=VELCHE(1,IELMIN)
  460.          MELVAL=IELVAY(ION)
  461.          IELMIN=MIN(IEL,MELVAL.VELCHE(/2))
  462.          Y1=VELCHE(1,IELMIN)
  463.          IF (IDIM.EQ.3) THEN
  464.             MELVAL=IELVAZ(ION)
  465.             IELMIN=MIN(IEL,MELVAL.VELCHE(/2))
  466.             Z1=VELCHE(1,IELMIN)
  467.          ENDIF
  468.       ENDIF
  469.  
  470. *     Label 1100 => noeud IP suivant de l'element courant
  471.       DO 1100 IN=1,IPT1.NUM(/1)
  472.          IP=ICPR(IPT1.NUM(IN,IEL))
  473.  
  474. *        Label 1110 => voisin suivant via le noeud IP
  475.          DO 1110 KK=JMEM(IP)+1,JMEM(IP+1)
  476.             JON=LMEMO(KK)
  477.             JEL=IMEMO(KK)
  478.             INDIC=LISIND(JON)
  479.  
  480.  
  481. *           TESTS SUR L'ELEMENT VOISIN (JON;JEL) : SI L'UN DES TESTS
  482. *           ECHOUE, ALORS CET ELEMENT N'APPARTIENT PAS A CETTE ZONE
  483. *           ********************************************************
  484.  
  485. *           1) CONDITION SINE QUA NONE : IL N'A PAS DEJA ETE ATTRIBUE
  486. *              A UNE AUTRE ZONE
  487. *              ======================================================
  488.             IF (INDIC(JEL).NE.0) GOTO 1110
  489.  
  490. *           2) OPTION 'FACE' (UNIQUEMENT POUR LES MAILLAGES DE SURFACES)
  491. *              =========================================================
  492.             IF (KFA.GT.0) THEN
  493.                IF (NBSOUS.NE.0) THEN
  494.                   IPT3=LISOUS(JON)
  495.                ELSE
  496.                   IPT3=MELEME
  497.                ENDIF
  498.  
  499. *              a) Verification que les elements ont au moins 1 autre
  500. *                 noeud que IP en commun (attention : on ne verifie
  501. *                 pas qu'ils appartiennent a une meme arete)
  502.                DO 1150 I1=1,IPT1.NUM(/1)
  503.                   IP1=ICPR(IPT1.NUM(I1,IEL))
  504.                   IF (IP1.EQ.IP) GOTO 1150
  505.                   DO 1160 I2=1,IPT3.NUM(/1)
  506.                      IP2=ICPR(IPT3.NUM(I2,JEL))
  507.                      IF (IP1.EQ.IP2) GOTO 1170
  508. 1160              CONTINUE
  509. 1150           CONTINUE
  510.                GOTO 1110
  511.  
  512. *              b) Verification qu'il n'y a que 2 elements qui
  513. *                 contiennent les noeuds IP et IP1
  514. 1170           CONTINUE
  515.                NL=0
  516.                DO 1180 K1=JMEM(IP)+1,JMEM(IP+1)
  517.                   IF (K1.EQ.KK) GOTO 1180
  518.                   I1=LMEMO(K1)
  519.                   J1=IMEMO(K1)
  520.                   DO 1190 K2=JMEM(IP1)+1,JMEM(IP1+1)
  521.                      I2=LMEMO(K2)
  522.                      J2=IMEMO(K2)
  523.                      IF (I1.EQ.I2.AND.J1.EQ.J2) NL=NL+1
  524. 1190              CONTINUE
  525. 1180           CONTINUE
  526.                IF (NL.NE.1) GOTO 1110
  527.             ENDIF
  528.  
  529. *           3) OPTION 'LIGN' (UNIQUEMENT POUR LES MAILLAGES DE LIGNES)
  530. *              VERIFICATION QU'IL N'Y A QUE 2 ELEMENTS QUI CONTIENNENT
  531. *              LE NOEUD IP
  532. *              =======================================================
  533.             IF (KLI.GT.0) THEN
  534.                IF (JMEM(IP+1)-JMEM(IP).NE.2) GOTO 1110
  535.             ENDIF
  536.  
  537. *           4) OPTION 'ANGL' (POUR LES MAILLAGES DE LIGNES ET/OU DE
  538. *              SURFACE) : VERIFICATION QUE L'ANGLE ENTRE 2 ELEMENTS
  539. *              VOISINS EST INFERIEUR A UNE VALEUR SEUIL
  540. *              ====================================================
  541.             IF (KAN.GT.0) THEN
  542.  
  543. *              Vecteur directeur de l'element 2
  544. *              (vecteur directeur de l'element 1 sorti de la boucle)
  545.                MELVAL=IELVAX(JON)
  546.                JELMIN=MIN(JEL,MELVAL.VELCHE(/2))
  547.                X2=VELCHE(1,JELMIN)
  548.                MELVAL=IELVAY(JON)
  549.                JELMIN=MIN(JEL,MELVAL.VELCHE(/2))
  550.                Y2=VELCHE(1,JELMIN)
  551.  
  552. *              Produit scalaire et norme
  553.                XN1=X1*X1+Y1*Y1
  554.                XN2=X2*X2+Y2*Y2
  555.                CA=X1*X2+Y1*Y2
  556.  
  557. *              Prise en compte 3eme direction le cas echeant
  558.                IF (IDIM.EQ.3) THEN
  559.                   MELVAL=IELVAZ(JON)
  560.                   JELMIN=MIN(JEL,MELVAL.VELCHE(/2))
  561.                   Z2=VELCHE(1,JELMIN)
  562.  
  563.                   XN1=XN1+Z1*Z1
  564.                   XN2=XN2+Z2*Z2
  565.                   CA=CA+Z1*Z2
  566.                ENDIF
  567.  
  568. *              Determination de l'angle en degres entre les 2 vecteurs
  569.                CA=CA/((XN1*XN2)**0.5)
  570.                IF (CA.GT.1.D0) CA=1.D0
  571.                IF (CA.LT.-1.D0) CA=-1.D0
  572.  
  573.                IF ((ITQ.GT.0.OR.(ITQ.LE.0.AND.CA.NE.-1.D0)).AND.
  574.      &             (CA.LT.ANG)) GOTO 1110
  575.             ENDIF
  576.  
  577. *           5) OPTION 'MAIL' : VERIFICATION QUE L'INTERFACE COMMUNE
  578. *              ENTRE 2 ELEMENTS VOISINS N'APPARTIENT PAS A MEL2
  579. *              ====================================================
  580.             IF (KMA.GT.0) THEN
  581.  
  582. *              Test rapide grace au noeud commun deja connu (IP et IPP
  583. *              sont les numeros locaux du meme noeud dans MEL1 et MEL2)
  584.                IPP=ICPR2(IPT1.NUM(IN,IEL))
  585.                IF (IPP.EQ.0) GOTO 999
  586.  
  587. *              Nb. de noeuds a l'interface entre ION/IEL et JON/JEL
  588. *              (IMPOSSIBLE A SAVOIR A PRIORI => EXEMPLE : CUB8/PY5)
  589.                IF (NBSOUS.NE.0) THEN
  590.                   IPT3=LISOUS(JON)
  591.                ELSE
  592.                   IPT3=MELEME
  593.                ENDIF
  594.                NBNIN=0
  595.                DO 1200 I1=1,IPT1.NUM(/1)
  596.                   IP1=ICPR(IPT1.NUM(I1,IEL))
  597.                   IF (IP1.EQ.IP) GOTO 1200
  598.                   DO I2=1,IPT3.NUM(/1)
  599.                      IP2=ICPR(IPT3.NUM(I2,JEL))
  600.                      IF (IP1.EQ.IP2) THEN
  601. *                       On a trouve un noeud de l'interface, mais s'il
  602. *                       n'est pas dans MEL2 => inutile d'aller plus loin
  603.                         IF (ICPR2(IPT3.NUM(I2,JEL)).EQ.0) GOTO 999
  604. *                       Sinon on le memorise et on en cherche d'autres
  605.                         NBNIN=NBNIN+1
  606.                         LISIN(NBNIN)=IPT3.NUM(I2,JEL)
  607.                         GOTO 1200
  608.                      ENDIF
  609.                   ENDDO
  610. 1200           CONTINUE
  611.  
  612. *              A ce stade, on connait tous les noeuds a l'interface
  613. *              entre les 2 elements voisins, et on sait qu'ils sont
  614. *              tous dans MEL2 => IL RESTE A VERIFIER QU'ILS SONT DANS
  615. *              UN MEME ELEMENT DE MEL2 (la liste des possibilites
  616. *              est reduite grace aux tableaux JMEM2/IMEMO2/LMEMO2)
  617.                DO 1210 K2=JMEM2(IPP)+1,JMEM2(IPP+1)
  618.                   KON=LMEMO2(K2)
  619.                   KEL=IMEMO2(K2)
  620.                   IF (NBSOU2.NE.0) THEN
  621.                      IPT5=IPT2.LISOUS(KON)
  622.                   ELSE
  623.                      IPT5=IPT2
  624.                   ENDIF
  625.  
  626. *                 Inutile de tester cet element de MEL2 s'il n'a pas
  627. *                 assez de noeuds...
  628.                   IF (NBNNE(IPT5.ITYPEL).LT.NBNIN+1) GOTO 1210
  629.  
  630. *                 Test de tous les noeuds de l'interface
  631.                   DO 1220 K1=1,NBNIN
  632.                      INO=LISIN(K1)
  633.                      IF (ICPR2(INO).EQ.IPP) GOTO 1220
  634.                      DO K3=1,IPT5.NUM(/1)
  635.                         IF (INO.EQ.IPT5.NUM(K3,KEL)) GOTO 1220
  636.                      ENDDO
  637.                      GOTO 1210
  638. 1220              CONTINUE
  639.  
  640. *                 => L'INTERFACE ENTRE LES DEUX VOISINS EST INCLUSE
  641. *                    DANS UN ELEMENT DE MEL2
  642.                   GOTO 1110
  643.  
  644. 1210           CONTINUE
  645.  
  646.             ENDIF
  647.  
  648. *           => TOUS LES TESTS SONT PASSES : ON AJOUTE L'ELEMENT AUX
  649. *              LISTES LISCO1/LISCO2 ET ON LUI ATTRIBUE LA ZONE COURANTE
  650. 999         CONTINUE
  652.             ILRMP=ILRMP+1
  653.             LISCO1(ILRMP)=JON
  654.             LISCO2(ILRMP)=JEL
  655.  
  656. 1110     CONTINUE
  657. 1100  CONTINUE
  658.  
  659. *     Tous les voisins de l'element courant ont ete testes : on va
  660. *     regarder s'il reste des elements dans LISCO1/LISCO2
  661.       ILEXT=ILEXT+1
  662.       GOTO 1120
  663.  
  664. *     ON A FINI, ON ENCHAINE SUR LA ZONE SUIVANTE
  665. 1130  CONTINUE
  666.       GOTO 1000
  667.  
  668.  
  669. *     +---------------------------------------------------------------+
  670. *     |                                                               |
  671. *     |   C R E A T I O N   D E   L ' O B J E T   D E   S O R T I E   |
  672. *     |                                                               |
  673. *     +---------------------------------------------------------------+
  674.  
  675. 1500  CONTINUE
  676.       NBCOMP=NBCOMP-1
  677.  
  678. *     CREATION DE L'OBJET TABLE DE SORTIE
  679. *     ***********************************
  681.       IF (KESCL.GT.0) M=M+2
  682.       SEGINI,MTABLE
  683.       ITAB=MTABLE
  684.       IF (KESCL.GT.0) THEN
  685.          CALL ECCTAB(ITAB,'MOT',0,0.D0,'SOUSTYPE',.TRUE.,
  686.      &                     0,'MOT',0,0.D0,'ESCLAVE',.TRUE.,0)
  687.          CALL ECCTAB(ITAB,'MOT',0,0.D0,'CREATEUR',.TRUE.,
  688.      &                     0,'MOT',0,0.D0,'PART',.TRUE.,0)
  689.       ENDIF
  690.  
  691.  
  692. *     CREATION DES MAILLAGES DES DIFFERENTES ZONES TROUVEES
  693. *     *****************************************************
  694.       DO 2000 ICOMP=1,NBCOMP
  695.          ISS=0
  696.  
  697.          NBNN=0
  698.          NBELEM=0
  699.          NBSOUS=0
  700.          NBREF=0
  701.          SEGINI,IPT7
  702.  
  703. *        Boucle sur les LISOUS du maillage initial
  704.          DO 2010 IS=1,NBS
  705.             IF (LISOUS(/1).NE.0) IPT1=LISOUS(IS)
  706.             INDIC=LISIND(IS)
  707.  
  708. *           Decompte du nombre d'elements du LISOUS initial appartenant
  709. *           a la composante courante
  710.             IELEM=0
  711.             DO 2001 I=1,INDIC(/1)
  712.                IF (INDIC(I).EQ.ICOMP) IELEM=IELEM+1
  713. 2001        CONTINUE
  714.  
  715. *           Si besoin, on cree et on remplit le LISOUS correspondant
  716. *           dans IPT7
  717.             IF (IELEM.NE.0) THEN
  718.                ISS=ISS+1
  719.                NBNN=IPT1.NUM(/1)
  720.                NBELEM=IELEM
  721.                NBSOUS=0
  722.                NBREF=0
  723.                SEGINI,IPT3
  724.                IPT3.ITYPEL=IPT1.ITYPEL
  725.                IELEM=0
  726.                DO 2005 I=1,INDIC(/1)
  727.                   IF (INDIC(I).NE.ICOMP) GOTO 2005
  728.                   IELEM=IELEM+1
  729.                   DO 2006 J=1,IPT1.NUM(/1)
  730.                      IPT3.NUM(J,IELEM)=IPT1.NUM(J,I)
  731. 2006              CONTINUE
  732.                   IPT3.ICOLOR(IELEM)=IPT1.ICOLOR(I)
  733. 2005           CONTINUE
  734.                NBNN=0
  735.                NBELEM=0
  736.                NBSOUS=IPT7.LISOUS(/1)+1
  737.                NBREF=0
  738.                SEGADJ,IPT7
  739.                IPT7.LISOUS(NBSOUS)=IPT3
  740.                SEGDES,IPT3
  741.             ENDIF
  742. 2010     CONTINUE
  743. *
  744. *        S'il n'y a qu'un seul LISOUS, on modifie la structure de IPT7
  745.          IF (IPT7.LISOUS(/1).EQ.1) THEN
  746.             IPT=IPT7.LISOUS(1)
  747.             SEGSUP,IPT7
  748.             IPT7=IPT
  749.          ENDIF
  750.          SEGDES,IPT7
  751.  
  752. *        Ajout de IPT7 a l'indice ICOMP de MTABLE
  753.          CALL ECCTAB(MTABLE,'ENTIER  ',ICOMP,CC,' ',.TRUE.,I,
  754.      &                      'MAILLAGE',I,XX,' ',.TRUE.,IPT7)
  755. 2000  CONTINUE
  756.  
  757. *     FIN DE LA SUBROUTINE : UN PEU DE MENAGE...
  758. *     ******************************************
  759.       SEGSUP,JMEM,LMEMO,IMEMO,LISCO1,LISCO2
  760.       DO I=1,LISIND(/1)
  761.          INDIC=LISIND(I)
  762.          SEGSUP,INDIC
  763.       ENDDO
  764.       SEGSUP,LISIND
  765.       IF (KMA.NE.0) THEN
  766.          SEGSUP,ICPR2,JMEM2,LMEMO2,IMEMO2,LISIN
  767.          IF (NBSOU2.GT.0) THEN
  768.             DO I=1,IPT2.LISOUS(/1)
  769.                IPT5=IPT2.LISOUS(I)
  770.                SEGDES,IPT5
  771.             ENDDO
  772.          ENDIF
  773.          SEGDES,IPT2
  774.       ENDIF
  775. 9999  CONTINUE
  776.       SEGSUP,ICPR
  777.       NBSOUS=LISOUS(/1)
  778.       IF (NBSOUS.GT.0) THEN
  779.          DO I=1,NBSOUS
  780.             IPT1=LISOUS(I)
  781.             SEGDES,IPT1
  782.          ENDDO
  783.       ENDIF
  784.       SEGDES,MELEME
  785.       IF (KAN.GT.0) THEN
  786.          DO K=1,N1
  787.             MELVAL=IELVAX(K)
  788.             SEGSUP,MELVAL
  789.             MELVAL=IELVAY(K)
  790.             SEGSUP,MELVAL
  791.             IF (IDIM.EQ.3) THEN
  792.                MELVAL=IELVAZ(K)
  793.                SEGSUP,MELVAL
  794.             ENDIF
  795.          ENDDO
  796.          SEGSUP,MIELVA
  797.       ENDIF
  798.  
  799.       END
  800.  
  801.  

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