Télécharger la plaquette Cast3M


fr en

Télécharger ijet.eso

Retour à la liste

Numérotation des lignes : 

ijet

  1. C IJET      SOURCE    CB215821  20/11/25    13:30:08     10792          
  2. C IJET     SOURCE    PM        01/12/11    21:44:48     4041
  3.       SUBROUTINE IJET
  4. CA   Ph. MIRANDA
  5. CD   30/11/2001
  6. CC--------------------------------------------------------------------------
  7. CC
  8. CC     Le sous-programme IJET calcule par une methode analytique exacte
  9. CC     l'intersection entre des segments (connus par les coordonnees de
  10. CC     leurs extremites) et les facettes triangulaires d'un maillage
  11. CC     (connues par les coordonnees de leurs 3 extremites et leurs
  12. CC     normales). De plus, ce programme permet de deduire les valeurs
  13. CC     d'un champ appuye sur les noeuds du maillage aux points
  14. CC     d'intersection des segments par interpolation.
  15. CC
  16. CC     IJET reprend les fonctionnalites de la procedure GIBIANE @INTSEC
  17. CC     inspiree de @INTERC concue par R.Mitteau (CEA/DSM/DRFC) et
  18. CC     B.Riou (CSSI) en oct 1998.
  19. CC
  20. CC     version initiale : Ph.Miranda (CSSI)    en nov 2001
  21. CC     modification     : T.Charras  (CEA/DM2S) en oct 2003
  22. CC     modification     : A.Moal     (CSSI)    en jan 2004
  23. CC
  24. CC     Syntaxe generale :
  25. CC
  26. CC     CHDIST MINTER CHFN CHDEP = IJET CHOLD CHNEW TOL TAB1 ;
  27. CC
  28. CC--------------------------------------------------------------------------
  29. CE     En entree :
  30. CE
  31. CE     CHOLD : coordonnees des points extremites initiales des segments
  32. CE            (type CHPOINT a 3 composantes)
  33. CE     CHNEW : coordonnees des points extremites finales des segments
  34. CE            (type CHPOINT a 3 composantes appuye sur le meme support
  35. CE             que CH_OLD)
  36. CE     TOL   : tolerance (type FLOTTANT)
  37. CE     TAB1  : Table qui doit contenir les parametres suivants en entree
  38. CE
  39. CE      <chamx1 (type MCHAML) coordonnee x du premier noeud de l'element
  40. CE      <chamy1 (type MCHAML) coordonnee y du premier noeud de l'element
  41. CE      <chamz1 (type MCHAML) coordonnee z du premier noeud de l'element
  42. CE
  43. CE      <chamx2 (type MCHAML) coordonnee x du deuxieme noeud de l'element
  44. CE      <chamy2 (type MCHAML) coordonnee y du deuxieme noeud de l'element
  45. CE      <chamz2 (type MCHAML) coordonnee z du deuxieme noeud de l'element
  46. CE
  47. CE      <chamx3 (type MCHAML) coordonnee x du troisieme noeud de l'element
  48. CE      <chamy3 (type MCHAML) coordonnee y du troisieme noeud de l'element
  49. CE      <chamz3 (type MCHAML) coordonnee z du troisieme noeud de l'element
  50. CE
  51. CE      <chamf1 (type MCHAML) flux au premier noeud de l'element
  52. CE      <chamf2 (type MCHAML) flux au second  noeud de l'element
  53. CE      <chamf3 (type MCHAML) flux au troisieme noeud de l'element
  54. CE
  55. CE      <cosx (type MCHAML) cosinus directeur selon x
  56. CE      <cosy (type MCHAML) cosinus directeur selon y
  57. CE      <cosz (type MCHAML) cosinus directeur selon z
  58. CC
  59. CS      En sortie :
  60. CS
  61. CS      CHDIST : CHPOINT appuye sur le support de CHOLD et contenant
  62. CS               le pas si pas d'intersection et la distance entre
  63. CS               point initial du segment et le point d'intersection
  64. CS               a la facette sinon (type CHPOINT)
  65. CS      MINTER : contient les noeuds intersectes du support de CHOLD.
  66. CS               (type MAILLAGE)
  67. CS      CHFN   : CHPOINT appuye sur le support de CHOLD et contenant le
  68. CS               le champ de flux interpole si intersection et 0. sinon.
  69. CS               (type CHPOINT a une composante)
  70. CS      CHDEP  : champ de deplacement des points interceptes
  71. CS               (type CHPOINT a 3 composantes)
  72. CS
  73. CC      Remarque : il est necessaire de travailler avec les champs par
  74. CC                 element du maillage interceptant et non pas avec les
  75. CC                 coordonnées des points de l'objet maillage car ces
  76. CC                 champs peuvent avoir subi des changements de repere
  77. CC
  78. CC--------------------------------------------------------------------------
  79.       IMPLICIT INTEGER(I-N)
  80.       IMPLICIT REAL*8 (A-H,O-Z)
  81. C
  82.  
  83. -INC PPARAM
  84. -INC CCOPTIO
  85. -INC SMCHPOI
  86. -INC SMCHAML
  87. -INC SMELEME
  88. -INC SMCOORD
  89. -INC CCGEOME
  90. -INC CCREEL
  91. -INC SMTABLE
  92. C
  93.       SEGMENT MTRAV1
  94.         REAL*8 XPN(NZ1,NBP),YPN(NZ1,NBP),ZPN(NZ1,NBP)
  95.         REAL*8 XPN1(NZ1,NBP),YPN1(NZ1,NBP),ZPN1(NZ1,NBP)
  96.         REAL*8 CHPAS(NZ1,NBP),CHPAS0(NZ1,NBP),CHFN(NZ1,NBP)
  97.         REAL*8 CHDEX(NZ1,NBP),CHDEY(NZ1,NBP),CHDEZ(NZ1,NBP)
  98.         LOGICAL INTER(NZ1,NBP)
  99.       ENDSEGMENT
  100. C
  101.       SEGMENT MTRAV2
  102.         REAL*8 COORFAC(NZ2,NEL,3,3),TFN(NZ2,NEL,3)
  103.         REAL*8 TXN(NZ2,NEL),TYN(NZ2,NEL),TZN(NZ2,NEL)
  104.       ENDSEGMENT
  105. C
  106.       POINTEUR CHOLD.MCHPOI,CHNEW.MCHPOI
  107.       POINTEUR CHX1.MCHPOI,CHX2.MCHPOI,CHX3.MCHPOI
  108.       POINTEUR CHY1.MCHPOI,CHY2.MCHPOI,CHY3.MCHPOI
  109.       POINTEUR CHZ1.MCHPOI,CHZ2.MCHPOI,CHZ3.MCHPOI
  110.       POINTEUR TAB1.MTABLE
  111.       POINTEUR ICHOLD.MELEME
  112.       POINTEUR CHCOSX.MCHAML,CHCOSY.MCHAML,CHCOSZ.MCHAML
  113.       POINTEUR CHF1.MCHAML,CHF2.MCHAML,CHF3.MCHAML
  114. C
  115.       CHARACTER*8 MTYPR
  116.       CHARACTER*4 REP,CHARRE
  117.       CHARACTER*9 MOTCLE(12)
  118.       logical login,logre
  119. C
  120. C  LECTURE DES OBJETS
  121. C
  122.       CALL LIROBJ('CHPOINT',CHOLD,1,IRETOU)
  123.       IF (IERR .NE. 0) RETURN
  124.       CALL LIROBJ('CHPOINT',CHNEW,1,IRETOU)
  125.       IF (IERR .NE. 0) RETURN
  126.       CALL LIRREE(TOL1,1,IRETOU)
  127.       IF (IERR .NE. 0) RETURN
  128.       CALL LIROBJ('TABLE   ',TAB1,1,IRETOU)
  129.       IF (IERR .NE. 0) RETURN
  130. C
  131. C  Dimensionnement des tableaux du segment MTRAV1
  132. C
  133.       MCHPOI=CHOLD
  134.       SEGACT MCHPOI
  135.       NZ1=IPCHP(/1)
  136.       DO 50 I=1,IPCHP(/1)
  137.          MSOUPO=IPCHP(I)
  138.          SEGACT MSOUPO
  139.          MELEME=IGEOC
  140.          SEGACT MELEME
  141.          IF (I .EQ. 1) NBP=NUM(/2)
  142.          IF (NUM(/2) .GT. NBP) NBP=NUM(/2)
  143.          SEGDES MELEME
  144.          SEGDES MSOUPO
  145. 50    CONTINUE
  146.       SEGINI MTRAV1
  147. C
  148. C
  149. C   EXTRACTION DES COORDONNEES DES POINTS Pn et Pn+1 a partir des CHPOINT
  150. C       CHOLD et CHNEW
  151. C
  152.       DO 100 I=1,IPCHP(/1)
  153.          MSOUPO=IPCHP(I)
  154.          SEGACT MSOUPO
  155.          ICHOLD=IGEOC
  156.          MPOVAL=IPOVAL
  157.          SEGACT MPOVAL
  158.          DO 200 J=1,NOCOMP(/2)
  159.              IF (NOCOMP(J).EQ. 'X   ')then
  160.               do 331 k=1,vpocha(/1)
  161.                 XPN(I,K)=VPOCHA(K,J)
  162.   331         continue
  163.              elseif(NOCOMP(J).EQ. 'Y   ') then
  164.               do 332 k=1,vpocha(/1)
  165.                 YPN(I,K)=VPOCHA(K,J)
  166.   332         continue
  167.              elseif(NOCOMP(J).EQ. 'Z   ')then
  168.               do 333 k=1,vpocha(/1)
  169.                 ZPN(I,K)=VPOCHA(K,J)
  170.   333         continue
  171.              endif
  172. 200      CONTINUE
  173.          SEGDES MPOVAL
  174.          SEGDES MELEME
  175.          SEGDES MSOUPO
  176. 100   CONTINUE
  177.       SEGDES MCHPOI
  178. C
  179.       MCHPOI=CHNEW
  180.       SEGACT MCHPOI
  181.       DO 110 I=1,IPCHP(/1)
  182.          MSOUPO=IPCHP(I)
  183.          SEGACT MSOUPO
  184.          MPOVAL=IPOVAL
  185.          SEGACT MPOVAL
  186.          DO 210 J=1,NOCOMP(/2)
  187.                IF (NOCOMP(J).EQ. 'X   ')then
  188.                    DO 310 K=1,VPOCHA(/1)
  189.                     XPN1(I,K)=VPOCHA(K,J)
  190.   310              CONTINUE
  191.                ELSEIF (NOCOMP(J).EQ. 'Y   ')then
  192.                    DO 311 K=1,VPOCHA(/1)
  193.                     YPN1(I,K)=VPOCHA(K,J)
  194.   311              CONTINUE
  195.                ELSEIF (NOCOMP(J).EQ. 'Z   ')then
  196.                    DO 312 K=1,VPOCHA(/1)
  197.                     ZPN1(I,K)=VPOCHA(K,J)
  198.   312              CONTINUE
  199.                ENDIF
  200. 210      CONTINUE
  201.          SEGDES MPOVAL
  202.          SEGDES MSOUPO
  203. 110   CONTINUE
  204.       SEGDES MCHPOI
  205. C
  206.          MOTCLE(1) = '<COSX'
  207.          MOTCLE(2) = '<COSY'
  208.          MOTCLE(3) = '<COSZ'
  209.          MOTCLE(4) = '<CHAMX1  '
  210.          MOTCLE(5) = '<CHAMY1  '
  211.          MOTCLE(6) = '<CHAMZ1  '
  212.          MOTCLE(7) = '<CHAMX2  '
  213.          MOTCLE(8) = '<CHAMY2  '
  214.          MOTCLE(9) = '<CHAMZ2  '
  215.          MOTCLE(10) = '<CHAMX3  '
  216.          MOTCLE(11) = '<CHAMY3  '
  217.          MOTCLE(12) = '<CHAMZ3  '
  218. C
  219.       SEGACT TAB1
  220.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(1),LOGIN,IOBIN,
  221.      .                 'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHCOSX)
  222.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(2),LOGIN,IOBIN,
  223.      .                 'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHCOSY)
  224.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(3),LOGIN,IOBIN,
  225.      .                 'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHCOSZ)
  226. C
  227.       MCHELM=CHCOSX
  228.       SEGACT MCHELM
  229. C
  230. C  Dimensionnement des tableaux du segment MTRAV2
  231. C
  232.       NZ2 = ICHAML(/1)
  233.       DO 60 I=1,ICHAML(/1)
  234.          MCHAML=ICHAML(I)
  235.          SEGACT MCHAML
  236.          DO 70 IC=1,IELVAL(/1)
  237.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  238.                 MELVAL=IELVAL(IC)
  239.                 SEGACT MELVAL
  240.                 IF (I .EQ. 1) NEL=VELCHE(/2)
  241.                 IF (VELCHE(/2) .GT. NEL) NEL=VELCHE(/2)
  242.                 SEGDES MELVAL
  243.              ENDIF
  244. 70      CONTINUE
  245.         SEGDES MCHAML
  246. 60    CONTINUE
  247.       SEGINI MTRAV2
  248. C
  249. C  Extraction des trois cosinus directeurs (CHAMELEM) a partir de TAB1
  250. C  Les elements du maillage associe sont ici des POI1
  251. C
  252.       DO 120 I=1,ICHAML(/1)
  253.          MCHAML=ICHAML(I)
  254.          SEGACT MCHAML
  255.          DO 220 IC=1,IELVAL(/1)
  256.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  257.                 MELVAL=IELVAL(IC)
  258.                 SEGACT MELVAL
  259.                 DO 320 IEL=1,VELCHE(/2)
  260.                    TXN(I,IEL) = VELCHE(1,IEL)
  261. 320             CONTINUE
  262.                 SEGDES MELVAL
  263.              ENDIF
  264. 220      CONTINUE
  265.          SEGDES MCHAML
  266. 120   CONTINUE
  267.       SEGDES MCHELM
  268. C
  269.       MCHELM=CHCOSY
  270.       SEGACT MCHELM
  271.       DO 130 I=1,ICHAML(/1)
  272.          MCHAML=ICHAML(I)
  273.          SEGACT MCHAML
  274.          DO 230 IC=1,IELVAL(/1)
  275.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  276.                 MELVAL=IELVAL(IC)
  277.                 SEGACT MELVAL
  278.                 DO 330 IEL=1,VELCHE(/2)
  279.                    TYN(I,IEL) = VELCHE(1,IEL)
  280. 330             CONTINUE
  281.                 SEGDES MELVAL
  282.              ENDIF
  283. 230      CONTINUE
  284.          SEGDES MCHAML
  285. 130   CONTINUE
  286.       SEGDES MCHELM
  287. C
  288.       MCHELM=CHCOSZ
  289.       SEGACT MCHELM
  290.       DO 140 I=1,ICHAML(/1)
  291.          MCHAML=ICHAML(I)
  292.          SEGACT MCHAML
  293.          DO 240 IC=1,IELVAL(/1)
  294.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  295.                 MELVAL=IELVAL(IC)
  296.                 SEGACT MELVAL
  297.                 DO 340 IEL=1,VELCHE(/2)
  298.                    TZN(I,IEL) = VELCHE(1,IEL)
  299. 340             CONTINUE
  300.                 SEGDES MELVAL
  301.              ENDIF
  302. 240      CONTINUE
  303.          SEGDES MCHAML
  304. 140   CONTINUE
  305.       SEGDES MCHELM
  306. C
  307. C   Extraction des coordonnees des 3 points de chaque facette
  308. C
  309.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(4),LOGIN,
  310.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHX1)
  311.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(5),LOGIN,
  312.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHY1)
  313.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(6),LOGIN,
  314.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHZ1)
  315. C
  316.       MCHELM=CHX1
  317.       SEGACT MCHELM
  318.       DO 141 I=1,ICHAML(/1)
  319.          MCHAML=ICHAML(I)
  320.          SEGACT MCHAML
  321.          DO 241 IC=1,IELVAL(/1)
  322.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  323.                 MELVAL=IELVAL(IC)
  324.                 SEGACT MELVAL
  325.                 DO 341 IEL=1,VELCHE(/2)
  326.                    COORFAC(I,IEL,1,1) = VELCHE(1,IEL)
  327. 341             CONTINUE
  328.                 SEGDES MELVAL
  329.              ENDIF
  330. 241      CONTINUE
  331.          SEGDES MCHAML
  332. 141   CONTINUE
  333.       SEGDES MCHELM
  334. C
  335.       MCHELM=CHY1
  336.       SEGACT MCHELM
  337.       DO 142 I=1,ICHAML(/1)
  338.          MCHAML=ICHAML(I)
  339.          SEGACT MCHAML
  340.          DO 242 IC=1,IELVAL(/1)
  341.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  342.                 MELVAL=IELVAL(IC)
  343.                 SEGACT MELVAL
  344.                 DO 342 IEL=1,VELCHE(/2)
  345.                    COORFAC(I,IEL,1,2) = VELCHE(1,IEL)
  346. 342             CONTINUE
  347.                 SEGDES MELVAL
  348.              ENDIF
  349. 242      CONTINUE
  350.          SEGDES MCHAML
  351. 142   CONTINUE
  352.       SEGDES MCHELM
  353. C
  354.       MCHELM=CHZ1
  355.       SEGACT MCHELM
  356.       DO 143 I=1,ICHAML(/1)
  357.          MCHAML=ICHAML(I)
  358.          SEGACT MCHAML
  359.          DO 243 IC=1,IELVAL(/1)
  360.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  361.                 MELVAL=IELVAL(IC)
  362.                 SEGACT MELVAL
  363.                 DO 343 IEL=1,VELCHE(/2)
  364.                    COORFAC(I,IEL,1,3) = VELCHE(1,IEL)
  365. 343             CONTINUE
  366.                 SEGDES MELVAL
  367.              ENDIF
  368. 243      CONTINUE
  369.          SEGDES MCHAML
  370. 143   CONTINUE
  371.       SEGDES MCHELM
  372. C
  373.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(7),LOGIN,
  374.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHX2)
  375.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(8),LOGIN,
  376.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHY2)
  377.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(9),LOGIN,
  378.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHZ2)
  379. C
  380.       MCHELM=CHX2
  381.       SEGACT MCHELM
  382.       DO 144 I=1,ICHAML(/1)
  383.          MCHAML=ICHAML(I)
  384.          SEGACT MCHAML
  385.          DO 244 IC=1,IELVAL(/1)
  386.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  387.                 MELVAL=IELVAL(IC)
  388.                 SEGACT MELVAL
  389.                 DO 344 IEL=1,VELCHE(/2)
  390.                     COORFAC(I,IEL,2,1) = VELCHE(1,IEL)
  391. 344             CONTINUE
  392.                 SEGDES MELVAL
  393.              ENDIF
  394. 244      CONTINUE
  395.          SEGDES MCHAML
  396. 144   CONTINUE
  397.       SEGDES MCHELM
  398. C
  399.       MCHELM=CHY2
  400.       SEGACT MCHELM
  401.       DO 145 I=1,ICHAML(/1)
  402.          MCHAML=ICHAML(I)
  403.          SEGACT MCHAML
  404.          DO 245 IC=1,IELVAL(/1)
  405.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  406.                 MELVAL=IELVAL(IC)
  407.                 SEGACT MELVAL
  408.                 DO 345 IEL=1,VELCHE(/2)
  409.                    COORFAC(I,IEL,2,2) = VELCHE(1,IEL)
  410. 345             CONTINUE
  411.                 SEGDES MELVAL
  412.              ENDIF
  413. 245      CONTINUE
  414.          SEGDES MCHAML
  415. 145   CONTINUE
  416.       SEGDES MCHELM
  417. C
  418.       MCHELM=CHZ2
  419.       SEGACT MCHELM
  420.       DO 146 I=1,ICHAML(/1)
  421.          MCHAML=ICHAML(I)
  422.          SEGACT MCHAML
  423.          DO 246 IC=1,IELVAL(/1)
  424.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  425.                 MELVAL=IELVAL(IC)
  426.                 SEGACT MELVAL
  427.                 DO 346 IEL=1,VELCHE(/2)
  428.                    COORFAC(I,IEL,2,3) = VELCHE(1,IEL)
  429. 346             CONTINUE
  430.                 SEGDES MELVAL
  431.              ENDIF
  432. 246      CONTINUE
  433.          SEGDES MCHAML
  434. 146   CONTINUE
  435.       SEGDES MCHELM
  436. C
  437.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(10),LOGIN,
  438.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHX3)
  439.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(11),LOGIN,
  440.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHY3)
  441.       CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,MOTCLE(12),LOGIN,
  442.      .           IOBIN,'MCHAML  ',IVALRE,XVALRE,CHARRE,LOGRE,CHZ3)
  443. C
  444.       MCHELM=CHX3
  445.       SEGACT MCHELM
  446.       DO 147 I=1,ICHAML(/1)
  447.          MCHAML=ICHAML(I)
  448.          SEGACT MCHAML
  449.          DO 247 IC=1,IELVAL(/1)
  450.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  451.                 MELVAL=IELVAL(IC)
  452.                 SEGACT MELVAL
  453.                 DO 347 IEL=1,VELCHE(/2)
  454.                    COORFAC(I,IEL,3,1) = VELCHE(1,IEL)
  455. 347             CONTINUE
  456.                 SEGDES MELVAL
  457.              ENDIF
  458. 247      CONTINUE
  459.          SEGDES MCHAML
  460. 147   CONTINUE
  461.       SEGDES MCHELM
  462. C
  463.       MCHELM=CHY3
  464.       SEGACT MCHELM
  465.       DO 148 I=1,ICHAML(/1)
  466.          MCHAML=ICHAML(I)
  467.          SEGACT MCHAML
  468.          DO 248 IC=1,IELVAL(/1)
  469.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  470.                 MELVAL=IELVAL(IC)
  471.                 SEGACT MELVAL
  472.                 DO 348 IEL=1,VELCHE(/2)
  473.                    COORFAC(I,IEL,3,2) = VELCHE(1,IEL)
  474. 348             CONTINUE
  475.                 SEGDES MELVAL
  476.              ENDIF
  477. 248      CONTINUE
  478.          SEGDES MCHAML
  479. 148   CONTINUE
  480.       SEGDES MCHELM
  481. C
  482.       MCHELM=CHZ3
  483.       SEGACT MCHELM
  484.       DO 149 I=1,ICHAML(/1)
  485.          MCHAML=ICHAML(I)
  486.          SEGACT MCHAML
  487.          DO 249 IC=1,IELVAL(/1)
  488.              IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  489.                 MELVAL=IELVAL(IC)
  490.                 SEGACT MELVAL
  491.                 DO 349 IEL=1,VELCHE(/2)
  492.                    COORFAC(I,IEL,3,3) = VELCHE(1,IEL)
  493. 349             CONTINUE
  494.                 SEGDES MELVAL
  495.              ENDIF
  496. 249      CONTINUE
  497.          SEGDES MCHAML
  498. 149   CONTINUE
  499.       SEGDES MCHELM
  500. C
  501. C     Lecture des flux
  502. C
  503.          MTYPR = '        '
  504.        CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,'<CHAMF1',LOGIN,
  505.      .              IOBIN,MTYPR,IVALRE,XVALRE,CHARRE,LOGRE,CHF1)
  506.          IF (MTYPR .EQ. 'MCHAML  ') THEN
  507.             MCHELM=CHF1
  508.             SEGACT MCHELM
  509.             DO 150 I=1,ICHAML(/1)
  510.                MCHAML=ICHAML(I)
  511.                SEGACT MCHAML
  512.                DO 250 IC=1,IELVAL(/1)
  513.                   IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  514.                      MELVAL=IELVAL(IC)
  515.                      SEGACT MELVAL
  516.                      DO 350 IEL=1,VELCHE(/2)
  517.                         TFN(I,IEL,1) = VELCHE(1,IEL)
  518. 350                  CONTINUE
  519.                      SEGDES MELVAL
  520.                   ENDIF
  521. 250            CONTINUE
  522.                SEGDES MCHAML
  523. 150         CONTINUE
  524.             SEGDES MCHELM
  525.          ENDIF
  526. C
  527.          MTYPR = '        '
  528.        CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,'<CHAMF2',LOGIN,
  529.      .           IOBIN,MTYPR,IVALRE,XVALRE,CHARRE,LOGRE,CHF2)
  530.          IF (MTYPR .EQ. 'MCHAML  ') THEN
  531.            MCHELM=CHF2
  532.            SEGACT MCHELM
  533.            DO 151 I=1,ICHAML(/1)
  534.               MCHAML=ICHAML(I)
  535.               SEGACT MCHAML
  536.               DO 251 IC=1,IELVAL(/1)
  537.                  IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  538.                     MELVAL=IELVAL(IC)
  539.                     SEGACT MELVAL
  540.                     DO 351 IEL=1,VELCHE(/2)
  541.                        TFN(I,IEL,2) = VELCHE(1,IEL)
  542. 351                 CONTINUE
  543.                     SEGDES MELVAL
  544.                  ENDIF
  545. 251           CONTINUE
  546.               SEGDES MCHAML
  547. 151        CONTINUE
  548.            SEGDES MCHELM
  549.          ENDIF
  550. C
  551.          MTYPR = '        '
  552.        CALL ACCTAB(TAB1,'MOT     ',IVALIN,XVALIN,'<CHAMF3',LOGIN,
  553.      .           IOBIN,MTYPR,IVALRE,XVALRE,CHARRE,LOGRE,CHF3)
  554. C
  555.          IF (MTYPR .EQ. 'MCHAML  ') THEN
  556.             MCHELM=CHF3
  557.             SEGACT MCHELM
  558.             DO 152 I=1,ICHAML(/1)
  559.                MCHAML=ICHAML(I)
  560.                SEGACT MCHAML
  561.                DO 252 IC=1,IELVAL(/1)
  562.                   IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  563.                      MELVAL=IELVAL(IC)
  564.                      SEGACT MELVAL
  565.                      DO 352 IEL=1,VELCHE(/2)
  566.                         TFN(I,IEL,3) = VELCHE(1,IEL)
  567. 352                  CONTINUE
  568.                      SEGDES MELVAL
  569.                   ENDIF
  570. 252            CONTINUE
  571.                SEGDES MCHAML
  572. 152         CONTINUE
  573.             SEGDES MCHELM
  574.          ENDIF
  575. C
  576.       SEGDES TAB1
  577. C
  578. C  Initialisation du tableau des indicateurs d'interception
  579. C
  580.       DO 302 I=1,NZ1
  581.          DO 301 K=1,NBP
  582.             INTER(I,K)=.FALSE.
  583. 301      CONTINUE
  584. 302   CONTINUE
  585. C
  586. C   Boucle sur les facettes du maillage OMBRANT
  587. C
  588.       MCHELM=CHX1
  589.       SEGACT MCHELM
  590.       DO 400 I=1,ICHAML(/1)
  591.          MCHAML=ICHAML(I)
  592.          SEGACT MCHAML
  593.          DO 450 IC=1,IELVAL(/1)
  594.             IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  595.                MELVAL=IELVAL(IC)
  596.                SEGACT MELVAL
  597.                DO 500 J=1,VELCHE(/2)
  598.                   XA1 =  COORFAC(I,J,1,1)
  599.                   YA1 =  COORFAC(I,J,1,2)
  600.                   ZA1 =  COORFAC(I,J,1,3)
  601. C
  602.                   XB1 =  COORFAC(I,J,2,1)
  603.                   YB1 =  COORFAC(I,J,2,2)
  604.                   ZB1 =  COORFAC(I,J,2,3)
  605. C
  606.                   XC1 =  COORFAC(I,J,3,1)
  607.                   YC1 =  COORFAC(I,J,3,2)
  608.                   ZC1 =  COORFAC(I,J,3,3)
  609. C
  610. C  Boucle sur les sous-zones du maillage OMBRE
  611. C
  612.                   DO 550 IZ=1,NZ1
  613. C
  614. C  Boucle sur les lignes de champ (ou sur les elements POI1 de OMBRE)
  615. C
  616.                      DO 600 K=1,NBP
  617.                         IF (.NOT. INTER(IZ,K)) THEN
  618. C
  619. C  Produits scalaires
  620. C
  621.                            PS1 = ((XPN(IZ,K)-XA1)*TXN(I,J)) +
  622.      .                           ((YPN(IZ,K)-YA1)*TYN(I,J)) +
  623.      .                           ((ZPN(IZ,K)-ZA1)*TZN(I,J))
  624.                            PS2 = ((XPN1(IZ,K)-XA1)*TXN(I,J)) +
  625.      .                           ((YPN1(IZ,K)-YA1)*TYN(I,J)) +
  626.      .                           ((ZPN1(IZ,K)-ZA1)*TZN(I,J))
  627. C             -> -->
  628. C   Calcul de PnPn+1 :
  629. C
  630.                            DEX = XPN1(IZ,K) - XPN(IZ,K)
  631.                            DEY = YPN1(IZ,K) - YPN(IZ,K)
  632.                            DEZ = ZPN1(IZ,K) - ZPN(IZ,K)
  633. C
  634. C   calcul de la distance de connection par defaut
  635. C
  636.              CHPAS0(IZ,K) = (DEX*DEX + DEY*DEY + DEZ*DEZ)**.5
  637. C
  638. C   Test d'intersection avec le plan de la facette
  639. C
  640.                   IF ( (PS1*PS2 .LT. (-1*TOL1*TOL1))
  641.      .             .OR. (PS2 .EQ. 0.) )               THEN
  642. C
  643. C                   -> -->   ->
  644. C   Calcul de AD1 = PnPn+1 . N  :
  645. C
  646.                 AD1 = DEX*TXN(I,J) + DEY*TYN(I,J) + DEZ*TZN(I,J)
  647. C
  648. C   Test de detection des lignes de champ paralleles a la facette
  649. C
  650.                            IF (ABS(AD1) .GT. TOL1) THEN
  651. C
  652. C
  653. C    Calcul des points d'intersection M entre les lignes de champ et
  654. C                                             le plan de la facette
  655. C
  656.                            XM = ( (XA1*DEX*TXN(I,J)) -
  657.      .                          ((((YPN(IZ,K)-YA1)*DEX)-
  658.      .                             (XPN(IZ,K)*DEY))*TYN(I,J)) -
  659.      .                          ((((ZPN(IZ,K)-ZA1)*DEX)-
  660.      .                       (XPN(IZ,K)*DEZ))*TZN(I,J)) ) / AD1
  661. C
  662.                            YM = ( (YA1*DEY*TYN(I,J)) -
  663.      .                          ((((ZPN(IZ,K)-ZA1)*DEY)-
  664.      .                             (YPN(IZ,K)*DEZ))*TZN(I,J)) -
  665.      .                          ((((XPN(IZ,K)-XA1)*DEY)-
  666.      .                       (YPN(IZ,K)*DEX))*TXN(I,J)) ) / AD1
  667. C
  668.                            ZM = ( (ZA1*DEZ*TZN(I,J)) -
  669.      .                          ((((XPN(IZ,K)-XA1)*DEZ)-
  670.      .                             (ZPN(IZ,K)*DEX))*TXN(I,J)) -
  671.      .                          ((((YPN(IZ,K)-YA1)*DEZ)-
  672.      .                       (ZPN(IZ,K)*DEY))*TYN(I,J)) ) / AD1
  673. C
  674. C   Calcul des vecteurs AM, BM, CM
  675. C
  676.                            XAM = XM - XA1
  677.                            XBM = XM - XB1
  678.                            XCM = XM - XC1
  679. C
  680. C
  681.                            YAM = YM - YA1
  682.                            YBM = YM - YB1
  683.                            YCM = YM - YC1
  684. C
  685.                            ZAM = ZM - ZA1
  686.                            ZBM = ZM - ZB1
  687.                            ZCM = ZM - ZC1
  688. C
  689. C   Denominateurs
  690. C
  691.                            DX = (YCM*ZBM)-(ZCM*YBM)+(ZAM*YBM)-
  692.      .                          (ZAM*YCM)-(YAM*ZBM)+(YAM*ZCM)
  693. C
  694.                            DY = (ZCM*XBM)-(XCM*ZBM)+(XAM*ZBM)-
  695.      .                          (XAM*ZCM)-(ZAM*XBM)+(ZAM*XCM)
  696. C
  697.                            DZ = (YCM*XBM)-(XCM*YBM)+(XAM*YBM)-
  698.      .                          (XAM*YCM)-(YAM*XBM)+(YAM*XCM)
  699. C
  700.                            IF (ABS(DZ) .GT. TOL1) THEN
  701.                                D1 = DZ
  702.                                A1 = YCM*XBM-XCM*YBM
  703.                                B1 = YAM*XCM-XAM*YCM
  704.                                C1 = YBM*XAM-XBM*YAM
  705.                            ELSEIF (ABS(DY) .GT. TOL1) THEN
  706.                                D1 = DY
  707.                                A1 = ZCM*XBM-XCM*ZBM
  708.                                B1 = ZAM*XCM-XAM*ZCM
  709.                                C1 = ZBM*XAM-XBM*ZAM
  710.                            ELSE
  711.                                D1 = DX
  712.                                A1 = YCM*ZBM-ZCM*YBM
  713.                                B1 = YAM*ZCM-ZAM*YCM
  714.                                C1 = YBM*ZAM-ZBM*YAM
  715.                            ENDIF
  716. C
  717. C    Coordonnees barycentriques
  718. C
  719.                            ALPHA = A1 / D1
  720.                            BETA =  B1 / D1
  721.                            GAMMA = C1 / D1
  722. C
  723.                            IF ( (ALPHA .GE. 0.) .AND.
  724.      .                          (BETA  .GE. 0.) .AND.
  725.      .                          (GAMMA .GE. 0.)      ) THEN
  726. C
  727.                                     INTER(IZ,K)=.TRUE.
  728. C
  729. C
  730. C   Calcul du Flux Normalise
  731. C
  732.                         CHFN(IZ,K) =  TFN(I,J,1)*ALPHA +
  733.      .                                TFN(I,J,2)*BETA  +
  734.      .                                TFN(I,J,3)*GAMMA
  735. C
  736. C   Calcul des champs de deplacement :
  737. C
  738.                              CHDEX(IZ,K) = XPN(IZ,K)-XM
  739.                              CHDEY(IZ,K) = YPN(IZ,K)-YM
  740.                              CHDEZ(IZ,K) = ZPN(IZ,K)-ZM
  741.  
  742. C
  743. C   Calcul de d(Pn,M) :
  744. C
  745.                          CHPAS(IZ,K) = ((XPN(IZ,K)-XM)**2. +
  746.      .          (YPN(IZ,K)-YM)**2. + (ZPN(IZ,K)-ZM)**2.)**.5
  747.                                  ENDIF
  748.                               ENDIF
  749.                            ENDIF
  750.                         ENDIF
  751. 600                  CONTINUE
  752. 550               CONTINUE
  753. 500            CONTINUE
  754.                SEGDES MELVAL
  755.             ENDIF
  756. 450      CONTINUE
  757.          SEGDES MCHAML
  758. 400   CONTINUE
  759.       SEGDES MCHELM
  760. C
  761. C--------------------------------------------------------
  762. C
  763. C  Creation du maillage MINTER des points d'intersection
  764. C
  765.       NBNN = 1
  766.       NBELEM = NBP
  767.       NBSOUS = 0
  768.       NBREF = 0
  769.       SEGINI MELEME
  770.       MINTER  = MELEME
  771.       ITYPEL = 1
  772. C
  773. C   Construction de MINTER
  774. C
  775.       IPT1 = ICHOLD
  776.       SEGACT IPT1
  777.       NINTER = 0
  778.       DO 701 I=1,NZ1
  779.          DO 700 K=1,NBP
  780.             IF (INTER(I,K)) THEN
  781.                NINTER = NINTER+1
  782.                NUM(1,NINTER) = IPT1.NUM(1,K)
  783.             ENDIF
  784. 700      CONTINUE
  785. 701   CONTINUE
  786.       NBELEM = NINTER
  787.       SEGADJ MELEME
  788.       SEGDES MELEME,IPT1
  789. C--------------------------------------------------------
  790. C
  791. C  Creation du CHPOINT CHDIST
  792. C
  793. C segment MCHPOI
  794.       NAT = 1
  795.       NSOUPO = NZ1
  796.       SEGINI MCHPOI
  797.       ICHDIS = MCHPOI
  798.       IFOPOI = 2
  799.       JATTRI(1) = 2
  800. C segment MSOUPO
  801.       NC = 1
  802.       DO 901 I=1,NZ1
  803.          SEGINI MSOUPO
  804.          IPCHP(I)=MSOUPO
  805.          NOCOMP(1)= 'SCAL'
  806. C segment MELEME : maillage associe a CHDIST
  807.          NBNN = 1
  808.          NBELEM = NBP
  809.          NBSOUS = 0
  810.          NBREF = 0
  811. C
  812. C   LE CHAMPOINT CHDIST S'APPUIE SUR LE MEME MAILLAGE QUE CHOLD
  813. C
  814.          IGEOC=ICHOLD
  815. C
  816. C   Construction de MPOVAL
  817. C
  818.          N = NBP
  819.          SEGINI MPOVAL
  820.          IPOVAL=MPOVAL
  821. C
  822. C   On met la valeur de la distance de connection dans le chpoint CHDIST
  823. C       Ce point n'existera plus dans le prochain maillage MINTER
  824. C
  825.          DO 900 K=1,NBP
  826.             IF (INTER(I,K)) THEN
  827.                 VPOCHA(K,1)=CHPAS(I,K)
  828.             ELSE
  829.                 VPOCHA(K,1)=CHPAS0(I,K)
  830.             ENDIF
  831. 900      CONTINUE
  832. C
  833.          SEGDES MPOVAL
  834.          SEGDES MSOUPO
  835. 901   CONTINUE
  836.       SEGDES MCHPOI
  837. C
  838. C--------------------------------------------------------
  839. C
  840. C  Creation du CHPOINT CHFN
  841. C
  842. C segment MCHPOI
  843.          NAT = 1
  844.          NSOUPO = NZ1
  845.          SEGINI MCHPOI
  846.          ICHFN = MCHPOI
  847.          IFOPOI = 2
  848.          JATTRI(1)=2
  849. C segment MSOUPO
  850.          NC = 1
  851.          DO 911 I=1,NZ1
  852.             SEGINI MSOUPO
  853.             IPCHP(I)=MSOUPO
  854.             NOCOMP(1)= 'SCAL'
  855. C segment MELEME : maillage associe a CHDIST
  856.             NBNN = 1
  857.             NBELEM = NBP
  858.             NBSOUS = 0
  859.             NBREF = 0
  860. C
  861. C   LE CHAMPOINT CHFN S'APPUIE SUR LE MEME MAILLAGE QUE CHOLD
  862. C
  863.             IGEOC=ICHOLD
  864. C
  865. C   Construction de MPOVAL
  866. C
  867.             N = NBP
  868.             SEGINI MPOVAL
  869.             IPOVAL=MPOVAL
  870. C
  871. C   On met la valeur de la distance de connection dans le chpoint CHDIST
  872. C       Ce point n'existera plus dans le prochain maillage MINTER
  873. C
  874.             DO 910 K=1,NBP
  875.                IF (INTER(I,K)) THEN
  876.                    VPOCHA(K,1)=CHFN(I,K)
  877.                ELSE
  878.                    VPOCHA(K,1)=0.
  879.                ENDIF
  880. 910         CONTINUE
  881. C
  882. C
  883. C   Desactivation des segments
  884. C
  885.             SEGDES MPOVAL
  886.             SEGDES MSOUPO
  887. 911      CONTINUE
  888.          SEGDES MCHPOI
  889. C--------------------------------------------------------
  890. C
  891. C  Creation du CHPOINT CHDEP
  892. C
  893. C segment MCHPOI
  894.          NAT = 1
  895.          NSOUPO = 1
  896.          SEGINI MCHPOI
  897.          ICHDEP = MCHPOI
  898.          IFOPOI = 2
  899.          JATTRI(1) = 1
  900. C segment MSOUPO
  901.          NC = 3
  902.          SEGINI MSOUPO
  903.          IPCHP(1)=MSOUPO
  904.          NOCOMP(1)= 'UX  '
  905.          NOCOMP(2)= 'UY  '
  906.          NOCOMP(3)= 'UZ  '
  907. C segment MELEME : maillage associe a CHDIST
  908.          NBNN = 1
  909.          NBELEM = NBP
  910.          NBSOUS = 0
  911.          NBREF = 0
  912. C
  913. C   LE CHAMPOINT CHDEP S'APPUIE SUR LE MEME MAILLAGE QUE MINTER
  914. C
  915.          IGEOC=MINTER
  916. C
  917. C   Construction de MPOVAL
  918. C
  919.          N = NINTER
  920.          SEGINI MPOVAL
  921.          IPOVAL=MPOVAL
  922. C
  923.          NINTER = 0
  924. C
  925. C   On met la valeur du champ de deplacement dans le chpoint CHDEP
  926. C
  927.          DO 921 I=1,NZ1
  928.             DO 920 K=1,NBP
  929.                IF (INTER(I,K)) THEN
  930.                    NINTER = NINTER+1
  931.                    VPOCHA(NINTER,1)= 0. - CHDEX(I,K)
  932.                    VPOCHA(NINTER,2)= 0. - CHDEY(I,K)
  933.                    VPOCHA(NINTER,3)= 0. - CHDEZ(I,K)
  934.                ENDIF
  935. 920         CONTINUE
  936. 921      CONTINUE
  937. C
  938. C
  939. C   Desactivation des segments
  940. C
  941.          SEGDES MPOVAL
  942.          SEGDES MSOUPO
  943.          SEGDES MCHPOI
  944. C--------------------------------------------------------
  945. C
  946. C   Ecriture des nouveaux objets en sortie
  947. C
  948.       CALL ECROBJ('CHPOINT ',ICHDEP)
  949.       CALL ECROBJ('CHPOINT ',ICHFN)
  950. C
  951.       CALL ECROBJ('MAILLAGE',MINTER)
  952. C
  953.       CALL ECROBJ('CHPOINT ',ICHDIS)
  954. C
  955.       SEGSUP MTRAV1
  956.       SEGSUP MTRAV2
  957. C
  958.       RETURN
  959.       END
  960.  
  961.  
  962.  
  963.  
  964.  
  965.  
  966.  
  967.  
  968.  
  969.  
  970.  

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