Télécharger la plaquette Cast3M


fr en

Télécharger kdom2.eso

Retour à la liste

Numérotation des lignes : 

kdom2

  1. C KDOM2     SOURCE    KK2000    14/04/10    21:15:11     8032
  2.       SUBROUTINE KDOM2(MELEMQ)
  3. C
  4. C************************************************************************
  5. C
  6. C PROJET            :  CASTEM 2000
  7. C
  8. C NOM               :  KDOM2
  9. C
  10. C DESCRIPTION       :  Subroutine called by KDOM1
  11. C                      In the case of EULER model, we change the
  12. C                      position of the QUAF points
  13. C
  14. C LANGAGE           :  FORTRAN 77 + ESOPE 2000 (avec estensions CISI)
  15. C
  16. C AUTEUR            :  A. BECCANTINI, DRN/DMT/SEMT/LTMF
  17. C
  18. C************************************************************************
  19. C
  20. C
  21.  
  22.       IMPLICIT INTEGER(I-N)
  23.       IMPLICIT REAL*8(A-H,O-Z)
  24.  
  25.       INTEGER MELEMQ, NBSOUS, ISOUS, NBEL, IEL, NN1, NN2, NN3, I1
  26.      &     ,NN4,NN5,NN6,NN7,NN8,NN9,NN10,NN11,NN12,NN13,NN14,NN15
  27.      &     ,NN16,NN17,NN18,NN19,NN20,NN21,NN22,NN23,NN24,NN25,NN26
  28.      &     ,NN27
  29.       REAL*8 XCEN(3),VOL,X1(3),X2(3),X3(3),X4(3),X5(3),X6(3)
  30.      &     ,VOL1,VOL2,VOL3,VOL4,VOL5,VOL6
  31.       LOGICAL LOSEG3, LOTRI7, LOQUA9, LOTE15, LOPY19, LOPR21, LOCU27
  32. C
  33.  
  34. -INC PPARAM
  35. -INC CCOPTIO
  36. -INC SMELEME
  37. -INC SMCOORD
  38. C
  39. C Elements allowed are SEG3, TRI7, QUA9, TE15, 'PY19', 'PR21', 'CU27'
  40. C ITYPEL                3      7    11     35    36      34      33
  41. C They can be referred just once
  42. C
  43.       LOSEG3=.FALSE.
  44.       LOTRI7=.FALSE.
  45.       LOQUA9=.FALSE.
  46.       LOTE15=.FALSE.
  47.       LOPY19=.FALSE.
  48.       LOPR21=.FALSE.
  49.       LOCU27=.FALSE.
  50. C
  51.       MELEME=MELEMQ
  52.       SEGACT MELEME
  53.       NBSOUS=MELEME.LISOUS(/1)
  54.       IF(NBSOUS .EQ. 0) NBSOUS=1
  55.       DO ISOUS=1,NBSOUS,1
  56.          IF(NBSOUS .NE. 1)THEN
  57.             IPT1=MELEME.LISOUS(ISOUS)
  58.             SEGACT IPT1
  59.          ELSE
  60.             IPT1=MELEME
  61.          ENDIF
  62.          NBEL=IPT1.NUM(/2)
  63. C
  64.          IF(IPT1.ITYPEL .EQ. 3)THEN
  65.             IF(LOSEG3)THEN
  66. C              Elt already referred
  67.                WRITE(IOIMP,*) 'Subroutine kdom2'
  68.                WRITE(IOIMP,*) 'Mesh type not recognized'
  69.                CALL ERREUR(5)
  70.             ENDIF
  71.             LOSEG3=.TRUE.
  72. C
  73. C********** SEG3
  74. C
  75.             DO IEL=1,NBEL,1
  76.                NN1=IPT1.NUM(1,IEL)
  77.                NN2=IPT1.NUM(2,IEL)
  78.                NN3=IPT1.NUM(3,IEL)
  79.                DO I1=1,IDIM
  80.                   XCOOR((NN2-1)*(IDIM+1)+I1)=0.5D0*
  81.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  82.      &                 XCOOR((NN1-1)*(IDIM+1)+I1))
  83.                ENDDO
  84.             ENDDO
  85.          ELSEIF(IPT1.ITYPEL .EQ. 7)THEN
  86. C
  87. C********** TRI7
  88. C
  89.             IF(LOTRI7)THEN
  90. C              Elt already referred
  91.                WRITE(IOIMP,*) 'Subroutine kdom2'
  92.                WRITE(IOIMP,*) 'Mesh type not recognized'
  93.                CALL ERREUR(5)
  94.             ENDIF
  95.             LOTRI7=.TRUE.
  96. C
  97.             DO IEL=1,NBEL,1
  98.                NN1=IPT1.NUM(1,IEL)
  99.                NN2=IPT1.NUM(2,IEL)
  100.                NN3=IPT1.NUM(3,IEL)
  101.                NN4=IPT1.NUM(4,IEL)
  102.                NN5=IPT1.NUM(5,IEL)
  103.                NN6=IPT1.NUM(6,IEL)
  104.                NN7=IPT1.NUM(7,IEL)
  105.                DO I1=1,IDIM
  106.                   XCOOR((NN2-1)*(IDIM+1)+I1)=0.5D0*
  107.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  108.      &                 XCOOR((NN3-1)*(IDIM+1)+I1))
  109.                   XCOOR((NN4-1)*(IDIM+1)+I1)=0.5D0*
  110.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  111.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  112.                   XCOOR((NN6-1)*(IDIM+1)+I1)=0.5D0*
  113.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  114.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  115.                   XCOOR((NN7-1)*(IDIM+1)+I1)=
  116.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  117.      &                 XCOOR((NN3-1)*(IDIM+1)+I1)+
  118.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))/3.0D0
  119.                ENDDO
  120.             ENDDO
  121.          ELSEIF(IPT1.ITYPEL .EQ. 11)THEN
  122. C
  123. C********** QUA9
  124. C
  125.             IF(LOQUA9)THEN
  126. C              Elt already referred
  127.                WRITE(IOIMP,*) 'Subroutine kdom2'
  128.                WRITE(IOIMP,*) 'Mesh type not recognized'
  129.                CALL ERREUR(5)
  130.             ENDIF
  131.             LOQUA9=.TRUE.
  132. C
  133.             DO IEL=1,NBEL,1
  134.                NN1=IPT1.NUM(1,IEL)
  135.                NN2=IPT1.NUM(2,IEL)
  136.                NN3=IPT1.NUM(3,IEL)
  137.                NN4=IPT1.NUM(4,IEL)
  138.                NN5=IPT1.NUM(5,IEL)
  139.                NN6=IPT1.NUM(6,IEL)
  140.                NN7=IPT1.NUM(7,IEL)
  141.                NN8=IPT1.NUM(8,IEL)
  142.                NN9=IPT1.NUM(9,IEL)
  143. C
  144. C**************** Centre de l'elt
  145. C
  146.                CALL KDOM7(NN1,NN3,NN5,NN7,XCEN)
  147.                IF(IERR.NE.0)GOTO 9999
  148.                DO I1=1,IDIM
  149.                   XCOOR((NN2-1)*(IDIM+1)+I1)=0.5D0*
  150.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  151.      &                 XCOOR((NN3-1)*(IDIM+1)+I1))
  152.                   XCOOR((NN4-1)*(IDIM+1)+I1)=0.5D0*
  153.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  154.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  155.                   XCOOR((NN6-1)*(IDIM+1)+I1)=0.5D0*
  156.      &                 (XCOOR((NN7-1)*(IDIM+1)+I1)+
  157.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  158.                   XCOOR((NN8-1)*(IDIM+1)+I1)=0.5D0*
  159.      &                 (XCOOR((NN7-1)*(IDIM+1)+I1)+
  160.      &                 XCOOR((NN1-1)*(IDIM+1)+I1))
  161.                   XCOOR((NN9-1)*(IDIM+1)+I1)=XCEN(I1)
  162.                ENDDO
  163.             ENDDO
  164.          ELSEIF(IPT1.ITYPEL .EQ. 35)THEN
  165. C
  166. C********** TE15
  167. C
  168.             IF(LOTE15)THEN
  169. C              Elt already referred
  170.                WRITE(IOIMP,*) 'Subroutine kdom2'
  171.                WRITE(IOIMP,*) 'Mesh type not recognized'
  172.                CALL ERREUR(5)
  173.             ENDIF
  174.             LOTE15=.TRUE.
  175. C
  176.             DO IEL=1,NBEL,1
  177.                NN1=IPT1.NUM(1,IEL)
  178.                NN2=IPT1.NUM(2,IEL)
  179.                NN3=IPT1.NUM(3,IEL)
  180.                NN4=IPT1.NUM(4,IEL)
  181.                NN5=IPT1.NUM(5,IEL)
  182.                NN6=IPT1.NUM(6,IEL)
  183.                NN7=IPT1.NUM(7,IEL)
  184.                NN8=IPT1.NUM(8,IEL)
  185.                NN9=IPT1.NUM(9,IEL)
  186.                NN10=IPT1.NUM(10,IEL)
  187.                NN11=IPT1.NUM(11,IEL)
  188.                NN12=IPT1.NUM(12,IEL)
  189.                NN13=IPT1.NUM(13,IEL)
  190.                NN14=IPT1.NUM(14,IEL)
  191.                NN15=IPT1.NUM(15,IEL)
  192.                DO I1=1,IDIM
  193.                   XCOOR((NN2-1)*(IDIM+1)+I1)=0.5D0*
  194.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  195.      &                 XCOOR((NN3-1)*(IDIM+1)+I1))
  196.                   XCOOR((NN4-1)*(IDIM+1)+I1)=0.5D0*
  197.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  198.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  199.                   XCOOR((NN6-1)*(IDIM+1)+I1)=0.5D0*
  200.      &                 (XCOOR((NN5-1)*(IDIM+1)+I1)+
  201.      &                 XCOOR((NN1-1)*(IDIM+1)+I1))
  202.                   XCOOR((NN7-1)*(IDIM+1)+I1)=0.5D0*
  203.      &                 (XCOOR((NN10-1)*(IDIM+1)+I1)+
  204.      &                 XCOOR((NN1-1)*(IDIM+1)+I1))
  205.                   XCOOR((NN8-1)*(IDIM+1)+I1)=0.5D0*
  206.      &                 (XCOOR((NN10-1)*(IDIM+1)+I1)+
  207.      &                 XCOOR((NN3-1)*(IDIM+1)+I1))
  208.                   XCOOR((NN9-1)*(IDIM+1)+I1)=0.5D0*
  209.      &                 (XCOOR((NN10-1)*(IDIM+1)+I1)+
  210.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  211.                   XCOOR((NN11-1)*(IDIM+1)+I1)=
  212.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  213.      &                 XCOOR((NN3-1)*(IDIM+1)+I1)+
  214.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))/3.0D0
  215.                   XCOOR((NN12-1)*(IDIM+1)+I1)=
  216.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  217.      &                 XCOOR((NN3-1)*(IDIM+1)+I1)+
  218.      &                 XCOOR((NN10-1)*(IDIM+1)+I1))/3.0D0
  219.                   XCOOR((NN13-1)*(IDIM+1)+I1)=
  220.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  221.      &                 XCOOR((NN5-1)*(IDIM+1)+I1)+
  222.      &                 XCOOR((NN10-1)*(IDIM+1)+I1))/3.0D0
  223.                   XCOOR((NN14-1)*(IDIM+1)+I1)=
  224.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  225.      &                 XCOOR((NN5-1)*(IDIM+1)+I1)+
  226.      &                 XCOOR((NN10-1)*(IDIM+1)+I1))/3.0D0
  227.                   XCOOR((NN15-1)*(IDIM+1)+I1)=
  228.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  229.      &                 XCOOR((NN3-1)*(IDIM+1)+I1)+
  230.      &                 XCOOR((NN5-1)*(IDIM+1)+I1)+
  231.      &                 XCOOR((NN10-1)*(IDIM+1)+I1))/4.0D0
  232.                ENDDO
  233.             ENDDO
  234.          ELSEIF(IPT1.ITYPEL .EQ. 36)THEN
  235. C
  236. C********** PY19
  237. C
  238.             IF(LOPY19)THEN
  239. C              Elt already referred
  240.                WRITE(IOIMP,*) 'Subroutine kdom2'
  241.                WRITE(IOIMP,*) 'Mesh type not recognized'
  242.                CALL ERREUR(5)
  243.             ENDIF
  244.             LOPY19=.TRUE.
  245. C
  246.             DO IEL=1,NBEL,1
  247.                NN1=IPT1.NUM(1,IEL)
  248.                NN2=IPT1.NUM(2,IEL)
  249.                NN3=IPT1.NUM(3,IEL)
  250.                NN4=IPT1.NUM(4,IEL)
  251.                NN5=IPT1.NUM(5,IEL)
  252.                NN6=IPT1.NUM(6,IEL)
  253.                NN7=IPT1.NUM(7,IEL)
  254.                NN8=IPT1.NUM(8,IEL)
  255.                NN9=IPT1.NUM(9,IEL)
  256.                NN10=IPT1.NUM(10,IEL)
  257.                NN11=IPT1.NUM(11,IEL)
  258.                NN12=IPT1.NUM(12,IEL)
  259.                NN13=IPT1.NUM(13,IEL)
  260.                NN14=IPT1.NUM(14,IEL)
  261.                NN15=IPT1.NUM(15,IEL)
  262.                NN16=IPT1.NUM(16,IEL)
  263.                NN17=IPT1.NUM(17,IEL)
  264.                NN18=IPT1.NUM(18,IEL)
  265.                NN19=IPT1.NUM(19,IEL)
  266. C              Computation of the center position
  267.                CALL KDOM7(NN1,NN3,NN5,NN7,XCEN)
  268.                IF(IERR.NE.0)GOTO 9999
  269.                DO I1=1,IDIM
  270.                   XCOOR((NN2-1)*(IDIM+1)+I1)=0.5D0*
  271.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  272.      &                 XCOOR((NN3-1)*(IDIM+1)+I1))
  273.                   XCOOR((NN4-1)*(IDIM+1)+I1)=0.5D0*
  274.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  275.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  276.                   XCOOR((NN6-1)*(IDIM+1)+I1)=0.5D0*
  277.      &                 (XCOOR((NN5-1)*(IDIM+1)+I1)+
  278.      &                 XCOOR((NN7-1)*(IDIM+1)+I1))
  279.                   XCOOR((NN8-1)*(IDIM+1)+I1)=0.5D0*
  280.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  281.      &                 XCOOR((NN7-1)*(IDIM+1)+I1))
  282.                   XCOOR((NN9-1)*(IDIM+1)+I1)=0.5D0*
  283.      &                 (XCOOR((NN13-1)*(IDIM+1)+I1)+
  284.      &                 XCOOR((NN1-1)*(IDIM+1)+I1))
  285.                   XCOOR((NN10-1)*(IDIM+1)+I1)=0.5D0*
  286.      &                 (XCOOR((NN13-1)*(IDIM+1)+I1)+
  287.      &                 XCOOR((NN3-1)*(IDIM+1)+I1))
  288.                   XCOOR((NN11-1)*(IDIM+1)+I1)=0.5D0*
  289.      &                 (XCOOR((NN13-1)*(IDIM+1)+I1)+
  290.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  291.                   XCOOR((NN12-1)*(IDIM+1)+I1)=0.5D0*
  292.      &                 (XCOOR((NN13-1)*(IDIM+1)+I1)+
  293.      &                 XCOOR((NN7-1)*(IDIM+1)+I1))
  294.                   XCOOR((NN14-1)*(IDIM+1)+I1)=XCEN(I1)
  295.                   XCOOR((NN15-1)*(IDIM+1)+I1)=
  296.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  297.      &                 XCOOR((NN3-1)*(IDIM+1)+I1)+
  298.      &                 XCOOR((NN13-1)*(IDIM+1)+I1))/3.0D0
  299.                   XCOOR((NN16-1)*(IDIM+1)+I1)=
  300.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  301.      &                 XCOOR((NN5-1)*(IDIM+1)+I1)+
  302.      &                 XCOOR((NN13-1)*(IDIM+1)+I1))/3.0D0
  303.                   XCOOR((NN17-1)*(IDIM+1)+I1)=
  304.      &                 (XCOOR((NN5-1)*(IDIM+1)+I1)+
  305.      &                 XCOOR((NN7-1)*(IDIM+1)+I1)+
  306.      &                 XCOOR((NN13-1)*(IDIM+1)+I1))/3.0D0
  307.                   XCOOR((NN18-1)*(IDIM+1)+I1)=
  308.      &                (XCOOR((NN1-1)*(IDIM+1)+I1)+
  309.      &                 XCOOR((NN7-1)*(IDIM+1)+I1)+
  310.      &                 XCOOR((NN13-1)*(IDIM+1)+I1))/3.0D0
  311.                ENDDO
  312. C
  313. C************* To compute the center, we divide the pyramid into
  314. C              4 tetrahedra (see kdom3)
  315. C
  316.                CALL KDOM3(NN1,NN3,NN5,NN7,NN14,NN13,XCEN,VOL)
  317.                DO I1=1,3,1
  318.                   XCOOR((NN19-1)*(IDIM+1)+I1)=XCEN(I1)
  319.                ENDDO
  320.             ENDDO
  321.          ELSEIF(IPT1.ITYPEL .EQ. 34)THEN
  322. C
  323. C********** PR21
  324. C
  325.             IF(LOPR21)THEN
  326. C              Elt already referred
  327.                WRITE(IOIMP,*) 'Subroutine kdom2'
  328.                WRITE(IOIMP,*) 'Mesh type not recognized'
  329.                CALL ERREUR(5)
  330.             ENDIF
  331.             LOPR21=.TRUE.
  332. C
  333.             DO IEL=1,NBEL,1
  334.                NN1=IPT1.NUM(1,IEL)
  335.                NN2=IPT1.NUM(2,IEL)
  336.                NN3=IPT1.NUM(3,IEL)
  337.                NN4=IPT1.NUM(4,IEL)
  338.                NN5=IPT1.NUM(5,IEL)
  339.                NN6=IPT1.NUM(6,IEL)
  340.                NN7=IPT1.NUM(7,IEL)
  341.                NN8=IPT1.NUM(8,IEL)
  342.                NN9=IPT1.NUM(9,IEL)
  343.                NN10=IPT1.NUM(10,IEL)
  344.                NN11=IPT1.NUM(11,IEL)
  345.                NN12=IPT1.NUM(12,IEL)
  346.                NN13=IPT1.NUM(13,IEL)
  347.                NN14=IPT1.NUM(14,IEL)
  348.                NN15=IPT1.NUM(15,IEL)
  349.                NN16=IPT1.NUM(16,IEL)
  350.                NN17=IPT1.NUM(17,IEL)
  351.                NN18=IPT1.NUM(18,IEL)
  352.                NN19=IPT1.NUM(19,IEL)
  353.                NN20=IPT1.NUM(20,IEL)
  354.                NN21=IPT1.NUM(21,IEL)
  355.                DO I1=1,IDIM
  356.                   XCOOR((NN2-1)*(IDIM+1)+I1)=0.5D0*
  357.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  358.      &                 XCOOR((NN3-1)*(IDIM+1)+I1))
  359.                   XCOOR((NN4-1)*(IDIM+1)+I1)=0.5D0*
  360.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  361.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  362.                   XCOOR((NN6-1)*(IDIM+1)+I1)=0.5D0*
  363.      &                 (XCOOR((NN5-1)*(IDIM+1)+I1)+
  364.      &                 XCOOR((NN1-1)*(IDIM+1)+I1))
  365.                   XCOOR((NN7-1)*(IDIM+1)+I1)=0.5D0*
  366.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  367.      &                 XCOOR((NN10-1)*(IDIM+1)+I1))
  368.                   XCOOR((NN8-1)*(IDIM+1)+I1)=0.5D0*
  369.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  370.      &                 XCOOR((NN12-1)*(IDIM+1)+I1))
  371.                   XCOOR((NN9-1)*(IDIM+1)+I1)=0.5D0*
  372.      &                 (XCOOR((NN5-1)*(IDIM+1)+I1)+
  373.      &                 XCOOR((NN14-1)*(IDIM+1)+I1))
  374.                   XCOOR((NN11-1)*(IDIM+1)+I1)=0.5D0*
  375.      &                 (XCOOR((NN10-1)*(IDIM+1)+I1)+
  376.      &                 XCOOR((NN12-1)*(IDIM+1)+I1))
  377.                   XCOOR((NN13-1)*(IDIM+1)+I1)=0.5D0*
  378.      &                 (XCOOR((NN12-1)*(IDIM+1)+I1)+
  379.      &                 XCOOR((NN14-1)*(IDIM+1)+I1))
  380.                   XCOOR((NN15-1)*(IDIM+1)+I1)=0.5D0*
  381.      &                 (XCOOR((NN10-1)*(IDIM+1)+I1)+
  382.      &                 XCOOR((NN14-1)*(IDIM+1)+I1))
  383.                   XCOOR((NN19-1)*(IDIM+1)+I1)=
  384.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  385.      &                 XCOOR((NN3-1)*(IDIM+1)+I1)+
  386.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))/3.0D0
  387.                   XCOOR((NN20-1)*(IDIM+1)+I1)=
  388.      &                 (XCOOR((NN10-1)*(IDIM+1)+I1)+
  389.      &                 XCOOR((NN12-1)*(IDIM+1)+I1)+
  390.      &                 XCOOR((NN14-1)*(IDIM+1)+I1))/3.0D0
  391.                ENDDO
  392. C              Computation of the QUA center positions
  393.                CALL KDOM7(NN1,NN10,NN12,NN3,XCEN)
  394.                IF(IERR.NE.0)GOTO 9999
  395.                DO I1=1,IDIM,1
  396.                   XCOOR((NN16-1)*(IDIM+1)+I1)=XCEN(I1)
  397.                ENDDO
  398.                CALL KDOM7(NN5,NN3,NN12,NN14,XCEN)
  399.                IF(IERR.NE.0)GOTO 9999
  400.                DO I1=1,IDIM,1
  401.                   XCOOR((NN17-1)*(IDIM+1)+I1)=XCEN(I1)
  402.                ENDDO
  403.                CALL KDOM7(NN10,NN1,NN5,NN14,XCEN)
  404.                IF(IERR.NE.0)GOTO 9999
  405.                DO I1=1,IDIM,1
  406.                   XCOOR((NN18-1)*(IDIM+1)+I1)=XCEN(I1)
  407.                ENDDO
  408. C
  409. C************* To compute the center, we divide the PR21
  410. C              into 2 tetrahedra and 3 pyramids
  411. C              (1,3,5,PC), (10,14,12,PC), (1,10,12,3,16,PC),
  412. C              (3,12,14,5,17,PC),(1,5,14,10,18,PC)
  413. C                         |
  414. C                         center of 3,12,14,5
  415. C              The results does not depend of the choice of
  416. C              PC.
  417. C
  418. C
  419. C
  420.                DO I1=1,3,1
  421.                   XCOOR((NN21-1)*(IDIM+1)+I1)=
  422.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  423.      &                 XCOOR((NN3-1)*(IDIM+1)+I1)+
  424.      &                 XCOOR((NN5-1)*(IDIM+1)+I1)+
  425.      &                 XCOOR((NN10-1)*(IDIM+1)+I1)+
  426.      &                 XCOOR((NN12-1)*(IDIM+1)+I1)+
  427.      &                 XCOOR((NN14-1)*(IDIM+1)+I1))/6.0D0
  428.  
  429.                ENDDO
  430. C
  431.                CALL KDOM3(NN10,NN12,NN3,NN1,NN16,NN21,X1,VOL1)
  432.                CALL KDOM3(NN3,NN12,NN14,NN5,NN17,NN21,X2,VOL2)
  433.                CALL KDOM3(NN1,NN5,NN14,NN10,NN18,NN21,X3,VOL3)
  434.                CALL KDOM4(NN1,NN3,NN5,NN21,X4,VOL4)
  435.                CALL KDOM4(NN10,NN14,NN12,NN21,X5,VOL5)
  436.                VOL=VOL1+VOL2+VOL3+VOL4+VOL5
  437. C
  438.                DO I1=1,3,1
  439.                   XCOOR((NN21-1)*(IDIM+1)+I1)=(VOL1*X1(I1)+VOL2*X2(I1)
  440.      &                 +VOL3*X3(I1)+VOL4*X4(I1)+VOL5*X5(I1))/VOL
  441.                ENDDO
  442. C
  443. C********** We check whether the center of gravity of the PR21 is
  444. C           inside its domain
  445. C
  446.                CALL KDOM8(NN1,NN2,NN3,NN4,NN5,NN6,NN7,NN8,NN9,NN10,
  447.      &              NN11,NN12,NN13,NN14,NN15,NN16,NN17,NN18,NN19,NN20,
  448.      &              NN21)
  449.                IF(IERR.NE.0)GOTO 9999
  450.             ENDDO
  451.             IF(NBSOUS .NE. 1) SEGDES IPT1
  452.          ELSEIF(IPT1.ITYPEL .EQ. 33)THEN
  453. C
  454. C********** CU27
  455. C
  456.             IF(LOCU27)THEN
  457. C              Elt already referred
  458.                WRITE(IOIMP,*) 'Subroutine kdom2'
  459.                WRITE(IOIMP,*) 'Mesh type not recognized'
  460.                CALL ERREUR(5)
  461.             ENDIF
  462.             LOCU27=.TRUE.
  463. C
  464.             DO IEL=1,NBEL,1
  465.                NN1=IPT1.NUM(1,IEL)
  466.                NN2=IPT1.NUM(2,IEL)
  467.                NN3=IPT1.NUM(3,IEL)
  468.                NN4=IPT1.NUM(4,IEL)
  469.                NN5=IPT1.NUM(5,IEL)
  470.                NN6=IPT1.NUM(6,IEL)
  471.                NN7=IPT1.NUM(7,IEL)
  472.                NN8=IPT1.NUM(8,IEL)
  473.                NN9=IPT1.NUM(9,IEL)
  474.                NN10=IPT1.NUM(10,IEL)
  475.                NN11=IPT1.NUM(11,IEL)
  476.                NN12=IPT1.NUM(12,IEL)
  477.                NN13=IPT1.NUM(13,IEL)
  478.                NN14=IPT1.NUM(14,IEL)
  479.                NN15=IPT1.NUM(15,IEL)
  480.                NN16=IPT1.NUM(16,IEL)
  481.                NN17=IPT1.NUM(17,IEL)
  482.                NN18=IPT1.NUM(18,IEL)
  483.                NN19=IPT1.NUM(19,IEL)
  484.                NN20=IPT1.NUM(20,IEL)
  485.                NN21=IPT1.NUM(21,IEL)
  486.                NN22=IPT1.NUM(22,IEL)
  487.                NN23=IPT1.NUM(23,IEL)
  488.                NN24=IPT1.NUM(24,IEL)
  489.                NN25=IPT1.NUM(25,IEL)
  490.                NN26=IPT1.NUM(26,IEL)
  491.                NN27=IPT1.NUM(27,IEL)
  492.                DO I1=1,IDIM
  493.                   XCOOR((NN2-1)*(IDIM+1)+I1)=0.5D0*
  494.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  495.      &                 XCOOR((NN3-1)*(IDIM+1)+I1))
  496.                   XCOOR((NN4-1)*(IDIM+1)+I1)=0.5D0*
  497.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  498.      &                 XCOOR((NN5-1)*(IDIM+1)+I1))
  499.                   XCOOR((NN6-1)*(IDIM+1)+I1)=0.5D0*
  500.      &                 (XCOOR((NN5-1)*(IDIM+1)+I1)+
  501.      &                 XCOOR((NN7-1)*(IDIM+1)+I1))
  502.                   XCOOR((NN8-1)*(IDIM+1)+I1)=0.5D0*
  503.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  504.      &                 XCOOR((NN7-1)*(IDIM+1)+I1))
  505.                   XCOOR((NN9-1)*(IDIM+1)+I1)=0.5D0*
  506.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  507.      &                 XCOOR((NN13-1)*(IDIM+1)+I1))
  508.                   XCOOR((NN10-1)*(IDIM+1)+I1)=0.5D0*
  509.      &                 (XCOOR((NN3-1)*(IDIM+1)+I1)+
  510.      &                 XCOOR((NN15-1)*(IDIM+1)+I1))
  511.                   XCOOR((NN11-1)*(IDIM+1)+I1)=0.5D0*
  512.      &                 (XCOOR((NN5-1)*(IDIM+1)+I1)+
  513.      &                 XCOOR((NN17-1)*(IDIM+1)+I1))
  514.                   XCOOR((NN12-1)*(IDIM+1)+I1)=0.5D0*
  515.      &                 (XCOOR((NN7-1)*(IDIM+1)+I1)+
  516.      &                 XCOOR((NN19-1)*(IDIM+1)+I1))
  517.                   XCOOR((NN14-1)*(IDIM+1)+I1)=0.5D0*
  518.      &                 (XCOOR((NN13-1)*(IDIM+1)+I1)+
  519.      &                 XCOOR((NN15-1)*(IDIM+1)+I1))
  520.                   XCOOR((NN16-1)*(IDIM+1)+I1)=0.5D0*
  521.      &                 (XCOOR((NN17-1)*(IDIM+1)+I1)+
  522.      &                 XCOOR((NN15-1)*(IDIM+1)+I1))
  523.                   XCOOR((NN18-1)*(IDIM+1)+I1)=0.5D0*
  524.      &                 (XCOOR((NN17-1)*(IDIM+1)+I1)+
  525.      &                 XCOOR((NN19-1)*(IDIM+1)+I1))
  526.                   XCOOR((NN20-1)*(IDIM+1)+I1)=0.5D0*
  527.      &                 (XCOOR((NN19-1)*(IDIM+1)+I1)+
  528.      &                 XCOOR((NN13-1)*(IDIM+1)+I1))
  529.                ENDDO
  530. C              Computation of the QUA center positions
  531.                CALL KDOM7(NN1,NN13,NN15,NN3,XCEN)
  532.                IF(IERR.NE.0) GOTO 9999
  533.                DO I1=1,IDIM,1
  534.                   XCOOR((NN21-1)*(IDIM+1)+I1)=XCEN(I1)
  535.                ENDDO
  536.                CALL KDOM7(NN3,NN15,NN17,NN5,XCEN)
  537.                IF(IERR.NE.0)GOTO 9999
  538.                DO I1=1,IDIM,1
  539.                   XCOOR((NN22-1)*(IDIM+1)+I1)=XCEN(I1)
  540.                ENDDO
  541.                CALL KDOM7(NN5,NN17,NN19,NN7,XCEN)
  542.                IF(IERR.NE.0)GOTO 9999
  543.                DO I1=1,IDIM,1
  544.                   XCOOR((NN23-1)*(IDIM+1)+I1)=XCEN(I1)
  545.                ENDDO
  546.                CALL KDOM7(NN7,NN19,NN13,NN1,XCEN)
  547.                IF(IERR.NE.0)GOTO 9999
  548.                DO I1=1,IDIM,1
  549.                   XCOOR((NN24-1)*(IDIM+1)+I1)=XCEN(I1)
  550.                ENDDO
  551.                CALL KDOM7(NN1,NN3,NN5,NN7,XCEN)
  552.                IF(IERR.NE.0)GOTO 9999
  553.                DO I1=1,IDIM,1
  554.                   XCOOR((NN25-1)*(IDIM+1)+I1)=XCEN(I1)
  555.                ENDDO
  556.                CALL KDOM7(NN13,NN19,NN17,NN15,XCEN)
  557.                IF(IERR.NE.0)GOTO 9999
  558.                DO I1=1,IDIM,1
  559.                   XCOOR((NN26-1)*(IDIM+1)+I1)=XCEN(I1)
  560.                ENDDO
  561. C
  562. C************* To compute the center, we divide the cube
  563. C              into 6 pyramids.
  564. C              (1,3,5,7,25,PC)
  565. C              (13,19,17,15,26,PC)
  566. C              (1,13,15,3,21,PC)
  567. C              (3,15,17,5,22,PC)
  568. C              (5,17,19,7,23,PC)
  569. C              (1,7,19,13,24,PC)
  570.  
  571. C              The results does not depend of the choice of
  572. C              PC.
  573. C
  574.                DO I1=1,3,1
  575.                   XCOOR((NN27-1)*(IDIM+1)+I1)=
  576.      &                 (XCOOR((NN1-1)*(IDIM+1)+I1)+
  577.      &                 XCOOR((NN3-1)*(IDIM+1)+I1)+
  578.      &                 XCOOR((NN5-1)*(IDIM+1)+I1)+
  579.      &                 XCOOR((NN7-1)*(IDIM+1)+I1)+
  580.      &                 XCOOR((NN13-1)*(IDIM+1)+I1)+
  581.      &                 XCOOR((NN15-1)*(IDIM+1)+I1)+
  582.      &                 XCOOR((NN17-1)*(IDIM+1)+I1)+
  583.      &                 XCOOR((NN19-1)*(IDIM+1)+I1))/8.0D0
  584.  
  585.                ENDDO
  586. C
  587.                CALL KDOM3(NN1,NN3,NN5,NN7,NN25,NN27,X1,VOL1)
  588.                CALL KDOM3(NN13,NN19,NN17,NN15,NN26,NN27,X2,VOL2)
  589.                CALL KDOM3(NN1,NN13,NN15,NN3,NN21,NN27,X3,VOL3)
  590.                CALL KDOM3(NN3,NN15,NN17,NN5,NN22,NN27,X4,VOL4)
  591.                CALL KDOM3(NN5,NN17,NN19,NN7,NN23,NN27,X5,VOL5)
  592.                CALL KDOM3(NN1,NN7,NN19,NN13,NN24,NN27,X6,VOL6)
  593.                VOL=VOL1+VOL2+VOL3+VOL4+VOL5+VOL6
  594.                DO I1=1,3,1
  595.                   XCOOR((NN27-1)*(IDIM+1)+I1)=(VOL1*X1(I1)+VOL2*X2(I1)
  596.      &                 +VOL3*X3(I1)+VOL4*X4(I1)+VOL5*X5(I1)+VOL6*X6(I1))
  597.      &                 /VOL
  598.                ENDDO
  599. C
  600. C********** We check whether the center of gravity of the CU27 is
  601. C           inside its domain
  602. C
  603.                CALL KDOM9(NN1,NN2,NN3,NN4,NN5,NN6,NN7,NN8,NN9,NN10,
  604.      &              NN11,NN12,NN13,NN14,NN15,NN16,NN17,NN18,NN19,NN20,
  605.      &              NN21,NN22,NN23,NN24,NN25,NN26,NN27)
  606.                IF(IERR.NE.0)GOTO 9999
  607.  
  608.             ENDDO
  609.          ELSE
  610.             WRITE(IOIMP,*) 'Euler model'
  611.             CALL ERREUR(21)
  612.          ENDIF
  613.          IF(NBSOUS .NE. 1) SEGDES IPT1
  614.       ENDDO
  615. C
  616.       SEGDES MELEME
  617. C
  618.  9999 RETURN
  619. C
  620.       END
  621.  
  622.  
  623.  
  624.  

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