Télécharger la plaquette Cast3M


fr en

Télécharger nor4d3.eso

Retour à la liste

Numérotation des lignes : 

nor4d3

  1. C NOR4D3    SOURCE    OF166741  24/10/03    21:15:27     12022          
  2.       SUBROUTINE NOR4D3(
  3.      &            MELEFA,MELEFL,MLECEN,MELEFP,MLESOM,MPONOR,
  4.      &            MPOSUR,MELTFA,MLEFA,MLEFA2,MPOTEN,MPOCHP,MLENCL,
  5.      &            MPOVCL,ICHTE,ICHCL,ICHCO,MPOVCO,IOP,
  6.      &            IPO2,SCMB,INDLI,VAL1,VAL2,IND22,IND2,IND,
  7.      &            IPO3,TAB,MPOGRA,MELVA1,MELVA2,
  8.      &            NSOMM,NBMAX,NBFAC,NBCOT,MCHEL2,MCHAM2)
  9. C
  10. C PROJET            :  CASTEM 2000
  11. C
  12. C NOM               :  NORV4
  13. C
  14. C DESCRIPTION       : Appelle par NORV1
  15. C
  16. C LANGAGE           :  FORTRAN 77 + ESOPE 2000 (avec extensions CISI)
  17. C
  18. C AUTEUR            :  C. LE POTIER, DM2S/SFME/MTMS
  19. C
  20. C************************************************************************
  21. C
  22.       IMPLICIT INTEGER(I-N)
  23.       IMPLICIT REAL*8(A-H,O-Z)
  24.  
  25. -INC PPARAM
  26. -INC CCOPTIO
  27. -INC SMLENTI
  28. -INC SMELEME
  29. -INC SMCHPOI
  30. -INC SMCOORD
  31. -INC SMLREEL
  32. -INC SMCHAML
  33.  
  34.       POINTEUR MELEFL.MELEME, MELEFP.MELEME, MELEFA.MELEME,
  35.      &         MELTFA.MELEME,MELEP2.MELEME
  36.       POINTEUR MPOSUR.MPOVAL, MPONOR.MPOVAL,
  37.      & MPOCHP.MPOVAL, MPOVCL.MPOVAL, MPGSOM.MPOVAL, MPVOSO.MPOVAL,
  38.      & MPOGRA.MPOVAL,MPOTEN.MPOVAL,MPOVCO.MPOVAL
  39.       POINTEUR MLENCL.MLENTI, MLECEN.MLENTI, MLESOM.MLENTI,
  40.      &         MLEFA.MLENTI,MLEFA2.MLENTI
  41.       INTEGER NBNN,NBREF,NBMAX
  42.  
  43. C**** Variable de SMLENTI, SMCHPOI
  44. C
  45.       INTEGER JG, N, NC,  NSOUPO, NAT, NBSOUS, NBNO,NBELEM
  46. C
  47. C**** Les includes
  48. C
  49.       INTEGER I1,ICOMP,ICOMGR,IGEOM
  50.      &     ,IOP1,ICEN,ISOMM,IFAC,IFACEL,IFACEP,INORM
  51.      &     ,ISURF,IMAIL,ICHPO,ICHCL,ICHGRA
  52.      &     ,NTOT,NSOMM,NCOMP,NFAC,NCEN
  53.      &     ,NLCF,NGCF,NGCF1,NGCF2,NGCG,NGCD,NLCG,NLCD,NGS1,NGS2
  54.      &     ,NLS1,NLS2,NLFCL
  55.      &     ,ISOUS,IELEM,INOEUD,ICELL,NCON,NFIN,INDGA,INDDR
  56.       INTEGER ICEN2
  57.       REAL*8 SCNX,SCNY,SCNZ,SURF,VOL,VAL,VALX,VALY,VALZ,XG,XD,XF,XS1,XS2
  58.      &     ,YG,YD,YF,YS1,YS2,ZG,ZD,ZF,ZS1,ZS2,
  59.      &     PSCA,XNORM,VECX,VECY,PSCAGX,PSCAGY,PSCAGZ,
  60.      &     PSCADX,PSCADY,PSCADZ,K11G,K22G,K21G,K31G,K32G,K33G,
  61.      &     K11D,K22D,K21D,K31D,K32D,K33D,VXG1,VXG2,
  62.      &     VXAU,VYAU,VZAU,VXD1,VXD2,VYG1,VYG2,VYD1,VYD2,VZG1,VZG2,VZD1,VZD2,
  63.      &     TRG1,TRG2,TRG3,
  64.      &     TRD1,TRD2,TRD3,TRG,TRD,AUX,AUY,AUZ,AUXMA,THETA,COEFDD
  65.       REAL*8 XLONG,AG1,AG2,AD1,AD2,PSCAG1,PSCAG2,PSCAD1,PSCAD2,
  66.      & COEF1,COEF2,COEF3,COEF4,SCN1X,SCN1Y,SCN1Z,VX,VY,VZ,COEF1X,COEF2X,
  67.      & COEF1Y,COEF2Y,COEF1Z,COEF2Z,CX,CY,CZ,ANCX,ANCY,ANCZ,
  68.      & DIFFX,DIFFY,DIFFZ,XLONGG,XLONGD
  69.      & VALG,COEF,GX,GY,GZ,UN,EXPR1,EXPR2
  70.  
  71. c      REAL*8 VECXG1(2),VECYG1(2)
  72. c      REAL*8 VECXG2(2),VECYG2(2)
  73. c      REAL*8 VECXD1(2),VECYD1(2)
  74. c      REAL*8 VECXD2(2),VECYD2(2)
  75.  
  76.       REAL*8 VECXG(4,4),VECYG(4,4),VECZG(4,4)
  77.       REAL*8 VECXD(4,4),VECYD(4,4),VECZD(4,4)
  78.       REAL*8 VOLUG(4),SURFAGX(4),SURFAGY(4),SURFAGZ(4),COEFG(4)
  79.       REAL*8 VOLUD(4),COEFD(4)
  80.       REAL*8 PX(3,4),PY(3,4),PZ(3,4),TRGAUX(4)
  81.       REAL*8 XPRO,TRACE
  82.       INTEGER NGS(4),NLS(4)
  83. C      INTEGER XS(4),YS(4),ZS(4)
  84.       INTEGER NLOCFG(4,4),NLOCFD(4,4)
  85.       INTEGER IGNS(4)
  86.       REAL*8 EPS
  87.       INTEGER ICRIT
  88.       CHARACTER*8 TYPE
  89. C
  90. C      CHARACTER*(4) NOMCOM(3)
  91. C      DATA NOMCOM  /'P1DX','P1DY','P1DZ'/
  92.  
  93.       INTEGER NDIM
  94.       SEGMENT MMAT1
  95.       REAL*8 PM(NDIM,NDIM),PM1(NDIM,NDIM),XSOL(NDIM)
  96.       INTEGER IC(NDIM)
  97.       ENDSEGMENT
  98.  
  99.       INTEGER K1,K2
  100.       SEGMENT INDICE
  101.       INTEGER NUME(K1,K2)
  102.       ENDSEGMENT
  103.       POINTEUR IND.INDICE,IND2.INDICE,IND22.INDICE
  104.  
  105.       SEGMENT MATRICE
  106.       REAL*8 MAT(K1,K2)
  107.       ENDSEGMENT
  108.       POINTEUR VAL1.MATRICE,VAL2.MATRICE,SCMB.MATRICE
  109.  
  110.       INTEGER K3
  111.       SEGMENT POINT2
  112.       INTEGER POINT(K3)
  113.       ENDSEGMENT
  114.       POINTEUR IPO2.POINT2
  115.  
  116.       SEGMENT MATRICE2
  117.       REAL*8 MAT2(K1,K2)
  118.       ENDSEGMENT
  119.       POINTEUR MATR1.MATRICE2,MATR2.MATRICE2
  120.  
  121.       SEGMENT POINT3
  122.       INTEGER POINT33(K3)
  123.       ENDSEGMENT
  124.       POINTEUR IPO3.POINT3
  125.  
  126.       SEGMENT INDICE3
  127.       INTEGER NU(K1,K2)
  128.       ENDSEGMENT
  129.       POINTEUR INDIC.INDICE3
  130.  
  131.       SEGMENT REP
  132.       INTEGER  ID(K3)
  133.       ENDSEGMENT
  134.       POINTEUR TAB.REP,INDLI.REP
  135.  
  136.       INTEGER K5
  137.       SEGMENT NBFAC
  138.       INTEGER NBFACEL(K5)
  139.       INTEGER IMELEM(K5)
  140.       ENDSEGMENT
  141.  
  142.       INTEGER K6
  143.       SEGMENT NBCOT
  144.       INTEGER NBCOTE(K6)
  145.       INTEGER IMECOTE(K6)
  146.       ENDSEGMENT
  147.  
  148.       INTEGER K7,K8
  149.       SEGMENT MISZERO
  150.       INTEGER  TABL(K7)
  151.       INTEGER  TABL2(K7)
  152.       INTEGER  TABL1(K8),IPOS(K8),ICOURANT(K8)
  153.       REAL*8 XMAX(K7)
  154.       ENDSEGMENT
  155.       POINTEUR ITAB.MISZERO
  156.  
  157.  
  158. c       WRITE(6,*) 'COUCOU NORV4'
  159.  
  160.  
  161.        NMAI1 = 0
  162.        NMAI2 = 0
  163.        NMAI3 = 0
  164. C       NMAI4 = 0
  165.  
  166. C      Initialisation sinon NAN...
  167.        DO II=1,4
  168.          SURFAGX(II)=0.D0
  169.          SURFAGY(II)=0.D0
  170.          SURFAGZ(II)=0.D0
  171.        ENDDO
  172.  
  173.        NBSO = MAX(1,MELTFA.LISOUS(/1))
  174. c       WRITE(6,*) 'NBSO MAILLE= ',NBSO
  175. c       WRITE(6,*) 'MELTFA= ',MELTFA
  176. C       IELTFA = MELTFA
  177.        IF (NBSO.EQ.1) THEN
  178.        K5 = MELTFA.NUM(/2)
  179.        ELSEIF (NBSO.EQ.2) THEN
  180.        IPT1 = MELTFA.LISOUS(1)
  181.        SEGACT IPT1
  182.        N1 = IPT1.NUM(/2)
  183.        NMAI1 = N1
  184.        SEGDES IPT1
  185.        IPT2 = MELTFA.LISOUS(2)
  186.        SEGACT IPT2
  187.        N2 = IPT2.NUM(/2)
  188.        NMAI2 = N2
  189.        SEGDES IPT2
  190.        K5 = N1 + N2
  191.        ELSEIF (NBSO.EQ.3) THEN
  192.        IPT1 = MELTFA.LISOUS(1)
  193.        SEGACT IPT1
  194.        N1 = IPT1.NUM(/2)
  195.        NMAI1 = N1
  196.        SEGDES IPT1
  197. c       WRITE(6,*) 'N1= ',N1
  198.        IPT2 = MELTFA.LISOUS(2)
  199.        SEGACT IPT2
  200.        N2 = IPT2.NUM(/2)
  201.        NMAI2 = N2
  202.        SEGDES IPT2
  203. c       WRITE(6,*) 'N2= ',N2
  204.        IPT3 = MELTFA.LISOUS(3)
  205.        SEGACT IPT3
  206.        N3 = IPT3.NUM(/2)
  207.        NMAI3 = N3
  208. c       WRITE(6,*) 'N3= ',N3
  209.        SEGDES IPT3
  210.        K5 = N1 + N2 + N3
  211.        ELSEIF (NBSO.EQ.4) THEN
  212.        IPT1 = MELTFA.LISOUS(1)
  213.        SEGACT IPT1
  214.        N1 = IPT1.NUM(/2)
  215.        NMAI1 = N1
  216.        SEGDES IPT1
  217.        IPT2 = MELTFA.LISOUS(2)
  218.        SEGACT IPT2
  219.        N2 = IPT2.NUM(/2)
  220.        NMAI2 = N2
  221.        SEGDES IPT2
  222.        IPT3 = MELTFA.LISOUS(3)
  223.        SEGACT IPT3
  224.        N3 = IPT3.NUM(/2)
  225.        NMAI3 = N3
  226.        SEGDES IPT3
  227.        IPT4 = MELTFA.LISOUS(4)
  228.        SEGACT IPT4
  229.        N4 = IPT4.NUM(/2)
  230. C       NMAI4 = N4
  231.        SEGDES IPT4
  232.        K5 = N1 + N2 + N3 + N4
  233.        ENDIF
  234. c       WRITE(6,*) 'K5= ',K5
  235.  
  236.  
  237.  
  238.        IF (NBSO.EQ.1) THEN
  239.        DO I = 1,K5
  240.            NTYPE = MELTFA.ITYPEL
  241. c           WRITE(6,*) 'NTYPE= ',NTYPE
  242.            IF (NTYPE .EQ. 16) THEN
  243.            NBFACEL(I) = 6
  244.            IMELEM(I) = MELTFA
  245.            ELSEIF (NTYPE .EQ. 25) THEN
  246.            NBFACEL(I) = 5
  247.            IMELEM(I) = MELTFA
  248.            ELSEIF (NTYPE .EQ. 23) THEN
  249.            NBFACEL(I) = 4
  250.            IMELEM(I) = MELTFA
  251.            ELSEIF (NTYPE .EQ. 9) THEN
  252.            NBFACEL(I) = 5
  253.            IMELEM(I) = MELTFA
  254.            ENDIF
  255. c       SEGDES MELTFA
  256.        ENDDO
  257.        ELSEIF (NBSO.EQ.2) THEN
  258.          IPT1 = MELTFA.LISOUS(1)
  259.          SEGACT IPT1
  260.          IPT2 = MELTFA.LISOUS(2)
  261.          SEGACT IPT2
  262.        DO I = 1,K5
  263.          N1 = IPT1.NUM(/2)
  264.          IF (I.LE.N1) THEN
  265.            NTYPE = IPT1.ITYPEL
  266.            IF (NTYPE .EQ. 16) THEN
  267.            NBFACEL(I) = 6
  268.            IMELEM(I) =  IPT1
  269.            ELSEIF (NTYPE .EQ. 25) THEN
  270.            NBFACEL(I) = 5
  271.            IMELEM(I) = IPT1
  272.            ELSEIF (NTYPE .EQ. 23) THEN
  273.            NBFACEL(I) = 4
  274.            IMELEM(I) = IPT1
  275.            ELSEIF (NTYPE .EQ. 9) THEN
  276.            NBFACEL(I) = 5
  277.            IMELEM(I) = IPT1
  278.            ENDIF
  279.          ELSE
  280.            NTYPE = IPT2.ITYPEL
  281.            IF (NTYPE .EQ. 16) THEN
  282.            NBFACEL(I) = 6
  283.            IMELEM(I) =  IPT2
  284.            ELSEIF (NTYPE .EQ. 25) THEN
  285.            NBFACEL(I) = 5
  286.            IMELEM(I) = IPT2
  287.            ELSEIF (NTYPE .EQ. 23) THEN
  288.            NBFACEL(I) = 4
  289.            IMELEM(I) = IPT2
  290.            ELSEIF (NTYPE .EQ. 9) THEN
  291.            NBFACEL(I) = 5
  292.            IMELEM(I) = IPT2
  293.            ENDIF
  294.         ENDIF
  295.       ENDDO
  296.       ELSEIF (NBSO.EQ.3) THEN
  297. c         WRITE(6,*) 'COUCOU'
  298.          IPT1 = MELTFA.LISOUS(1)
  299.          SEGACT IPT1
  300.            NTYPE = IPT1.ITYPEL
  301. c           WRITE(6,*) 'NTYPE= ',IPT1.ITYPEL
  302.          IPT2 = MELTFA.LISOUS(2)
  303.          SEGACT IPT2
  304.            NTYPE = IPT2.ITYPEL
  305. c           WRITE(6,*) 'NTYPE= ',IPT2.ITYPEL
  306.          IPT3 = MELTFA.LISOUS(3)
  307.          SEGACT IPT3
  308.            NTYPE = IPT3.ITYPEL
  309. c           WRITE(6,*) 'NTYPE= ',IPT3.ITYPEL
  310.          N1 = IPT1.NUM(/2)
  311.          N2 = IPT2.NUM(/2)
  312.          N3 = IPT3.NUM(/2)
  313.        DO I = 1,K5
  314.          IF (I.LE.N1) THEN
  315.            NTYPE = IPT1.ITYPEL
  316.            IF (NTYPE .EQ. 16) THEN
  317.            NBFACEL(I) = 6
  318.            IMELEM(I) =  IPT1
  319.            ELSEIF (NTYPE .EQ. 25) THEN
  320.            NBFACEL(I) = 5
  321.            IMELEM(I) = IPT1
  322.            ELSEIF (NTYPE .EQ. 23) THEN
  323.            NBFACEL(I) = 4
  324.            IMELEM(I) = IPT1
  325.            ELSEIF (NTYPE .EQ. 9) THEN
  326.            NBFACEL(I) = 5
  327.            IMELEM(I) = IPT1
  328.            ENDIF
  329.          ELSEIF (I.LE.(N1+N2)) THEN
  330.            NTYPE = IPT2.ITYPEL
  331.            IF (NTYPE .EQ. 16) THEN
  332.            NBFACEL(I) = 6
  333.            IMELEM(I) =  IPT2
  334.            ELSEIF (NTYPE .EQ. 25) THEN
  335.            NBFACEL(I) = 5
  336.            IMELEM(I) = IPT2
  337.            ELSEIF (NTYPE .EQ. 23) THEN
  338.            NBFACEL(I) = 4
  339.            IMELEM(I) = IPT2
  340.            ELSEIF (NTYPE .EQ. 9) THEN
  341.            NBFACEL(I) = 5
  342.            IMELEM(I) = IPT2
  343.            ENDIF
  344.          ELSE
  345.            NTYPE = IPT3.ITYPEL
  346.            IF (NTYPE .EQ. 16) THEN
  347.            NBFACEL(I) = 6
  348.            IMELEM(I) =  IPT3
  349.            ELSEIF (NTYPE .EQ. 25) THEN
  350.            NBFACEL(I) = 5
  351.            IMELEM(I) = IPT3
  352.            ELSEIF (NTYPE .EQ. 23) THEN
  353.            NBFACEL(I) = 4
  354.            IMELEM(I) = IPT3
  355.            ELSEIF (NTYPE .EQ. 9) THEN
  356.            NBFACEL(I) = 5
  357.            IMELEM(I) = IPT3
  358.            ENDIF
  359.          ENDIF
  360.         ENDDO
  361.       ELSEIF (NBSO.EQ.4) THEN
  362.          IPT1 = MELTFA.LISOUS(1)
  363.          SEGACT IPT1
  364.            NTYPE = IPT1.ITYPEL
  365. c           WRITE(6,*) 'NTYPE= ',IPT1.ITYPEL
  366.          IPT2 = MELTFA.LISOUS(2)
  367.          SEGACT IPT2
  368.            NTYPE = IPT2.ITYPEL
  369. c           WRITE(6,*) 'NTYPE= ',IPT2.ITYPEL
  370.          IPT3 = MELTFA.LISOUS(3)
  371.          SEGACT IPT3
  372.            NTYPE = IPT3.ITYPEL
  373. c           WRITE(6,*) 'NTYPE= ',IPT3.ITYPEL
  374.          IPT4 = MELTFA.LISOUS(4)
  375.          SEGACT IPT4
  376.            NTYPE = IPT4.ITYPEL
  377. c           WRITE(6,*) 'NTYPE= ',IPT4.ITYPEL
  378.          N1 = IPT1.NUM(/2)
  379.          N2 = IPT2.NUM(/2)
  380.          N3 = IPT3.NUM(/2)
  381.          N4 = IPT4.NUM(/2)
  382.        DO I = 1,K5
  383.          IF (I.LE.N1) THEN
  384.            NTYPE = IPT1.ITYPEL
  385.            IF (NTYPE .EQ. 16) THEN
  386.            NBFACEL(I) = 6
  387.            IMELEM(I) =  IPT1
  388.            ELSEIF (NTYPE .EQ. 25) THEN
  389.            NBFACEL(I) = 5
  390.            IMELEM(I) = IPT1
  391.            ELSEIF (NTYPE .EQ. 23) THEN
  392.            NBFACEL(I) = 4
  393.            IMELEM(I) = IPT1
  394.            ELSEIF (NTYPE .EQ. 9) THEN
  395.            NBFACEL(I) = 5
  396.            IMELEM(I) = IPT1
  397.            ENDIF
  398.          ELSEIF (I.LE.(N1+N2)) THEN
  399.            NTYPE = IPT2.ITYPEL
  400.            IF (NTYPE .EQ. 16) THEN
  401.            NBFACEL(I) = 6
  402.            IMELEM(I) =  IPT2
  403.            ELSEIF (NTYPE .EQ. 25) THEN
  404.            NBFACEL(I) = 5
  405.            IMELEM(I) = IPT2
  406.            ELSEIF (NTYPE .EQ. 23) THEN
  407.            NBFACEL(I) = 4
  408.            IMELEM(I) = IPT2
  409.            ELSEIF (NTYPE .EQ. 9) THEN
  410.            NBFACEL(I) = 5
  411.            IMELEM(I) = IPT2
  412.            ENDIF
  413.          ELSEIF (I.LE.(N1+N2+N3)) THEN
  414.            NTYPE = IPT3.ITYPEL
  415.            IF (NTYPE .EQ. 16) THEN
  416.            NBFACEL(I) = 6
  417.            IMELEM(I) =  IPT3
  418.            ELSEIF (NTYPE .EQ. 25) THEN
  419.            NBFACEL(I) = 5
  420.            IMELEM(I) = IPT3
  421.            ELSEIF (NTYPE .EQ. 23) THEN
  422.            NBFACEL(I) = 4
  423.            IMELEM(I) = IPT3
  424.            ELSEIF (NTYPE .EQ. 9) THEN
  425.            NBFACEL(I) = 5
  426.            IMELEM(I) = IPT3
  427.            ENDIF
  428.          ELSE
  429.            NTYPE = IPT4.ITYPEL
  430.            IF (NTYPE .EQ. 16) THEN
  431.            NBFACEL(I) = 6
  432.            IMELEM(I) =  IPT4
  433.            ELSEIF (NTYPE .EQ. 25) THEN
  434.            NBFACEL(I) = 5
  435.            IMELEM(I) = IPT4
  436.            ELSEIF (NTYPE .EQ. 23) THEN
  437.            NBFACEL(I) = 4
  438.            IMELEM(I) = IPT4
  439.            ELSEIF (NTYPE .EQ. 9) THEN
  440.            NBFACEL(I) = 5
  441.            IMELEM(I) = IPT4
  442.            ENDIF
  443.          ENDIF
  444.         ENDDO
  445.       ENDIF
  446.  
  447. c CAS OU LES FACES SONT DES TRIANGLES OU DES FACES
  448. C       MAUX = MELEFP
  449.        NFAI1 = 0
  450.        NBSOF = MAX(1,MELEFP.LISOUS(/1))
  451. c       WRITE(6,*) 'NBSOF FACE= ',NBSOF
  452. C       IELTFA = MELTFA
  453.        IF (NBSOF.EQ.1) THEN
  454.        K6 = MELEFP.NUM(/2)
  455.        ELSEIF (NBSOF.EQ.2) THEN
  456.        IPT5 = MELEFP.LISOUS(1)
  457.        SEGACT IPT5
  458.        N1 = IPT5.NUM(/2)
  459.        NFAI1 = N1
  460.        SEGDES IPT5
  461.        IPT6 = MELEFP.LISOUS(2)
  462.        SEGACT IPT6
  463.        N2 = IPT6.NUM(/2)
  464. C       NFAI2 = N2
  465.        SEGDES IPT6
  466.        K6 = N1 + N2
  467.        ENDIF
  468.  
  469. C ON EST ICI
  470.        IF (NBSOF.EQ.1) THEN
  471.        DO I = 1,K6
  472.            NTYPE = MELEFP.ITYPEL
  473. c           WRITE(6,*) 'NTYPE= ',NTYPE
  474.            IF (NTYPE .EQ. 5) THEN
  475.            NBCOTE(I) = 3
  476.            IMECOTE(I) = MELEFP
  477.            ELSE
  478.            NBCOTE(I) = 4
  479.            IMECOTE(I) = MELEFP
  480.            ENDIF
  481. c       SEGDES MELTFA
  482.        ENDDO
  483.        ELSEIF (NBSOF.EQ.2) THEN
  484. c       WRITE(6,*) 'POINT2'
  485.          IPT5 = MELEFP.LISOUS(1)
  486.          SEGACT IPT5
  487.          IPT6 = MELEFP.LISOUS(2)
  488.          SEGACT IPT6
  489. c          WRITE(6,*) 'IPT5= ',IPT5.ITYPEL
  490. c          WRITE(6,*) 'IPT6= ',IPT6.ITYPEL
  491.        DO I = 1,K6
  492.          N1 = IPT5.NUM(/2)
  493. C MISE A JOUR DE MLEFA.LECT
  494.          IF (I.LE.N1) THEN
  495.          N0 = IPT5.NUM(/1)
  496.          NGFAUX = IPT5.NUM(N0,I)
  497.          MLEFA2.LECT(NGFAUX) = I
  498. c           WRITE(6,*)  'NGFAUX = ',NGFAUX,
  499. c     &                 'MLEFA2=',MLEFA2.LECT(NGFAUX)
  500.            IF (IPT5.ITYPEL .EQ. 5) THEN
  501.            NBCOTE(I) = 3
  502.            IMECOTE(I) =  IPT5
  503.            ELSE
  504.            NBCOTE(I) = 4
  505.            IMECOTE(I) = IPT5
  506.            ENDIF
  507. c         SEGDES IPT5
  508.          ELSE
  509.          N0 = IPT6.NUM(/1)
  510.          NGFAUX = IPT6.NUM(N0,I-NFAI1)
  511.          MLEFA2.LECT(NGFAUX) = I
  512.            IF (IPT6.ITYPEL .EQ. 5) THEN
  513.            NBCOTE(I) = 3
  514.            IMECOTE(I) = IPT6
  515.            ELSE
  516.            NBCOTE(I) = 4
  517.            IMECOTE(I) =  IPT6
  518.            ENDIF
  519. c           SEGDES IPT6
  520.         ENDIF
  521.  
  522.         ENDDO
  523.         ENDIF
  524.  
  525.  
  526.  
  527.       NFAC=MELEFL.NUM(/2)
  528.  
  529.  
  530. c ON MAJORE SUPERIEUREMENT NBNN : ON LE REAJUSTERA PAR LA SUITE
  531.          NCON = ((2*NBMAX)) + 1
  532. c         NCON = ((3*NBMAX)/2) + 1
  533.          NBNN = NCON
  534.          NESSAI = NCON
  535.          NBNN = NESSAI
  536. c         WRITE(6,*) 'NBMAX= ',NBMAX
  537. c         WRITE(6,*) 'NBNN= ',NBNN
  538. c         WRITE(6,*) 'NFAC= ',NFAC
  539.  
  540.  
  541.  
  542. C DEFINITION DES PARAMETRES DU CHAMELEM DES COEFFICIENTS
  543.  
  544. c INITIALISATION DU CHAMELEM
  545.       N1=1
  546.       N2=1
  547. c      WRITE(6,*) 'N2= ',N2
  548.       N3=6
  549.       L1=8
  550.       SEGINI MCHELM
  551. C      ICOEFF = MCHELM
  552.       MCHELM.TITCHE='Gradient'
  553.       MCHELM.IFOCHE=IFOUR
  554. C
  555.       ISOUS=0
  556.       NBSOUS=0
  557.       NBREF=0
  558.       NBELEM = NFAC
  559.       ISOUS=ISOUS+1
  560.       SEGINI MELEME
  561. C     ITYPEL=32 -> 'POLY'
  562.       ITYPEL=32
  563.       MCHELM.IMACHE(ISOUS)=MELEME
  564.       SEGINI MCHAML
  565.       MCHELM.ICHAML(ISOUS)=MCHAML
  566.       MCHAML.NOMCHE(1)='SCAL'
  567.       MCHAML.TYPCHE(1)='REAL*8          '
  568.       N1PTEL=NESSAI
  569.       N1EL=NBELEM
  570.       N2PTEL=0
  571.       N2EL=0
  572.       SEGINI MELVA1
  573.       MCHAML.IELVAL(1)=MELVA1
  574.  
  575.  
  576.  
  577.  
  578. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  579. C CALCUL DE LA VITESSE EN CHAQUE FACE
  580. C      AA = 0.0
  581. C      BB = 0.0
  582. C INDICE QUI COMPTE LES COEFFICIENTS POUR CHAQUE FACE
  583.       NAUX2 = 0
  584.       NMOY = 0
  585.       DO NLCF= 1, NFAC, 1
  586.       NCON = 1
  587.  
  588.          NGCF=MELEFL.NUM(2,NLCF)
  589.          NGCG=MELEFL.NUM(1,NLCF)
  590.          NGCD=MELEFL.NUM(3,NLCF)
  591.          NLCG=MLECEN.LECT(NGCG)
  592.          NLCD=MLECEN.LECT(NGCD)
  593.  
  594.          SCNX=MPONOR.VPOCHA(NLCF,1)
  595.          SCNY=MPONOR.VPOCHA(NLCF,2)
  596.          SCNZ=MPONOR.VPOCHA(NLCF,3)
  597.          SCN1X = SCNX
  598.          SCN1Y = SCNY
  599.          SCN1Z = SCNZ
  600.          SURF=0.5D0*MPOSUR.VPOCHA(NLCF,1)
  601.          SCNX=SCNX*SURF
  602.          SCNY=SCNY*SURF
  603.          SCNZ=SCNZ*SURF
  604.  
  605. C        4=IDIM+1
  606.          ICELL=(4*(NGCG -1))+1
  607.          XG=MCOORD.XCOOR(ICELL)
  608.          YG=MCOORD.XCOOR(ICELL+1)
  609.          ZG=MCOORD.XCOOR(ICELL+2)
  610.          ICELL=(4*(NGCD -1))+1
  611.          XD=MCOORD.XCOOR(ICELL)
  612.          YD=MCOORD.XCOOR(ICELL+1)
  613.          ZD=MCOORD.XCOOR(ICELL+2)
  614.          ICELL=(4*(NGCF -1))+1
  615.          XF=MCOORD.XCOOR(ICELL)
  616.          YF=MCOORD.XCOOR(ICELL+1)
  617.          ZF=MCOORD.XCOOR(ICELL+2)
  618.  
  619.  
  620.          MELEME.NUM(1,NLCF)=NGCF
  621.          MELVA1.VELCHE(1,NLCF)=0.0D0
  622.  
  623.          DO J= 1,NBNN
  624.          MELVA1.VELCHE(J,NLCF) = 0.0D0
  625.          ENDDO
  626.  
  627. C MISE A ZERO DE NLOC
  628.             DO JA=1,4
  629.               DO IA=1,3
  630.               NLOCFG(IA,JA) = 0
  631.               NLOCFD(IA,JA) = 0
  632.               ENDDO
  633.             ENDDO
  634.  
  635.  
  636.          MELTFA =  IMELEM(NLCG)
  637.          NBF = NBFACEL(NLCG)
  638.          IF (NLCG.LE.NMAI1) THEN
  639.          NGAUX = NLCG
  640.          ELSEIF ((NLCG.GT.NMAI1).AND.(NLCG.LE.(NMAI1+NMAI2))) THEN
  641.          NGAUX = NLCG - NMAI1
  642.          ELSEIF ((NLCG.GT.(NMAI1+NMAI2)).AND.
  643.      &                  (NLCG.LE.(NMAI1+NMAI2+NMAI3))) THEN
  644.          NGAUX = NLCG - (NMAI1+NMAI2)
  645.          ELSEIF (NLCG.GT.(NMAI1+NMAI2+NMAI3)) THEN
  646.          NGAUX = NLCG - (NMAI1+NMAI2+NMAI3)
  647.          ENDIF
  648.  
  649.  
  650.  
  651. c ON TIENT COMPTE DU CHANGEMENT DE NUMEROTATION
  652.  
  653.  
  654.  
  655.          NLCF1 = MLEFA2.LECT(NGCF)
  656.          NBNO = NBCOTE(NLCF1)
  657.          MELEFP = IMECOTE(NLCF1)
  658.          IF (NLCF1.GT.NFAI1) THEN
  659.          NLCF1AUX = NLCF1 - NFAI1
  660.          ELSE
  661.          NLCF1AUX = NLCF1
  662.          ENDIF
  663.  
  664.  
  665.            DO IA=1,NBNO
  666.             NGS(IA) = MELEFP.NUM(IA,NLCF1AUX)
  667.             NLS(IA) = MLESOM.LECT(NGS(IA))
  668. C            ICELL=(4*(NGS(IA) -1))+1
  669. C            XS(IA)=nint(MCOORD.XCOOR(ICELL))
  670. C            YS(IA)=nint(MCOORD.XCOOR(ICELL+1))
  671. C            ZS(IA)=nint(MCOORD.XCOOR(ICELL+2))
  672.             ENDDO
  673.  
  674.  
  675.  
  676.  
  677. C ON REPERE LES VECTEURS PRINCIPAUX DE LA BASE
  678.  
  679.          DO JA = 1,NBNO
  680.          NGS(JA) = MELEFP.NUM(JA,NLCF1AUX)
  681. c         WRITE(6,*) 'NGAUX= ',NGAUX,'JA= ',JA,'NGS= ',NGS(JA)
  682. c         WRITE(6,*) 'NGCF= ',NGCF,'NLCF= ',NLCF
  683.  
  684.          ICOUR = 0
  685.          DO J = 1,NBF
  686.             N1 = MELTFA.NUM(J,NGAUX)
  687.             NL1 = MLEFA2.LECT(N1)
  688.              NBNO2 = NBCOTE(NL1)
  689.              MELEP2 = IMECOTE(NL1)
  690.              IF (NL1.GT.NFAI1) THEN
  691.              NL1AUX = NL1 - NFAI1
  692.              ELSE
  693.              NL1AUX = NL1
  694.              ENDIF
  695. c            WRITE(6,*) 'N1= ',N1,'NL1= ',NL1,'NL1AUX= ',NL1AUX
  696.  
  697.  
  698.              DO IA =1,NBNO2
  699.              NSOM1 = MELEP2.NUM(IA,NL1AUX)
  700. c             WRITE(6,*) 'NBNO2= ',NBNO2,'IA= ',IA,'NSOM1= ',NSOM1
  701.              IF (NSOM1.EQ.NGS(JA)) THEN
  702.  
  703.              ICELL=(4*(N1 -1))+1
  704.              XF=MCOORD.XCOOR(ICELL)
  705.              YF=MCOORD.XCOOR(ICELL+1)
  706.              ZF=MCOORD.XCOOR(ICELL+2)
  707.  
  708.              ICOUR = ICOUR + 1
  709.              VECXG(ICOUR,JA) = (XF - XG)
  710.              VECYG(ICOUR,JA) = (YF - YG)
  711.              VECZG(ICOUR,JA) = (ZF - ZG)
  712.              NLOCFG(ICOUR,JA) = N1
  713. C ON PERMUTE
  714.              IF (N1.EQ.NGCF) THEN
  715.              NAUX =  NLOCFG(1,JA)
  716.              VXAU = VECXG(1,JA)
  717.              VYAU = VECYG(1,JA)
  718.              VZAU = VECZG(1,JA)
  719.              VECXG(1,JA) = (XF - XG)
  720.              VECYG(1,JA) = (YF - YG)
  721.              VECZG(1,JA) = (ZF - ZG)
  722.              NLOCFG(1,JA) = N1
  723.              VECXG(ICOUR,JA) = VXAU
  724.              VECYG(ICOUR,JA) = VYAU
  725.              VECZG(ICOUR,JA) = VZAU
  726.              NLOCFG(ICOUR,JA) = NAUX
  727.              ENDIF
  728.              ENDIF
  729.             ENDDO
  730.          ENDDO
  731.          ENDDO
  732.  
  733.          MELTFA =  IMELEM(NLCD)
  734.          NBF = NBFACEL(NLCD)
  735.          IF (NLCD.LE.NMAI1) THEN
  736.          NDAUX = NLCD
  737.          ELSEIF ((NLCD.GT.NMAI1).AND.(NLCD.LE.(NMAI1+NMAI2))) THEN
  738.          NDAUX = NLCD - NMAI1
  739.          ELSEIF ((NLCD.GT.(NMAI1+NMAI2)).AND.
  740.      &                  (NLCD.LE.(NMAI1+NMAI2+NMAI3))) THEN
  741.          NDAUX = NLCD - (NMAI1+NMAI2)
  742.          ELSEIF (NLCD.GT.(NMAI1+NMAI2+NMAI3)) THEN
  743.          NDAUX = NLCD - (NMAI1+NMAI2+NMAI3)
  744.          ENDIF
  745.  
  746.  
  747. C ON REPERE LES VECTEURS PRINCIPAUX DE LA BASE
  748.          DO JA = 1,NBNO
  749.          NGS(JA) = MELEFP.NUM(JA,NLCF1AUX)
  750. c         WRITE(6,*) 'NDAUX= ',NDAUX,'JA= ',JA,'NGS= ',NGS(JA)
  751. c         WRITE(6,*) 'NGCF= ',NGCF,'NLCF= ',NLCF
  752.          ICOUR = 0
  753.          DO J = 1,NBF
  754.             N1 = MELTFA.NUM(J,NDAUX)
  755.             NL1 = MLEFA2.LECT(N1)
  756. c            WRITE(6,*) 'N1= ',N1,'NL1= ',NL1
  757.  
  758.              NBNO2 = NBCOTE(NL1)
  759.              MELEP2 = IMECOTE(NL1)
  760.              IF (NL1.GT.NFAI1) THEN
  761.              NL1AUX = NL1  - NFAI1
  762.              ELSE
  763.              NL1AUX = NL1
  764.              ENDIF
  765.  
  766.  
  767.              DO IA =1,NBNO2
  768.              NSOM1 = MELEP2.NUM(IA,NL1AUX)
  769. c             WRITE(6,*) 'NBNO2= ',NBNO2,'IA= ',IA,'NSOM1= ',NSOM1
  770.              IF (NSOM1.EQ.NGS(JA)) THEN
  771.  
  772.              ICELL=(4*(N1 -1))+1
  773.              XF=MCOORD.XCOOR(ICELL)
  774.              YF=MCOORD.XCOOR(ICELL+1)
  775.              ZF=MCOORD.XCOOR(ICELL+2)
  776.  
  777.              ICOUR = ICOUR + 1
  778.              VECXD(ICOUR,JA) = (XF - XD)
  779.              VECYD(ICOUR,JA) = (YF - YD)
  780.              VECZD(ICOUR,JA) = (ZF - ZD)
  781.              NLOCFD(ICOUR,JA) = N1
  782. C ON PERMUTE
  783.              IF (N1.EQ.NGCF) THEN
  784.              NAUX =  NLOCFD(1,JA)
  785.              VXAU = VECXD(1,JA)
  786.              VYAU = VECYD(1,JA)
  787.              VZAU = VECZD(1,JA)
  788.              VECXD(1,JA) = (XF - XD)
  789.              VECYD(1,JA) = (YF - YD)
  790.              VECZD(1,JA) = (ZF - ZD)
  791.              NLOCFD(1,JA) = N1
  792.              VECXD(ICOUR,JA) = VXAU
  793.              VECYD(ICOUR,JA) = VYAU
  794.              VECZD(ICOUR,JA) = VZAU
  795.              NLOCFD(ICOUR,JA) = NAUX
  796.              ENDIF
  797.              ENDIF
  798.             ENDDO
  799.          ENDDO
  800. c             WRITE(6,*) 'JA= ',JA
  801. c             WRITE(6,*) 'ICOUR= ',ICOUR
  802.          ENDDO
  803.  
  804.          MPOGRA.VPOCHA(NLCF,1) = 0.D0
  805.          DO JA = 1,NBNO
  806. c         WRITE(6,*) 'XPRO= ',XPRO
  807.  
  808.  
  809.           DO KA = 1,ICOUR
  810. C PRODUIT MIXTES
  811. C        PRODUIT VECTORIEL
  812.                 IF (KA.EQ.1) THEN
  813.                 PSCAGX = (VECYG(2,JA)*VECZG(3,JA)) -
  814.      &         (VECZG(2,JA)*VECYG(3,JA))
  815.                 PSCAGY = (VECZG(2,JA)*VECXG(3,JA)) -
  816.      &          (VECXG(2,JA)*VECZG(3,JA))
  817.                 PSCAGZ = (VECXG(2,JA)*VECYG(3,JA)) -
  818.      &          (VECYG(2,JA)*VECXG(3,JA))
  819.                 VOLUG(JA) = (VECXG(1,JA)* PSCAGX) +
  820.      &              (VECYG(1,JA)* PSCAGY) +
  821.      &              (VECZG(1,JA)* PSCAGZ)
  822.                 SURFAGX(KA) = 0.5D0* PSCAGX
  823.                 SURFAGY(KA) = 0.5D0* PSCAGY
  824.                 SURFAGZ(KA) = 0.5D0* PSCAGZ
  825.                 IF ( VOLUG(JA).GT.0) THEN
  826.                 SURFAGX(KA) =  - SURFAGX(KA)
  827.                 SURFAGY(KA) =  - SURFAGY(KA)
  828.                 SURFAGZ(KA) =  - SURFAGZ(KA)
  829.                 ENDIF
  830.                 VOLUG(JA) = 1.D0/6.D0*ABS(VOLUG(JA))
  831.                 ENDIF
  832.  
  833.                 IF (KA.EQ.2) THEN
  834.                 PSCAGX = (VECYG(3,JA)*VECZG(1,JA)) -
  835.      &          (VECZG(3,JA)*VECYG(1,JA))
  836.                 PSCAGY = (VECZG(3,JA)*VECXG(1,JA)) -
  837.      &          (VECXG(3,JA)*VECZG(1,JA))
  838.                 PSCAGZ = (VECXG(3,JA)*VECYG(1,JA)) -
  839.      &          (VECYG(3,JA)*VECXG(1,JA))
  840.                 SURFAGX(KA) = 0.5D0* PSCAGX
  841.                 SURFAGY(KA) = 0.5D0* PSCAGY
  842.                 SURFAGZ(KA) = 0.5D0* PSCAGZ
  843.                 PSCA = (VECXG(2,JA)* PSCAGX) + (VECYG(2,JA)* PSCAGY) +
  844.      &              (VECZG(2,JA)* PSCAGZ)
  845.                 IF ( PSCA.GT.0) THEN
  846.                 SURFAGX(KA) =  - SURFAGX(KA)
  847.                 SURFAGY(KA) =  - SURFAGY(KA)
  848.                 SURFAGZ(KA) =  - SURFAGZ(KA)
  849.                 ENDIF
  850.                 ENDIF
  851.  
  852.  
  853.                 IF (KA.EQ.3) THEN
  854.                 PSCAGX = (VECYG(1,JA)*VECZG(2,JA)) -
  855.      &          (VECZG(1,JA)*VECYG(2,JA))
  856.                 PSCAGY = (VECZG(1,JA)*VECXG(2,JA)) -
  857.      &          (VECXG(1,JA)*VECZG(2,JA))
  858.                 PSCAGZ = (VECXG(1,JA)*VECYG(2,JA)) -
  859.      &          (VECYG(1,JA)*VECXG(2,JA))
  860.  
  861.                 SURFAGX(KA) = 0.5D0* PSCAGX
  862.                 SURFAGY(KA) = 0.5D0* PSCAGY
  863.                 SURFAGZ(KA) = 0.5D0* PSCAGZ
  864.                 PSCA = (VECXG(3,JA)* PSCAGX) + (VECYG(3,JA)* PSCAGY) +
  865.      &              (VECZG(3,JA)* PSCAGZ)
  866.                 IF ( PSCA.GT.0) THEN
  867.                 SURFAGX(KA) =  - SURFAGX(KA)
  868.                 SURFAGY(KA) =  - SURFAGY(KA)
  869.                 SURFAGZ(KA) =  - SURFAGZ(KA)
  870.                 ENDIF
  871.                 ENDIF
  872.  
  873. c CALCUL DE MATRICE POUR LE NOEUD D INDICE NS1
  874.            IF (ICHTE.EQ.0) THEN
  875.            PX(KA,JA) = SURFAGX(KA)/(3.D0*VOLUG(JA))
  876.            PY(KA,JA) = SURFAGY(KA)/(3.D0*VOLUG(JA))
  877.            PZ(KA,JA) = SURFAGZ(KA)/(3.D0*VOLUG(JA))
  878.  
  879.            ELSE
  880.            IF (MPOTEN.VPOCHA(/2) .EQ.6) THEN
  881. C TENSEUR ANISOTROPE
  882.            K11G = MPOTEN.VPOCHA(NLCG,1)
  883.            K22G = MPOTEN.VPOCHA(NLCG,2)
  884.            K33G = MPOTEN.VPOCHA(NLCG,3)
  885.            K21G = MPOTEN.VPOCHA(NLCG,4)
  886.            K31G = MPOTEN.VPOCHA(NLCG,5)
  887.            K32G = MPOTEN.VPOCHA(NLCG,6)
  888.            ELSEIF (MPOTEN.VPOCHA(/2) .EQ.1) THEN
  889. C TENSEUR DIAGONAL
  890.            K11G = MPOTEN.VPOCHA(NLCG,1)
  891.            K22G = K11G
  892.            K33G = K11G
  893.            K21G = 0.0D0
  894.            K31G = 0.0D0
  895.            K32G = 0.0D0
  896.            ELSE
  897.            WRITE(6,*) 'TENSEUR NON PREVU'
  898.            STOP
  899.            ENDIF
  900.  
  901.            PSCAGX = (K11G*SURFAGX(KA)) + (K21G*SURFAGY(KA)) +
  902.      &     (K31G*SURFAGZ(KA))
  903.            PSCAGY = (K21G*SURFAGX(KA)) + (K22G*SURFAGY(KA)) +
  904.      &     (K32G*SURFAGZ(KA))
  905.            PSCAGZ = (K31G*SURFAGX(KA)) + (K32G*SURFAGY(KA)) +
  906.      &     (K33G*SURFAGZ(KA))
  907.            PX(KA,JA) = PSCAGX/(3.D0*VOLUG(JA))
  908.            PY(KA,JA) = PSCAGY/(3.D0*VOLUG(JA))
  909.            PZ(KA,JA) = PSCAGZ/(3.D0*VOLUG(JA))
  910.            ENDIF
  911.            ENDDO
  912.          ENDDO
  913.  
  914.           XPRO = 1.D0/NBNO
  915.           DO JA = 1,NBNO
  916.  
  917. C               MARQ = 0
  918.                DO I5 = 1,INDLI.ID(NLS(JA))
  919.                INDAUX = IND2.NUME(I5,NLS(JA))
  920.                IF (INDAUX.EQ.NLOCFG(2,JA)) THEN
  921.                IAFF = I5
  922.                IG2 = IAFF
  923.                GOTO 62
  924.                ENDIF
  925.                ENDDO
  926.  62            CONTINUE
  927.  
  928.                TRG2 =  SCMB.MAT(IAFF,NLS(JA))
  929. c               WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG2
  930.  
  931.                DO I5 = 1,INDLI.ID(NLS(JA))
  932.                INDAUX = IND2.NUME(I5,NLS(JA))
  933.                IF (INDAUX.EQ.NLOCFG(3,JA)) THEN
  934.                IAFF = I5
  935.                IG3 = IAFF
  936.                GOTO 629
  937.                ENDIF
  938.                ENDDO
  939.  629            CONTINUE
  940.                TRG3 =  SCMB.MAT(IAFF,NLS(JA))
  941. c               WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG3
  942.  
  943. C               MARQ = 0
  944.                DO I5 = 1,INDLI.ID(NLS(JA))
  945.                INDAUX = IND2.NUME(I5,NLS(JA))
  946.                IF (INDAUX.EQ.NGCF) THEN
  947.                IAFF = I5
  948.                IG = IAFF
  949.                GOTO 63
  950.                ENDIF
  951.                ENDDO
  952.  63              CONTINUE
  953.                TRG =  SCMB.MAT(IAFF,NLS(JA))
  954.                TRGAUX(JA) =  TRG
  955.                IGNS(JA) = IAFF
  956.  
  957. c               WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG
  958.  
  959.             VAL = MPOCHP.VPOCHA(NLCG,1)
  960. C            VALD =  MPOCHP.VPOCHA(NLCD,1)
  961. c ICI
  962.  
  963.             AUX =                     (XPRO*(
  964.      &                               (PX(1,JA) * (VAL - TRG))
  965.      &                             + (PX(2,JA) * (VAL - TRG2))
  966.      &                             + (PX(3,JA) * (VAL - TRG3))))
  967.             AUY =                             (XPRO*(
  968.      &                               (PY(1,JA) * (VAL - TRG))
  969.      &                             + (PY(2,JA) * (VAL - TRG2))
  970.      &                             + (PY(3,JA) * (VAL - TRG3))))
  971.             AUZ =                             (XPRO*(
  972.      &                               (PZ(1,JA) * (VAL - TRG))
  973.      &                             + (PZ(2,JA) * (VAL - TRG2))
  974.      &                             + (PZ(3,JA) * (VAL - TRG3))))
  975.             MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) +
  976.      &                               (AUX*SCN1X) + (AUY*SCN1Y) +
  977.      &                               (AUZ*SCN1Z)
  978.  
  979. c             IF (NLCF.EQ.791) THEN
  980. c            WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1)
  981. c            WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2)
  982. c            WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3)
  983. c            WRITE(6,*) 'PX= ',PX(1,JA),PX(2,JA),PX(3,JA)
  984. c            WRITE(6,*) 'PY= ',PY(1,JA),PY(2,JA),PY(3,JA)
  985. c            WRITE(6,*) 'PZ= ',PZ(1,JA),PZ(2,JA),PZ(3,JA)
  986. c            WRITE(6,*) 'TR ',TRG,TRG2,TRG3,'VAL',VAL
  987. c            WRITE(6,*) 'NLS=  ',NLS(JA),
  988. c     &      'IG= ',IG,'IG2= ',IG2,'IG3= ',IG3
  989. c             WRITE(6,*) 'NGCF= ',NGCF
  990. c             ENDIF
  991.  
  992.  
  993. C          ITROUVE = 0
  994.           INDIC = IPO3.POINT33(NLS(JA))
  995.           SEGACT INDIC *MOD
  996.           MATR1 = IPO2.POINT(NLS(JA))
  997.           SEGACT MATR1 *MOD
  998. c          NLS1 = NLS(JA)
  999. c          DO IAUX = 1,INDLI.ID(NLS1)
  1000. c          DO IAUX2 = 1,TAB.ID(NLS1)
  1001. c         WRITE(6,*) 'NLS1= ',NLS1,'IAUX= ',IAUX ,'IAUX2= ',
  1002. c     &                IAUX2,'VAUX',MATR1.MAT2(IAUX,IAUX2)
  1003. c     &                 ,'IND3= ',INDIC.NU(IAUX,IAUX2)
  1004. c          ENDDO
  1005. c          ENDDO
  1006.  
  1007.           DO ICOUR = 1,TAB.ID(NLS(JA))
  1008.              IA = ICOUR
  1009.                 J1 = INDIC.NU(IG,IA)
  1010.                 DO IAUX2 = 2,NCON
  1011.                 J2 = MELEME.NUM(IAUX2,NLCF)
  1012.                 IF (J1.EQ.J2) THEN
  1013.                 IAUX = IAUX2
  1014. C                ITROUVE = 1
  1015.                 GOTO 5119
  1016.                 ENDIF
  1017.               ENDDO
  1018. C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR
  1019.               NCON = NCON + 1
  1020.               IAUX = NCON
  1021.  5119         CONTINUE
  1022.  
  1023.  
  1024.                   CX = MATR1.MAT2(IG,IA)
  1025.                   CY = MATR1.MAT2(IG,IA)
  1026.                   CZ = MATR1.MAT2(IG,IA)
  1027.                   IF (J1.EQ.NGCG) THEN
  1028.                   CX =  CX - 1.D0
  1029.                   CY =  CY - 1.D0
  1030.                   CZ =  CZ - 1.D0
  1031.                   INDGA = IAUX
  1032.                   ENDIF
  1033.                   IF (J1.EQ.NGCD) THEN
  1034.                   INDDR = IAUX
  1035.                   ENDIF
  1036. c                MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1037. c     &            (XPRO*PX(1,JA)*CX)
  1038. c                MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) -
  1039. c     &            (XPRO*PY(1,JA)*CY)
  1040. c                MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) -
  1041. c     &            (XPRO*PZ(1,JA)*CZ)
  1042.                  AUX = (XPRO*PX(1,JA)*CX*SCN1X) +
  1043.      &                 (XPRO*PY(1,JA)*CY*SCN1Y) +
  1044.      &                 (XPRO*PZ(1,JA)*CZ*SCN1Z)
  1045.                 MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1046.      &          AUX
  1047.  
  1048.                  INDAUX = INDIC.NU(IG,IA)
  1049.                  MELEME.NUM(IAUX,NLCF)  = INDAUX
  1050.            ENDDO
  1051.  
  1052.  
  1053. C          ITROUVE = 0
  1054.           DO ICOUR = 1,TAB.ID(NLS(JA))
  1055.              IA = ICOUR
  1056.                 J1 = INDIC.NU(IG2,IA)
  1057.                 DO IAUX2 = 2,NCON
  1058.                 J2 = MELEME.NUM(IAUX2,NLCF)
  1059.                 IF (J1.EQ.J2) THEN
  1060.                 IAUX = IAUX2
  1061. C                ITROUVE = 1
  1062.                 GOTO 511
  1063.                 ENDIF
  1064.               ENDDO
  1065. C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR
  1066.               NCON = NCON + 1
  1067.               IAUX = NCON
  1068.  511          CONTINUE
  1069.  
  1070.                   CX = MATR1.MAT2(IG2,IA)
  1071.                   CY = MATR1.MAT2(IG2,IA)
  1072.                   CZ = MATR1.MAT2(IG2,IA)
  1073.                   IF (J1.EQ.NGCG) THEN
  1074.                   CX = CX - 1.D0
  1075.                   CY = CY - 1.D0
  1076.                   CZ = CZ - 1.D0
  1077.                   ENDIF
  1078. c                MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1079. c     &            (XPRO*PX(2,JA)*CX)
  1080. c                MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) -
  1081. c     &            (XPRO*PY(2,JA)*CY)
  1082. c                MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) -
  1083. c     &            (XPRO*PZ(2,JA)*CZ)
  1084.                  AUX = (XPRO*PX(2,JA)*CX*SCN1X) +
  1085.      &                 (XPRO*PY(2,JA)*CY*SCN1Y) +
  1086.      &                 (XPRO*PZ(2,JA)*CZ*SCN1Z)
  1087.                 MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1088.      &          AUX
  1089.                  INDAUX = INDIC.NU(IG2,IA)
  1090.                  MELEME.NUM(IAUX,NLCF)  = INDAUX
  1091.            ENDDO
  1092.  
  1093. C          ITROUVE = 0
  1094.           DO ICOUR = 1,TAB.ID(NLS(JA))
  1095.              IA = ICOUR
  1096.                 J1 = INDIC.NU(IG3,IA)
  1097.                 DO IAUX2 = 2,NCON
  1098.                 J2 = MELEME.NUM(IAUX2,NLCF)
  1099.                 IF (J1.EQ.J2) THEN
  1100.                 IAUX = IAUX2
  1101. C                ITROUVE = 1
  1102.                 GOTO 5118
  1103.                 ENDIF
  1104.               ENDDO
  1105. C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR
  1106.               NCON = NCON + 1
  1107.               IAUX = NCON
  1108.  5118          CONTINUE
  1109.  
  1110.                   CX = MATR1.MAT2(IG3,IA)
  1111.                   CY = MATR1.MAT2(IG3,IA)
  1112.                   CZ = MATR1.MAT2(IG3,IA)
  1113.                   IF (J1.EQ.NGCG) THEN
  1114.                   CX = CX - 1.D0
  1115.                   CY = CY - 1.D0
  1116.                   CZ = CZ - 1.D0
  1117.                   ENDIF
  1118. c                MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1119. c     &            (XPRO*PX(3,JA)*CX)
  1120. c                MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) -
  1121. c     &            (XPRO*PY(3,JA)*CY)
  1122. c                MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) -
  1123. c     &            (XPRO*PZ(3,JA)*CZ)
  1124.                  AUX = (XPRO*PX(3,JA)*CX*SCN1X) +
  1125.      &                 (XPRO*PY(3,JA)*CY*SCN1Y) +
  1126.      &                 (XPRO*PZ(3,JA)*CZ*SCN1Z)
  1127.                 MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1128.      &          AUX
  1129.                  INDAUX = INDIC.NU(IG3,IA)
  1130.                  MELEME.NUM(IAUX,NLCF)  = INDAUX
  1131.            ENDDO
  1132.  
  1133. c         WRITE(6,*) 'NLCF= ',NLCF,'JA= ',JA,'VOLUG(JA)  = ',VOLUG(JA),
  1134. c     &               'XPRO= ',XPRO
  1135.  
  1136. c FIN DE LA BOUCLE SUR LES NOEUDS
  1137.                  SEGDES MATR1
  1138.                  SEGDES INDIC
  1139.           ENDDO
  1140.  
  1141. c             MPOGRA.VPOCHA(NLCF,1) = 0.0D0
  1142. c             ISUP = NCON
  1143. c            DO J= ISUP+1,NBNN
  1144. c             DO J= 2,NBNN
  1145. c             MELVA1.VELCHE(J,NLCF) = 0.0D0
  1146. c             MELEME.NUM(J,NLCF) = MELEME.NUM(ISUP,NLCF)
  1147. c             ENDDO
  1148.              ISUP = NCON
  1149.              DO J= ISUP+1,NBNN
  1150.              MELVA1.VELCHE(J,NLCF) = 0.0D0
  1151.              MELEME.NUM(J,NLCF) = MELEME.NUM(ISUP,NLCF)
  1152.              ENDDO
  1153.  
  1154. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  1155. C  ON RAJOUTE LE CONVECTIF
  1156.  
  1157. C CALCUL PLUS PRECIS
  1158.  
  1159.  
  1160. c             IF (NLCD.NE.NLCG) THEN
  1161. c             MELVA1.VELCHE(INDDR,NLCF) = - MELVA1.VELCHE(INDGA,NLCF)
  1162. c             ELSE
  1163. c             DO J= 1,ISUP
  1164. c             ICENT =  MELEME.NUM(J,NLCF)
  1165. c             ICEN = MLECEN.LECT(ICENT)
  1166. c             DIFF = MELVA1.VELCHE(J,NLCF) + MELVA1.VELCHE(INDGA,NLCF)
  1167. c             DIFF = ABS(DIFF)
  1168. c             XX = ABS(MELVA1.VELCHE(INDGA,NLCF))
  1169. c             IF ((ICEN.EQ.0).AND.(DIFF.LT.(1e-5*XX))) THEN
  1170. c              MELVA1.VELCHE(J,NLCF) = - MELVA1.VELCHE(INDGA,NLCF)
  1171. c             ENDIF
  1172. c             ENDDO
  1173. c             ENDIF
  1174.  
  1175.  
  1176.              IF (ICHCO.GT.0) THEN
  1177. C BOUCLE POUR CALCUER INDGA,INDDR
  1178.              INDFR = 0
  1179.              DO J= 1,ISUP
  1180.              IF (MELEME.NUM(J,NLCF).EQ.NGCG) INDGA = J
  1181.              IF (MELEME.NUM(J,NLCF).EQ.NGCD) INDDR = J
  1182.              IF (MELEME.NUM(J,NLCF).EQ.NGCF) INDFR = J
  1183.              ENDDO
  1184.  
  1185.  
  1186.  
  1187. c             WRITE(6,*) 'NGCF= ',NGCF
  1188.              UN = MPOVCO.VPOCHA(NLCF,1)
  1189. c             WRITE(6,*) 'UN= ',UN
  1190. C OPTION CENTRE
  1191.              IF (IOP.EQ.2) THEN
  1192.              IF (NLCD.NE.NLCG) THEN
  1193.              VAL = 0.5D0*(MPOCHP.VPOCHA(NLCG,1) +
  1194.      &                        MPOCHP.VPOCHA(NLCD,1))*UN
  1195.              MPOGRA.VPOCHA(NLCF,1) =  MPOGRA.VPOCHA(NLCF,1) - VAL
  1196.              ELSE
  1197. C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION
  1198.  
  1199.               XPRO = 1.D0/NBNO
  1200.               TRACE = 0.0D0
  1201.               DO JA = 1,NBNO
  1202.               TRACE =TRACE + (XPRO*TRGAUX(JA))
  1203.               ENDDO
  1204. c             WRITE(6,*) 'NLCF= ',NLCF,'TRACE= ',TRACE
  1205.              VAL = TRACE*UN
  1206.              MPOGRA.VPOCHA(NLCF,1) =  MPOGRA.VPOCHA(NLCF,1) - VAL
  1207.              ENDIF
  1208.  
  1209. C ON COMPLETE MELVA1 POUR LE CONVECTIF
  1210.             IF (NLCD.NE.NLCG) THEN
  1211.             MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) -
  1212.      &                                  (0.5D0*UN)
  1213.             MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) -
  1214.      &                                  (0.5D0*UN)
  1215. C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDENCES DES TRACES
  1216. c POUR LES CONDITIONS MIXTES OU DE NEUMAN
  1217.  
  1218.  
  1219.             ELSE
  1220.  
  1221.             NLFCL = MLENCL.LECT(NGCF)
  1222.  
  1223. C ON RAJOUTE CECI POUR L OPTION GRADGEO
  1224.             IF (NLFCL.NE.0) THEN
  1225.             MELVA1.VELCHE(NCON+1,NLCF) =   - UN
  1226.             MELEME.NUM(NCON+1,NLCF) = NGCF
  1227.             ENDIF
  1228.  
  1229.             IF (NLFCL.EQ.0) THEN
  1230.             XPRO = 1.D0/NBNO
  1231.             DO JA = 1,NBNO
  1232.             INDIC = IPO3.POINT33(NLS(JA))
  1233.             SEGACT INDIC *MOD
  1234.             MATR1 = IPO2.POINT(NLS(JA))
  1235.             SEGACT MATR1 *MOD
  1236.  
  1237.             DO ICOUR = 1,TAB.ID(NLS(JA))
  1238.              IA = ICOUR
  1239.                 J1 = INDIC.NU(IGNS(JA),IA)
  1240.                 DO IAUX2 = 2,NCON
  1241.                 J2 = MELEME.NUM(IAUX2,NLCF)
  1242.                 IF (J1.EQ.J2) THEN
  1243.                 IAUX = IAUX2
  1244.                 GOTO 5169
  1245.                 ENDIF
  1246.               ENDDO
  1247.  5169         CONTINUE
  1248.  
  1249.  
  1250.                   CX = MATR1.MAT2(IGNS(JA),IA)
  1251.                 MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1252.      &          (UN*CX*XPRO)
  1253.             ENDDO
  1254.            SEGDES INDIC *MOD
  1255.            SEGDES MATR1 *MOD
  1256.            ENDDO
  1257.            ENDIF
  1258.            ENDIF
  1259.  
  1260. C OPTION UPWIND
  1261.             ELSEIF (IOP.EQ.1) THEN
  1262.              IF (NLCD.NE.NLCG) THEN
  1263.               IF (UN.GT.0.0D0) THEN
  1264.              VAL = (MPOCHP.VPOCHA(NLCG,1)*UN)
  1265.              ELSE
  1266.              VAL = (MPOCHP.VPOCHA(NLCD,1)*UN)
  1267.              ENDIF
  1268. c            WRITE(6,*) 'VAL= ',VAL
  1269.              MPOGRA.VPOCHA(NLCF,1) =  MPOGRA.VPOCHA(NLCF,1) - VAL
  1270.              ELSE
  1271. C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION
  1272. c ANCIENNE VERRUE
  1273.                IF (UN.GT.0.0D0) THEN
  1274.                VAL = (MPOCHP.VPOCHA(NLCG,1)*UN)
  1275.                ELSE
  1276.                XPRO = 1.D0/NBNO
  1277.                TRACE = 0.0D0
  1278.                DO JA = 1,NBNO
  1279.                TRACE =TRACE + (XPRO*TRGAUX(JA))
  1280.                ENDDO
  1281.  
  1282. c            WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF
  1283. c            WRITE(6,*) 'UN= ',UN,'TRACE=',TRACE
  1284.  
  1285.  
  1286.                VAL = TRACE*UN
  1287.                ENDIF
  1288.  
  1289.  
  1290.                MPOGRA.VPOCHA(NLCF,1) =  MPOGRA.VPOCHA(NLCF,1) - VAL
  1291.              ENDIF
  1292.  
  1293. C ON COMPLETE MELVA1 POUR LE CONVECTIF
  1294.             IF (NLCD.NE.NLCG) THEN
  1295. c            WRITE(6,*) 'UN= ',UN
  1296.                IF (UN.GT.0.0D0) THEN
  1297.                MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) -
  1298.      &                                     (UN)
  1299.                ELSE
  1300.                MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) -
  1301.      &                                     (UN)
  1302.                ENDIF
  1303. C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDENCES DES TRACES
  1304. c POUR LES CONDITIONS MIXTES OU DE NEUMAN
  1305.  
  1306.             ELSE
  1307. c ANCIENNE VERRUE
  1308.             IF (UN.GT.0.0D0) THEN
  1309.                MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) -
  1310.      &                                     (UN)
  1311.             ELSE
  1312.  
  1313.             NLFCL = MLENCL.LECT(NGCF)
  1314.  
  1315. C ON RAJOUTE CECI POUR L OPTION GRADGEO
  1316. C IL Y A PROBABLEMENT REDONDANCE
  1317.             IF (NLFCL.NE.0) THEN
  1318.             MELVA1.VELCHE(NCON+1,NLCF) =   - UN
  1319.             MELEME.NUM(NCON+1,NLCF) = NGCF
  1320.             ENDIF
  1321. * ENDIF
  1322.  
  1323.             IF (NLFCL.EQ.0) THEN
  1324.             XPRO = 1.D0/NBNO
  1325.             DO JA = 1,NBNO
  1326.             INDIC = IPO3.POINT33(NLS(JA))
  1327.             SEGACT INDIC *MOD
  1328.             MATR1 = IPO2.POINT(NLS(JA))
  1329.             SEGACT MATR1 *MOD
  1330.  
  1331.             DO ICOUR = 1,TAB.ID(NLS(JA))
  1332.              IA = ICOUR
  1333.                 J1 = INDIC.NU(IGNS(JA),IA)
  1334.                 DO IAUX2 = 2,NCON
  1335.                 J2 = MELEME.NUM(IAUX2,NLCF)
  1336.                 IF (J1.EQ.J2) THEN
  1337.                 IAUX = IAUX2
  1338.                 GOTO 5129
  1339.                 ENDIF
  1340.               ENDDO
  1341.  5129         CONTINUE
  1342.  
  1343.  
  1344.                   CX = MATR1.MAT2(IGNS(JA),IA)
  1345.                 MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1346.      &          (UN*CX*XPRO)
  1347.             ENDDO
  1348.  
  1349.            SEGDES INDIC *MOD
  1350.            SEGDES MATR1 *MOD
  1351.            ENDDO
  1352.            ENDIF
  1353.            ENDIF
  1354.  
  1355.  
  1356.            ENDIF
  1357.  
  1358.  
  1359. C OPTION UPWICENT
  1360.             ELSEIF (IOP.EQ.3) THEN
  1361.              IF (NLCD.NE.NLCG) THEN
  1362. c CALCUL DE THETA
  1363. c DANS MELVA1 :  INDDR EST LA NUM de NLCD
  1364. C COEFD = MELVA1(INDDR,NLCF)
  1365.              COEFDD =   MELVA1.VELCHE(INDDR,NLCF)
  1366.  
  1367.  
  1368.              THETA = 0.5D0
  1369.              EXPR1 =  - (THETA*UN) - COEFDD
  1370.              EXPR2 =  - ((1.D0 - THETA)*UN) + COEFDD
  1371. c             WRITE(6,*) 'EXPR1',EXPR1
  1372. c             WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF,
  1373. c     &       'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN
  1374. c             WRITE(6,*) 'EXPR2',EXPR2
  1375. c             WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF,
  1376. c     &       'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN
  1377.              IF (EXPR1.GT.0.0D0) THEN
  1378.                 IF (ABS(UN) .GT. 1e-20) THEN
  1379.                  THETA = - (0.5D0*COEFDD / UN)
  1380.                 ENDIF
  1381.              ENDIF
  1382.              IF (EXPR2.LT.0.0D0) THEN
  1383.                 IF (ABS(UN) .GT. 1e-20) THEN
  1384.                  THETA = 1.D0 - (0.5D0*COEFDD / UN)
  1385.                 ENDIF
  1386.              ENDIF
  1387.              THETA = MIN(1.D0,THETA)
  1388.              THETA = MAX(0.D0,THETA)
  1389.  
  1390.              EXPR1 = (THETA*UN) + COEFDD
  1391.              EXPR2 = ((1.D0 - THETA)*UN) - COEFDD
  1392. c             IF (EXPR1.LT.0.0D0) THEN
  1393. c             WRITE(6,*) 'EXPR1',EXPR1
  1394. c             WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF,
  1395. c     &       'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN
  1396. c             ENDIF
  1397. c             IF (EXPR2.GT.0.0D0) THEN
  1398. c             WRITE(6,*) 'EXPR2',EXPR2
  1399. c             WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF,
  1400. c     &       'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN
  1401. c             ENDIF
  1402. c             WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF,
  1403. c     &       'EXPR1= ',EXPR1,'EXPR2= ',EXPR2
  1404.  
  1405.  
  1406.              VAL = (MPOCHP.VPOCHA(NLCG,1)*UN*THETA) +
  1407.      &             (MPOCHP.VPOCHA(NLCD,1)*(1.D0-THETA)*UN)
  1408.              MPOGRA.VPOCHA(NLCF,1) =  MPOGRA.VPOCHA(NLCF,1) - VAL
  1409.              ELSE
  1410. C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION
  1411.              IF (INDFR.NE.0) THEN
  1412.              COEFDD =   MELVA1.VELCHE(INDFR,NLCF)
  1413.              ELSE
  1414.              COEFDD = 0.0D0
  1415.              ENDIF
  1416. C C'EST LE THETA LE MIEUX ADPATE POUR LES CONDITIONS AUX LIMITES
  1417.              THETA = 0.0D0
  1418.              EXPR1 =  - (THETA*UN) - COEFDD
  1419.              EXPR2 =  - ((1.D0 - THETA)*UN) + COEFDD
  1420.  
  1421.              IF (EXPR1.GT.0.0D0) THEN
  1422.                 IF (ABS(UN) .GT. 1e-20) THEN
  1423.                  THETA = - (0.5D0*COEFDD / UN)
  1424.                 ENDIF
  1425.              ENDIF
  1426.              IF (EXPR2.LT.0.0D0) THEN
  1427.                 IF (ABS(UN) .GT. 1e-20) THEN
  1428.                  THETA = 1.D0 - (0.5D0*COEFDD / UN)
  1429.                 ENDIF
  1430.              ENDIF
  1431.              THETA = MIN(1.D0,THETA)
  1432.              THETA = MAX(0.D0,THETA)
  1433. C ANCIENNE VERRUE
  1434. c             THETA = 0.0D0
  1435.  
  1436.               XPRO = 1.D0/NBNO
  1437.               TRACE = 0.0D0
  1438.               DO JA = 1,NBNO
  1439.               TRACE =TRACE + (XPRO*TRGAUX(JA))
  1440.               ENDDO
  1441.               VAL = TRACE*UN
  1442.              VAL = (MPOCHP.VPOCHA(NLCG,1)*UN*THETA) +
  1443.      &             (TRACE*(1.D0-THETA)*UN)
  1444.              MPOGRA.VPOCHA(NLCF,1) =  MPOGRA.VPOCHA(NLCF,1) - VAL
  1445.  
  1446. c            WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF,
  1447. c     &       'GR1= ',MPOGRA.VPOCHA(NLCF,1),'VAL = ',VAL,
  1448. c     &         'TRACE = ',TRACE,'UN= ',UN
  1449. c              MPOGRA.VPOCHA(NLCF,1) =  MPOGRA.VPOCHA(NLCF,1) - VAL
  1450. c            WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF,
  1451. c     &       'GR1= ',MPOGRA.VPOCHA(NLCF,1)
  1452.              ENDIF
  1453.  
  1454. C ON COMPLETE MELVA1 POUR LE CONVECTIF
  1455.             IF (NLCD.NE.NLCG) THEN
  1456.               MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) -
  1457.      &                                    (THETA*UN)
  1458.               MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) -
  1459.      &                                    ((1.D0-THETA)*UN)
  1460. C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDANCES DES TRACES
  1461. c POUR LES CONDITIONS MIXTES OU DE NEUMAN
  1462.  
  1463.             ELSE
  1464.  
  1465.              MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) -
  1466.      &                                    (THETA*UN)
  1467.             NLFCL = MLENCL.LECT(NGCF)
  1468. C ON EST ICI : CORRIGER
  1469. C ON RAJOUTE CECI POUR L OPTION GRADGEO
  1470. c            IF (NLFCL.NE.0) THEN
  1471. c            MELVA1.VELCHE(NCON+1,NLCF) =   - UN
  1472. c            MELEME.NUM(NCON+1,NLCF) = NGCF
  1473. c            ENDIF
  1474.  
  1475. c            IF (NLFCL.EQ.0) THEN
  1476.             XPRO = 1.D0/NBNO
  1477.             DO JA = 1,NBNO
  1478.             INDIC = IPO3.POINT33(NLS(JA))
  1479.             SEGACT INDIC *MOD
  1480.             MATR1 = IPO2.POINT(NLS(JA))
  1481.             SEGACT MATR1 *MOD
  1482.  
  1483.             DO ICOUR = 1,TAB.ID(NLS(JA))
  1484.              IA = ICOUR
  1485.                 J1 = INDIC.NU(IGNS(JA),IA)
  1486.                 DO IAUX2 = 2,NCON
  1487.                 J2 = MELEME.NUM(IAUX2,NLCF)
  1488.                 IF (J1.EQ.J2) THEN
  1489.                 IAUX = IAUX2
  1490.                 GOTO 5159
  1491.                 ENDIF
  1492.               ENDDO
  1493.  5159         CONTINUE
  1494.  
  1495.  
  1496.                   CX = MATR1.MAT2(IGNS(JA),IA)
  1497.                 MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) -
  1498.      &          ((1.D0-THETA)*UN*CX*XPRO)
  1499.             ENDDO
  1500.  
  1501.            SEGDES INDIC *MOD
  1502.            SEGDES MATR1 *MOD
  1503.            ENDDO
  1504. c           ENDIF
  1505.            ENDIF
  1506.  
  1507.            ENDIF
  1508.            ENDIF
  1509.  
  1510. c            DO J= 1,NBNN
  1511. c             PSCA = (MPOGRA.VPOCHA(NLCF,1)*SCN1X) +
  1512. c     &              (MPOGRA.VPOCHA(NLCF,2)*SCN1Y) +
  1513. c     &              (MPOGRA.VPOCHA(NLCF,3)*SCN1Z)
  1514. c             MPOGRA.VPOCHA(NLCF,1) = PSCA*SCN1X
  1515. c             MPOGRA.VPOCHA(NLCF,2) = PSCA*SCN1Y
  1516. c             MPOGRA.VPOCHA(NLCF,3) = PSCA*SCN1Z
  1517. c             ENDDO
  1518. c            IF (NLCF.EQ.791) THEN
  1519. c            WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1)
  1520. c            ENDIF
  1521. c            WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2)
  1522. c            WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3)
  1523.  
  1524.  
  1525. c             DO J= 1,NBNN
  1526. c             PSCA = (MELVA1.VELCHE(J,NLCF)*SCN1X) +
  1527. c     &              (MELVA2.VELCHE(J,NLCF)*SCN1Y) +
  1528. c     &              (MELVA3.VELCHE(J,NLCF)*SCN1Z)
  1529. c             MELVA1.VELCHE(J,NLCF) = PSCA*SCN1X
  1530. c             MELVA2.VELCHE(J,NLCF) = PSCA*SCN1Y
  1531. c             MELVA3.VELCHE(J,NLCF) = PSCA*SCN1Z
  1532. c             ENDDO
  1533. c             DO J= 1,NBNN
  1534. c            WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA1= ',
  1535. c     &     MELVA1.VELCHE(J,NLCF)
  1536. c            WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA2= ',
  1537. c     &     MELVA2.VELCHE(J,NLCF)
  1538. c            WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA3= ',
  1539. c     &     MELVA3.VELCHE(J,NLCF)
  1540. c              WRITE(6,*) 'MELEME=',MELEME.NUM(J,NLCF)
  1541. c             ENDDO
  1542. C            IF (1.EQ.0) THEN
  1543. C            WRITE(6,*) 'NGCF= ',NGCF
  1544. C            WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1)
  1545. C            WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2)
  1546. C            WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3)
  1547. C            VALD = MPOCHP.VPOCHA(NLCD,1)
  1548. C            WRITE(6,*) 'NLCF= ',NLCF,'VALD= ',VALD
  1549. C            WRITE(6,*) 'NLCF= ',NLCF,'VAL= ',VAL
  1550. C            WRITE(6,*) 'NLCF= ',NLCF,'TRG= ',TRG
  1551. C            WRITE(6,*) 'NLCF= ',NLCF,'TRG2= ',TRG2
  1552. C            WRITE(6,*) 'NLCF= ',NLCF,'TRG3= ',TRG3
  1553. C            WRITE(6,*) 'NLCF= ',NLCF,'TRG= ',TRGAUX
  1554. C            WRite(6,*) 'AG1=',AG1
  1555. C            WRite(6,*) 'AG2=',AG2
  1556. C            WRite(6,*) 'AD1=',AD1
  1557. C            WRite(6,*) 'AD2=',AD2
  1558. C            ENDIF
  1559.             NAUX2 = MAX(NAUX2,NCON)
  1560.             NMOY = NMOY + NCON
  1561.  
  1562. c             WRITE(6,*) 'INDGA= ',INDGA
  1563. c             WRITE(6,*) 'INDDR= ',INDDR
  1564. c             DO J= 1,NBNN
  1565. c             WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA1= ',
  1566. c     &       MELVA1.VELCHE(J,NLCF)
  1567. c             ENDDO
  1568.  
  1569.  
  1570.          ENDDO
  1571.          NMOY = NMOY/(NFAC*1.D0)
  1572.  
  1573.          IF (1.EQ.0) THEN
  1574.          DO NLCF = 1, NFAC, 1
  1575.          MPOGRA.VPOCHA(NLCF,1) = 0.D0
  1576.          SCNX=MPONOR.VPOCHA(NLCF,1)
  1577.          SCNY=MPONOR.VPOCHA(NLCF,2)
  1578.          SCNZ=MPONOR.VPOCHA(NLCF,3)
  1579.          SCN1X = SCNX
  1580.          SCN1Y = SCNY
  1581.          SCN1Z = SCNZ
  1582.             NGCF=MELEFL.NUM(2,NLCF)
  1583.              DO IVOI=2,MELEME.NUM(/1)
  1584.              ICENT =  MELEME.NUM(IVOI,NLCF)
  1585.              ICEN = MLECEN.LECT(ICENT)
  1586.              VAL = 0.0D0
  1587.              IF (ICEN.NE.0) THEN
  1588. c             WRITE(6,*) 'INTERIEUR'
  1589.              VAL = MPOCHP.VPOCHA(ICEN,1)
  1590.              ELSE
  1591.              ICENL = MLENCL.LECT(ICENT)
  1592.              IF (ICENL.GT.0) THEN
  1593. c              WRITE(6,*) 'DIRICHLET'
  1594.              VAL = MPOVCL.VPOCHA(ICENL,1)
  1595.              ENDIF
  1596.              ENDIF
  1597.  
  1598.              COEF1 =  MELVA1.VELCHE(IVOI,NLCF)
  1599.              MPOGRA.VPOCHA(NLCF,1)= MPOGRA.VPOCHA(NLCF,1) +
  1600.      &       (COEF1 * VAL)
  1601.  
  1602. c            WRITE(6,*) 'NLCF= ',NLCF,'VAL= ',VAL
  1603. c            WRITE(6,*) 'IVOI= ',IVOI,'MELEME= ', MELEME.NUM(IVOI,NLCF),
  1604. c     &                 'COEF1 = ',COEF1,'COEF2= ',COEF2,'COEF3= ',COEF3
  1605.              ENDDO
  1606. c             DO J= 1,NBNN
  1607. c             PSCA = (MPOGRA.VPOCHA(NLCF,1)*SCN1X) +
  1608. c     &              (MPOGRA.VPOCHA(NLCF,2)*SCN1Y) +
  1609. c     &              (MPOGRA.VPOCHA(NLCF,3)*SCN1Z)
  1610. c             MPOGRA.VPOCHA(NLCF,1) = PSCA*SCN1X
  1611. c             MPOGRA.VPOCHA(NLCF,2) = PSCA*SCN1Y
  1612. c             MPOGRA.VPOCHA(NLCF,3) = PSCA*SCN1Z
  1613. c             ENDDO
  1614. c            NGCF=MELEFL.NUM(2,NLCF)
  1615. c            WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF
  1616. c            WRITE(6,*) 'MPOGRA1= ', MPOGRA.VPOCHA(NLCF,1)
  1617. c            WRITE(6,*) 'MPOGRA2= ', MPOGRA.VPOCHA(NLCF,2)
  1618. c            WRITE(6,*) 'MPOGRA3= ', MPOGRA.VPOCHA(NLCF,3)
  1619.          ENDDO
  1620.              ENDIF
  1621.  
  1622.       IF (NBSO.EQ.2) THEN
  1623.       SEGDES  IPT1
  1624.       SEGDES IPT2
  1625.       ELSEIF (NBSO.EQ.3) THEN
  1626.       SEGDES IPT1
  1627.       SEGDES IPT2
  1628.       SEGDES IPT3
  1629.       ELSEIF (NBSO.EQ.4) THEN
  1630.       SEGDES IPT1
  1631.       SEGDES IPT2
  1632.       SEGDES IPT3
  1633.       SEGDES IPT4
  1634.       ENDIF
  1635.       IF (NBSOF.EQ.2) THEN
  1636.       SEGDES IPT5
  1637.       SEGDES IPT6
  1638.       ENDIF
  1639.  
  1640. C ON REJUSTE LE CHAMELEM
  1641.  
  1642.  
  1643. c         WRITE(6,*) 'NAUX2= ',NAUX2
  1644. c         WRITE(6,*) 'NMOY2= ',NMOY
  1645. c         IF (NAUX2.GT.NESSAI) THEN
  1646. c         WRITE(6,*) 'NAUX2 = ',NAUX2,'NESSAI = ',NESSAI
  1647. c         WRITE(6,*) 'NESSAI TROP PETIT'
  1648. c         STOP
  1649. c         ENDIF
  1650.          N1PTEL=NAUX2
  1651.          N1EL=NBELEM
  1652.          N2PTEL=0
  1653.          N2EL=0
  1654.          SEGADJ MELVA1
  1655.  
  1656.          NBNN = NAUX2
  1657.          SEGADJ MELEME
  1658.  
  1659.          K7 = NFAC
  1660.          K8 = NAUX2
  1661.          SEGINI ITAB
  1662.  
  1663. c ON SUPPRIME LES ZEROS INTERIEURS
  1664.          IND = 2
  1665.          DO NLCF = 1,NFAC
  1666.             IND = 2
  1667.             NFIN = NAUX2
  1668. C ON CALCULE D4ABORD LE MAXIMUM DE LA LIGNE
  1669.             AUXMA = 0
  1670.             DO J=2,NFIN
  1671.             ICEN = MLECEN.LECT(MELEME.NUM(J,NLCF))
  1672.  
  1673.             IF (ICEN.NE.0) THEN
  1674.             AUXMA = MAX(ABS(MELVA1.VELCHE(J,NLCF)),AUXMA)
  1675.             ENDIF
  1676.             ENDDO
  1677.             ITAB.XMAX(NLCF) = AUXMA
  1678.  
  1679.             DO  J=IND, NFIN
  1680.             ICEN = MLECEN.LECT(MELEME.NUM(J,NLCF))
  1681.             IF (ABS(MELVA1.VELCHE(J,NLCF)).Le.(1e-14*AUXMA)
  1682.      &          .AND.(ICEN.NE.0) ) THEN
  1683.  
  1684.  
  1685.                  DO K=1,NFIN-J
  1686.                  AUX = MELVA1.VELCHE(J+K,NLCF)
  1687.                  ICEN = MLECEN.LECT(MELEME.NUM(J+K,NLCF))
  1688. C ON DECALE DE K CRAN
  1689.                  IF (ABS(AUX).gt.(1e-14*AUXMA).OR.(ICEN.EQ.0)) THEN
  1690.                   DO I=J,NFIN - K
  1691.                    MELVA1.VELCHE(I,NLCF) = MELVA1.VELCHE(I+K,NLCF)
  1692.                    MELEME.NUM(I,NLCF) = MELEME.NUM(I+K,NLCF)
  1693.                   ENDDO
  1694. C MISE A ZERO DES TERMES DU BOUT
  1695.                   DO I=NFIN-K+1,NFIN
  1696.                    MELVA1.VELCHE(I,NLCF) = 0.0
  1697.                   ENDDO
  1698.                  GOTO 2000
  1699.                  ENDIF
  1700.                  ENDDO
  1701. 2000        CONTINUE
  1702.             ENDIF
  1703.             ENDDO
  1704.          ENDDO
  1705.  
  1706. c         IF (1.EQ.0) THEN
  1707. c ON RECONSTUIT UN CHAMELEM EN SUPPRIMANT LES 0
  1708.  
  1709.  
  1710.  
  1711.  
  1712.  
  1713. c TABLEAU CALCULANT LE NOMBRE DE VOISIN NON NUL POUR CHAQUE FACE
  1714.          NMOY = 0
  1715. C         INF = NAUX2
  1716.          DO NLCF = 1,NFAC
  1717.            IREC = 3
  1718.            DO J=NAUX2,1,-1
  1719.            IF (ABS(MELVA1.VELCHE(J,NLCF)).gt.
  1720.      &        (1e-14*ITAB.XMAX(NLCF))) THEN
  1721.            IREC = J
  1722.            GOTO 1111
  1723.            ENDIF
  1724.            ENDDO
  1725. 1111     CONTINUE
  1726.            ITAB.TABL(NLCF) = IREC
  1727.            NMOY = NMOY + IREC
  1728.         ENDDO
  1729.          NMOY = NMOY/(NFAC*1.D0)
  1730. c         WRITE(6,*) 'NEWMOY2= ',NMOY
  1731.  
  1732.  
  1733. C   TAB1(U) TABLEAU QUI CONTIENT LE NOMBRE DE FACE AYANT U VOISIN
  1734. C          NMAX = 0
  1735.           DO ICOUR = 1,NAUX2
  1736.           ITAB.TABL1(ICOUR) = 0
  1737.           ENDDO
  1738.  
  1739.           DO NLCF = 1,NFAC
  1740.           ICOUR =  ITAB.TABL(NLCF)
  1741.           ITAB.TABL1(ICOUR)=ITAB.TABL1(ICOUR) + 1
  1742.           ENDDO
  1743.  
  1744. C ON COMPTE LE NOMVRE DE SOUS DOMAINE
  1745.           NTSOUS = 0
  1746.           DO ICOUR = 1,NAUX2
  1747.           IF (ITAB.TABL1(ICOUR) .NE.0) NTSOUS = NTSOUS +1
  1748.           ENDDO
  1749. C IPOS INDICE DE LA PREMIERE FACE AYANT I VOISIN
  1750. C ICOUR INDICE COURANT INITIALISE A IPOS
  1751.  
  1752.  
  1753.           ITAB.IPOS(1) = 1
  1754.           DO I = 2,NAUX2
  1755.           ITAB.IPOS(I) = ITAB.IPOS(I-1) + ITAB.TABL1(I-1)
  1756.           ITAB.ICOURANT(I) = ITAB.IPOS(I)
  1757.           ENDDO
  1758.  
  1759.  
  1760. c TABL2 TABLEAU QUI RANGE DANS L4ODRES DES SOUS DOMAINES LES FACES NLCF
  1761.  
  1762.            DO I =1,NFAC
  1763.            IHELP = ITAB.TABL(I)
  1764.            IAUX = ITAB.ICOURANT(IHELP)
  1765.            ITAB.TABL2(IAUX) = I
  1766.            IAUX2 = ITAB.TABL(I)
  1767.            ITAB.ICOURANT(IAUX2) = ITAB.ICOURANT(IAUX2) + 1
  1768.            ENDDO
  1769.  
  1770.  
  1771. C**** Initialisation du MCHELM
  1772. C
  1773.       N1=NTSOUS
  1774.       N2=1
  1775.       N3=6
  1776.       L1=8
  1777.       SEGINI MCHEL2
  1778.       MCHEL2.TITCHE='Gradient'
  1779.       MCHEL2.IFOCHE=IFOUR
  1780. C
  1781.       ISOUS=0
  1782.       NBSOUS=0
  1783.       NBREF=0
  1784.       DO I1 = 1, NAUX2, 1
  1785.          NBELEM=ITAB.TABL1(I1)
  1786.          IF(NBELEM .GT. 0)THEN
  1787.             ISOUS=ISOUS+1
  1788.             NBNN=I1
  1789.             SEGINI IPT8
  1790. C           ITYPEL=32 -> 'POLY'
  1791.             ITYPEL=32
  1792.             MCHEL2.IMACHE(ISOUS)=IPT8
  1793.             MCHEL2.CONCHE(ISOUS)='    '
  1794.             MCHEL2.INFCHE(ISOUS,6)=1
  1795.             SEGINI MCHAM2
  1796.             MCHEL2.ICHAML(ISOUS)=MCHAM2
  1797.             MCHAM2.NOMCHE(1)='SCAL'
  1798.             MCHAM2.TYPCHE(1)='REAL*8          '
  1799.             N1PTEL=NBNN
  1800.             N1EL=NBELEM
  1801.             N2PTEL=0
  1802.             N2EL=0
  1803.             SEGINI MELVA2
  1804.             MCHAM2.IELVAL(1)=MELVA2
  1805.  
  1806.             DO I2=1,NBELEM,1
  1807.                DO I3=1,NBNN,1
  1808.                   IAUX = ITAB.IPOS(I1)
  1809.                   IP= ITAB.TABL2(I2 + IAUX - 1)
  1810. c             WRITE(6,*) 'IP= ',IP,'I1= ',I1,'I2=',I2,'I3=',I3
  1811.                   IPT8.NUM(I3,I2)= MELEME.NUM(I3,IP)
  1812.                   MELVA2.VELCHE(I3,I2)=MELVA1.VELCHE(I3,IP)
  1813.                ENDDO
  1814.             ENDDO
  1815. C
  1816.             SEGDES MCHAM2
  1817.             SEGDES IPT8
  1818.             SEGDES MELVA2
  1819.          ENDIF
  1820.       ENDDO
  1821.          SEGDES MCHEL2
  1822.          SEGDES ITAB
  1823. c        ENDIF
  1824.  
  1825. c         SEGDES MCHAML
  1826. c         SEGDES MELEME
  1827. c         SEGDES MELVA1
  1828. c         SEGDES MCHELM
  1829.  
  1830. C 3009
  1831.          SEGSUP MCHAML
  1832.          SEGSUP MELEME
  1833.          SEGSUP MELVA1
  1834.          SEGSUP MCHELM
  1835.  
  1836. C SUPRESSION DES SEGMENTS
  1837.          K3 = NSOMM
  1838.          DO I = 1,K3
  1839.          MATR1 = IPO2.POINT(I)
  1840.          SEGACT MATR1
  1841.          SEGSUP  MATR1
  1842.          ENDDO
  1843.  
  1844.          SEGSUP IPO3
  1845.          SEGSUP IPO2
  1846.          SEGSUP INDLI
  1847.          SEGSUP TAB
  1848.          SEGSUP IND2
  1849. c         SEGSUP IND
  1850. c         SEGSUP IND22
  1851. c         SEGSUP VAL1
  1852. c         SEGSUP VAL2
  1853.          SEGSUP SCMB
  1854. C 3009
  1855.          SEGSUP NBCOT
  1856.  
  1857.       RETURN
  1858.       END
  1859.  
  1860.  
  1861.  

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