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

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