C NOR4D3 SOURCE CB215821 20/11/25 13:34:57 10792 SUBROUTINE NOR4D3( & MELEFA,MELEFL,MLECEN,MELEFP,MLESOM,MPONOR, & MPOSUR,MELTFA,MLEFA,MLEFA2,MPOTEN,MPOCHP,MLENCL, & MPOVCL,ICHTE,ICHCL,ICHCO,MPOVCO,IOP, & IPO2,SCMB,INDLI,VAL1,VAL2,IND22,IND2,IND, & IPO3,TAB,MPOGRA,MELVA1,MELVA2, & NSOMM,NBMAX,NBFAC,NBCOT,MCHEL2,MCHAM2) C C PROJET : CASTEM 2000 C C NOM : NORV4 C C DESCRIPTION : Appelle par NORV1 C C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec extensions CISI) C C AUTEUR : C. LE POTIER, DM2S/SFME/MTMS C C************************************************************************ C IMPLICIT INTEGER(I-N) IMPLICIT REAL*8(A-H,O-Z) -INC PPARAM -INC CCOPTIO -INC SMLENTI -INC SMELEME -INC SMCHPOI -INC SMCOORD -INC SMLREEL -INC SMCHAML POINTEUR MELEFL.MELEME, MELEFP.MELEME, MELEFA.MELEME, & MELTFA.MELEME,MELEP2.MELEME POINTEUR MPOSUR.MPOVAL, MPONOR.MPOVAL, & MPOCHP.MPOVAL, MPOVCL.MPOVAL, MPGSOM.MPOVAL, MPVOSO.MPOVAL, & MPOGRA.MPOVAL,MPOTEN.MPOVAL,MPOVCO.MPOVAL POINTEUR MLENCL.MLENTI, MLECEN.MLENTI, MLESOM.MLENTI, & MLEFA.MLENTI,MLEFA2.MLENTI INTEGER NBNN,NBREF,NBMAX C**** Variable de SMLENTI, SMCHPOI C INTEGER JG, N, NC, NSOUPO, NAT, NBSOUS, NBNO,NBELEM C C**** Les includes C INTEGER I1,ICOMP,ICOMGR,IGEOM & ,IOP1,ICEN,ISOMM,IFAC,IFACEL,IFACEP,INORM & ,ISURF,IMAIL,ICHPO,ICHCL,ICHGRA & ,NTOT,NSOMM,NCOMP,NFAC,NCEN & ,NLCF,NGCF,NGCF1,NGCF2,NGCG,NGCD,NLCG,NLCD,NGS1,NGS2 & ,NLS1,NLS2,NLFCL & ,ISOUS,IELEM,INOEUD,ICELL,NCON,NFIN,INDGA,INDDR INTEGER ICEN2 REAL*8 SCNX,SCNY,SCNZ,SURF,VOL,VAL,VALX,VALY,VALZ,XG,XD,XF,XS1,XS2 & ,YG,YD,YF,YS1,YS2,ZG,ZD,ZF,ZS1,ZS2, & PSCA,XNORM,VECX,VECY,PSCAGX,PSCAGY,PSCAGZ, & PSCADX,PSCADY,PSCADZ,K11G,K22G,K21G,K31G,K32G,K33G, & K11D,K22D,K21D,K31D,K32D,K33D,VXG1,VXG2, & VXAU,VYAU,VZAU,VXD1,VXD2,VYG1,VYG2,VYD1,VYD2,VZG1,VZG2,VZD1,VZD2, & TRG1,TRG2,TRG3, & TRD1,TRD2,TRD3,TRG,TRD,AUX,AUY,AUZ,AUXMA,THETA,COEFDD REAL*8 XLONG,AG1,AG2,AD1,AD2,PSCAG1,PSCAG2,PSCAD1,PSCAD2, & COEF1,COEF2,COEF3,COEF4,SCN1X,SCN1Y,SCN1Z,VX,VY,VZ,COEF1X,COEF2X, & COEF1Y,COEF2Y,COEF1Z,COEF2Z,CX,CY,CZ,ANCX,ANCY,ANCZ, & DIFFX,DIFFY,DIFFZ,XLONGG,XLONGD & VALG,COEF,GX,GY,GZ,UN,EXPR1,EXPR2 c REAL*8 VECXG1(2),VECYG1(2) c REAL*8 VECXG2(2),VECYG2(2) c REAL*8 VECXD1(2),VECYD1(2) c REAL*8 VECXD2(2),VECYD2(2) REAL*8 VECXG(4,4),VECYG(4,4),VECZG(4,4) REAL*8 VECXD(4,4),VECYD(4,4),VECZD(4,4) REAL*8 VOLUG(4),SURFAGX(4),SURFAGY(4),SURFAGZ(4),COEFG(4) REAL*8 VOLUD(4),COEFD(4) REAL*8 PX(3,4),PY(3,4),PZ(3,4),TRGAUX(4) REAL*8 XPRO,TRACE INTEGER NGS(4),NLS(4) C INTEGER XS(4),YS(4),ZS(4) INTEGER NLOCFG(4,4),NLOCFD(4,4) INTEGER IGNS(4) REAL*8 EPS INTEGER ICRIT CHARACTER*8 TYPE C C CHARACTER*(4) NOMCOM(3) C DATA NOMCOM /'P1DX','P1DY','P1DZ'/ INTEGER NDIM SEGMENT MMAT1 REAL*8 PM(NDIM,NDIM),PM1(NDIM,NDIM),XSOL(NDIM) INTEGER IC(NDIM) ENDSEGMENT INTEGER K1,K2 SEGMENT INDICE INTEGER NUME(K1,K2) ENDSEGMENT POINTEUR IND.INDICE,IND2.INDICE,IND22.INDICE SEGMENT MATRICE REAL*8 MAT(K1,K2) ENDSEGMENT POINTEUR VAL1.MATRICE,VAL2.MATRICE,SCMB.MATRICE INTEGER K3 SEGMENT POINT2 INTEGER POINT(K3) ENDSEGMENT POINTEUR IPO2.POINT2 SEGMENT MATRICE2 REAL*8 MAT2(K1,K2) ENDSEGMENT POINTEUR MATR1.MATRICE2,MATR2.MATRICE2 SEGMENT POINT3 INTEGER POINT33(K3) ENDSEGMENT POINTEUR IPO3.POINT3 SEGMENT INDICE3 INTEGER NU(K1,K2) ENDSEGMENT POINTEUR INDIC.INDICE3 SEGMENT REP INTEGER ID(K3) ENDSEGMENT POINTEUR TAB.REP,INDLI.REP INTEGER K5 SEGMENT NBFAC INTEGER NBFACEL(K5) INTEGER IMELEM(K5) ENDSEGMENT INTEGER K6 SEGMENT NBCOT INTEGER NBCOTE(K6) INTEGER IMECOTE(K6) ENDSEGMENT INTEGER K7,K8 SEGMENT MISZERO INTEGER TABL(K7) INTEGER TABL2(K7) INTEGER TABL1(K8),IPOS(K8),ICOURANT(K8) REAL*8 XMAX(K7) ENDSEGMENT POINTEUR ITAB.MISZERO c WRITE(6,*) 'COUCOU NORV4' NMAI1 = 0 NMAI2 = 0 NMAI3 = 0 C NMAI4 = 0 C Initialisation sinon NAN... DO II=1,4 SURFAGX(II)=0.D0 SURFAGY(II)=0.D0 SURFAGZ(II)=0.D0 ENDDO NBSO = MAX(1,MELTFA.LISOUS(/1)) c WRITE(6,*) 'NBSO MAILLE= ',NBSO c WRITE(6,*) 'MELTFA= ',MELTFA C IELTFA = MELTFA IF (NBSO.EQ.1) THEN K5 = MELTFA.NUM(/2) ELSEIF (NBSO.EQ.2) THEN IPT1 = MELTFA.LISOUS(1) SEGACT IPT1 N1 = IPT1.NUM(/2) NMAI1 = N1 SEGDES IPT1 IPT2 = MELTFA.LISOUS(2) SEGACT IPT2 N2 = IPT2.NUM(/2) NMAI2 = N2 SEGDES IPT2 K5 = N1 + N2 ELSEIF (NBSO.EQ.3) THEN IPT1 = MELTFA.LISOUS(1) SEGACT IPT1 N1 = IPT1.NUM(/2) NMAI1 = N1 SEGDES IPT1 c WRITE(6,*) 'N1= ',N1 IPT2 = MELTFA.LISOUS(2) SEGACT IPT2 N2 = IPT2.NUM(/2) NMAI2 = N2 SEGDES IPT2 c WRITE(6,*) 'N2= ',N2 IPT3 = MELTFA.LISOUS(3) SEGACT IPT3 N3 = IPT3.NUM(/2) NMAI3 = N3 c WRITE(6,*) 'N3= ',N3 SEGDES IPT3 K5 = N1 + N2 + N3 ELSEIF (NBSO.EQ.4) THEN IPT1 = MELTFA.LISOUS(1) SEGACT IPT1 N1 = IPT1.NUM(/2) NMAI1 = N1 SEGDES IPT1 IPT2 = MELTFA.LISOUS(2) SEGACT IPT2 N2 = IPT2.NUM(/2) NMAI2 = N2 SEGDES IPT2 IPT3 = MELTFA.LISOUS(3) SEGACT IPT3 N3 = IPT3.NUM(/2) NMAI3 = N3 SEGDES IPT3 IPT4 = MELTFA.LISOUS(4) SEGACT IPT4 N4 = IPT4.NUM(/2) C NMAI4 = N4 SEGDES IPT4 K5 = N1 + N2 + N3 + N4 ENDIF c WRITE(6,*) 'K5= ',K5 IF (NBSO.EQ.1) THEN DO I = 1,K5 NTYPE = MELTFA.ITYPEL c WRITE(6,*) 'NTYPE= ',NTYPE IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = MELTFA ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = MELTFA ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = MELTFA ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = MELTFA ENDIF c SEGDES MELTFA ENDDO ELSEIF (NBSO.EQ.2) THEN IPT1 = MELTFA.LISOUS(1) SEGACT IPT1 IPT2 = MELTFA.LISOUS(2) SEGACT IPT2 DO I = 1,K5 N1 = IPT1.NUM(/2) IF (I.LE.N1) THEN NTYPE = IPT1.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT1 ENDIF ELSE NTYPE = IPT2.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT2 ENDIF ENDIF ENDDO ELSEIF (NBSO.EQ.3) THEN c WRITE(6,*) 'COUCOU' IPT1 = MELTFA.LISOUS(1) SEGACT IPT1 NTYPE = IPT1.ITYPEL c WRITE(6,*) 'NTYPE= ',IPT1.ITYPEL IPT2 = MELTFA.LISOUS(2) SEGACT IPT2 NTYPE = IPT2.ITYPEL c WRITE(6,*) 'NTYPE= ',IPT2.ITYPEL IPT3 = MELTFA.LISOUS(3) SEGACT IPT3 NTYPE = IPT3.ITYPEL c WRITE(6,*) 'NTYPE= ',IPT3.ITYPEL N1 = IPT1.NUM(/2) N2 = IPT2.NUM(/2) N3 = IPT3.NUM(/2) DO I = 1,K5 IF (I.LE.N1) THEN NTYPE = IPT1.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT1 ENDIF ELSEIF (I.LE.(N1+N2)) THEN NTYPE = IPT2.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT2 ENDIF ELSE NTYPE = IPT3.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT3 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT3 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT3 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT3 ENDIF ENDIF ENDDO ELSEIF (NBSO.EQ.4) THEN IPT1 = MELTFA.LISOUS(1) SEGACT IPT1 NTYPE = IPT1.ITYPEL c WRITE(6,*) 'NTYPE= ',IPT1.ITYPEL IPT2 = MELTFA.LISOUS(2) SEGACT IPT2 NTYPE = IPT2.ITYPEL c WRITE(6,*) 'NTYPE= ',IPT2.ITYPEL IPT3 = MELTFA.LISOUS(3) SEGACT IPT3 NTYPE = IPT3.ITYPEL c WRITE(6,*) 'NTYPE= ',IPT3.ITYPEL IPT4 = MELTFA.LISOUS(4) SEGACT IPT4 NTYPE = IPT4.ITYPEL c WRITE(6,*) 'NTYPE= ',IPT4.ITYPEL N1 = IPT1.NUM(/2) N2 = IPT2.NUM(/2) N3 = IPT3.NUM(/2) N4 = IPT4.NUM(/2) DO I = 1,K5 IF (I.LE.N1) THEN NTYPE = IPT1.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT1 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT1 ENDIF ELSEIF (I.LE.(N1+N2)) THEN NTYPE = IPT2.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT2 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT2 ENDIF ELSEIF (I.LE.(N1+N2+N3)) THEN NTYPE = IPT3.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT3 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT3 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT3 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT3 ENDIF ELSE NTYPE = IPT4.ITYPEL IF (NTYPE .EQ. 16) THEN NBFACEL(I) = 6 IMELEM(I) = IPT4 ELSEIF (NTYPE .EQ. 25) THEN NBFACEL(I) = 5 IMELEM(I) = IPT4 ELSEIF (NTYPE .EQ. 23) THEN NBFACEL(I) = 4 IMELEM(I) = IPT4 ELSEIF (NTYPE .EQ. 9) THEN NBFACEL(I) = 5 IMELEM(I) = IPT4 ENDIF ENDIF ENDDO ENDIF c CAS OU LES FACES SONT DES TRIANGLES OU DES FACES C MAUX = MELEFP NFAI1 = 0 NBSOF = MAX(1,MELEFP.LISOUS(/1)) c WRITE(6,*) 'NBSOF FACE= ',NBSOF C IELTFA = MELTFA IF (NBSOF.EQ.1) THEN K6 = MELEFP.NUM(/2) ELSEIF (NBSOF.EQ.2) THEN IPT5 = MELEFP.LISOUS(1) SEGACT IPT5 N1 = IPT5.NUM(/2) NFAI1 = N1 SEGDES IPT5 IPT6 = MELEFP.LISOUS(2) SEGACT IPT6 N2 = IPT6.NUM(/2) C NFAI2 = N2 SEGDES IPT6 K6 = N1 + N2 ENDIF C ON EST ICI IF (NBSOF.EQ.1) THEN DO I = 1,K6 NTYPE = MELEFP.ITYPEL c WRITE(6,*) 'NTYPE= ',NTYPE IF (NTYPE .EQ. 5) THEN NBCOTE(I) = 3 IMECOTE(I) = MELEFP ELSE NBCOTE(I) = 4 IMECOTE(I) = MELEFP ENDIF c SEGDES MELTFA ENDDO ELSEIF (NBSOF.EQ.2) THEN c WRITE(6,*) 'POINT2' IPT5 = MELEFP.LISOUS(1) SEGACT IPT5 IPT6 = MELEFP.LISOUS(2) SEGACT IPT6 c WRITE(6,*) 'IPT5= ',IPT5.ITYPEL c WRITE(6,*) 'IPT6= ',IPT6.ITYPEL DO I = 1,K6 N1 = IPT5.NUM(/2) C MISE A JOUR DE MLEFA.LECT IF (I.LE.N1) THEN N0 = IPT5.NUM(/1) NGFAUX = IPT5.NUM(N0,I) MLEFA2.LECT(NGFAUX) = I c WRITE(6,*) 'NGFAUX = ',NGFAUX, c & 'MLEFA2=',MLEFA2.LECT(NGFAUX) IF (IPT5.ITYPEL .EQ. 5) THEN NBCOTE(I) = 3 IMECOTE(I) = IPT5 ELSE NBCOTE(I) = 4 IMECOTE(I) = IPT5 ENDIF c SEGDES IPT5 ELSE N0 = IPT6.NUM(/1) NGFAUX = IPT6.NUM(N0,I-NFAI1) MLEFA2.LECT(NGFAUX) = I IF (IPT6.ITYPEL .EQ. 5) THEN NBCOTE(I) = 3 IMECOTE(I) = IPT6 ELSE NBCOTE(I) = 4 IMECOTE(I) = IPT6 ENDIF c SEGDES IPT6 ENDIF ENDDO ENDIF NFAC=MELEFL.NUM(/2) c ON MAJORE SUPERIEUREMENT NBNN : ON LE REAJUSTERA PAR LA SUITE NCON = ((2*NBMAX)) + 1 c NCON = ((3*NBMAX)/2) + 1 NBNN = NCON NESSAI = NCON NBNN = NESSAI c WRITE(6,*) 'NBMAX= ',NBMAX c WRITE(6,*) 'NBNN= ',NBNN c WRITE(6,*) 'NFAC= ',NFAC C DEFINITION DES PARAMETRES DU CHAMELEM DES COEFFICIENTS c INITIALISATION DU CHAMELEM N1=1 N2=1 c WRITE(6,*) 'N2= ',N2 N3=6 L1=8 SEGINI MCHELM C ICOEFF = MCHELM MCHELM.TITCHE='Gradient' MCHELM.IFOCHE=IFOUR C ISOUS=0 NBSOUS=0 NBREF=0 NBELEM = NFAC ISOUS=ISOUS+1 SEGINI MELEME C ITYPEL=32 -> 'POLY' ITYPEL=32 MCHELM.IMACHE(ISOUS)=MELEME SEGINI MCHAML MCHELM.ICHAML(ISOUS)=MCHAML MCHAML.NOMCHE(1)='SCAL' MCHAML.TYPCHE(1)='REAL*8 ' N1PTEL=NESSAI N1EL=NBELEM N2PTEL=0 N2EL=0 SEGINI MELVA1 MCHAML.IELVAL(1)=MELVA1 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C CALCUL DE LA VITESSE EN CHAQUE FACE C AA = 0.0 C BB = 0.0 C INDICE QUI COMPTE LES COEFFICIENTS POUR CHAQUE FACE NAUX2 = 0 NMOY = 0 DO NLCF= 1, NFAC, 1 NCON = 1 NGCF=MELEFL.NUM(2,NLCF) NGCG=MELEFL.NUM(1,NLCF) NGCD=MELEFL.NUM(3,NLCF) NLCG=MLECEN.LECT(NGCG) NLCD=MLECEN.LECT(NGCD) SCNX=MPONOR.VPOCHA(NLCF,1) SCNY=MPONOR.VPOCHA(NLCF,2) SCNZ=MPONOR.VPOCHA(NLCF,3) SCN1X = SCNX SCN1Y = SCNY SCN1Z = SCNZ SURF=0.5D0*MPOSUR.VPOCHA(NLCF,1) SCNX=SCNX*SURF SCNY=SCNY*SURF SCNZ=SCNZ*SURF C 4=IDIM+1 ICELL=(4*(NGCG -1))+1 XG=MCOORD.XCOOR(ICELL) YG=MCOORD.XCOOR(ICELL+1) ZG=MCOORD.XCOOR(ICELL+2) ICELL=(4*(NGCD -1))+1 XD=MCOORD.XCOOR(ICELL) YD=MCOORD.XCOOR(ICELL+1) ZD=MCOORD.XCOOR(ICELL+2) ICELL=(4*(NGCF -1))+1 XF=MCOORD.XCOOR(ICELL) YF=MCOORD.XCOOR(ICELL+1) ZF=MCOORD.XCOOR(ICELL+2) MELEME.NUM(1,NLCF)=NGCF MELVA1.VELCHE(1,NLCF)=0.0D0 DO J= 1,NBNN MELVA1.VELCHE(J,NLCF) = 0.0D0 ENDDO C MISE A ZERO DE NLOC DO JA=1,4 DO IA=1,3 NLOCFG(IA,JA) = 0 NLOCFD(IA,JA) = 0 ENDDO ENDDO MELTFA = IMELEM(NLCG) NBF = NBFACEL(NLCG) IF (NLCG.LE.NMAI1) THEN NGAUX = NLCG ELSEIF ((NLCG.GT.NMAI1).AND.(NLCG.LE.(NMAI1+NMAI2))) THEN NGAUX = NLCG - NMAI1 ELSEIF ((NLCG.GT.(NMAI1+NMAI2)).AND. & (NLCG.LE.(NMAI1+NMAI2+NMAI3))) THEN NGAUX = NLCG - (NMAI1+NMAI2) ELSEIF (NLCG.GT.(NMAI1+NMAI2+NMAI3)) THEN NGAUX = NLCG - (NMAI1+NMAI2+NMAI3) ENDIF c ON TIENT COMPTE DU CHANGEMENT DE NUMEROTATION NLCF1 = MLEFA2.LECT(NGCF) NBNO = NBCOTE(NLCF1) MELEFP = IMECOTE(NLCF1) IF (NLCF1.GT.NFAI1) THEN NLCF1AUX = NLCF1 - NFAI1 ELSE NLCF1AUX = NLCF1 ENDIF DO IA=1,NBNO NGS(IA) = MELEFP.NUM(IA,NLCF1AUX) NLS(IA) = MLESOM.LECT(NGS(IA)) C ICELL=(4*(NGS(IA) -1))+1 C XS(IA)=nint(MCOORD.XCOOR(ICELL)) C YS(IA)=nint(MCOORD.XCOOR(ICELL+1)) C ZS(IA)=nint(MCOORD.XCOOR(ICELL+2)) ENDDO C ON REPERE LES VECTEURS PRINCIPAUX DE LA BASE DO JA = 1,NBNO NGS(JA) = MELEFP.NUM(JA,NLCF1AUX) c WRITE(6,*) 'NGAUX= ',NGAUX,'JA= ',JA,'NGS= ',NGS(JA) c WRITE(6,*) 'NGCF= ',NGCF,'NLCF= ',NLCF ICOUR = 0 DO J = 1,NBF N1 = MELTFA.NUM(J,NGAUX) NL1 = MLEFA2.LECT(N1) NBNO2 = NBCOTE(NL1) MELEP2 = IMECOTE(NL1) IF (NL1.GT.NFAI1) THEN NL1AUX = NL1 - NFAI1 ELSE NL1AUX = NL1 ENDIF c WRITE(6,*) 'N1= ',N1,'NL1= ',NL1,'NL1AUX= ',NL1AUX DO IA =1,NBNO2 NSOM1 = MELEP2.NUM(IA,NL1AUX) c WRITE(6,*) 'NBNO2= ',NBNO2,'IA= ',IA,'NSOM1= ',NSOM1 IF (NSOM1.EQ.NGS(JA)) THEN ICELL=(4*(N1 -1))+1 XF=MCOORD.XCOOR(ICELL) YF=MCOORD.XCOOR(ICELL+1) ZF=MCOORD.XCOOR(ICELL+2) ICOUR = ICOUR + 1 VECXG(ICOUR,JA) = (XF - XG) VECYG(ICOUR,JA) = (YF - YG) VECZG(ICOUR,JA) = (ZF - ZG) NLOCFG(ICOUR,JA) = N1 C ON PERMUTE IF (N1.EQ.NGCF) THEN NAUX = NLOCFG(1,JA) VXAU = VECXG(1,JA) VYAU = VECYG(1,JA) VZAU = VECZG(1,JA) VECXG(1,JA) = (XF - XG) VECYG(1,JA) = (YF - YG) VECZG(1,JA) = (ZF - ZG) NLOCFG(1,JA) = N1 VECXG(ICOUR,JA) = VXAU VECYG(ICOUR,JA) = VYAU VECZG(ICOUR,JA) = VZAU NLOCFG(ICOUR,JA) = NAUX ENDIF ENDIF ENDDO ENDDO ENDDO MELTFA = IMELEM(NLCD) NBF = NBFACEL(NLCD) IF (NLCD.LE.NMAI1) THEN NDAUX = NLCD ELSEIF ((NLCD.GT.NMAI1).AND.(NLCD.LE.(NMAI1+NMAI2))) THEN NDAUX = NLCD - NMAI1 ELSEIF ((NLCD.GT.(NMAI1+NMAI2)).AND. & (NLCD.LE.(NMAI1+NMAI2+NMAI3))) THEN NDAUX = NLCD - (NMAI1+NMAI2) ELSEIF (NLCD.GT.(NMAI1+NMAI2+NMAI3)) THEN NDAUX = NLCD - (NMAI1+NMAI2+NMAI3) ENDIF C ON REPERE LES VECTEURS PRINCIPAUX DE LA BASE DO JA = 1,NBNO NGS(JA) = MELEFP.NUM(JA,NLCF1AUX) c WRITE(6,*) 'NDAUX= ',NDAUX,'JA= ',JA,'NGS= ',NGS(JA) c WRITE(6,*) 'NGCF= ',NGCF,'NLCF= ',NLCF ICOUR = 0 DO J = 1,NBF N1 = MELTFA.NUM(J,NDAUX) NL1 = MLEFA2.LECT(N1) c WRITE(6,*) 'N1= ',N1,'NL1= ',NL1 NBNO2 = NBCOTE(NL1) MELEP2 = IMECOTE(NL1) IF (NL1.GT.NFAI1) THEN NL1AUX = NL1 - NFAI1 ELSE NL1AUX = NL1 ENDIF DO IA =1,NBNO2 NSOM1 = MELEP2.NUM(IA,NL1AUX) c WRITE(6,*) 'NBNO2= ',NBNO2,'IA= ',IA,'NSOM1= ',NSOM1 IF (NSOM1.EQ.NGS(JA)) THEN ICELL=(4*(N1 -1))+1 XF=MCOORD.XCOOR(ICELL) YF=MCOORD.XCOOR(ICELL+1) ZF=MCOORD.XCOOR(ICELL+2) ICOUR = ICOUR + 1 VECXD(ICOUR,JA) = (XF - XD) VECYD(ICOUR,JA) = (YF - YD) VECZD(ICOUR,JA) = (ZF - ZD) NLOCFD(ICOUR,JA) = N1 C ON PERMUTE IF (N1.EQ.NGCF) THEN NAUX = NLOCFD(1,JA) VXAU = VECXD(1,JA) VYAU = VECYD(1,JA) VZAU = VECZD(1,JA) VECXD(1,JA) = (XF - XD) VECYD(1,JA) = (YF - YD) VECZD(1,JA) = (ZF - ZD) NLOCFD(1,JA) = N1 VECXD(ICOUR,JA) = VXAU VECYD(ICOUR,JA) = VYAU VECZD(ICOUR,JA) = VZAU NLOCFD(ICOUR,JA) = NAUX ENDIF ENDIF ENDDO ENDDO c WRITE(6,*) 'JA= ',JA c WRITE(6,*) 'ICOUR= ',ICOUR ENDDO MPOGRA.VPOCHA(NLCF,1) = 0.D0 DO JA = 1,NBNO c WRITE(6,*) 'XPRO= ',XPRO DO KA = 1,ICOUR C PRODUIT MIXTES C PRODUIT VECTORIEL IF (KA.EQ.1) THEN PSCAGX = (VECYG(2,JA)*VECZG(3,JA)) - & (VECZG(2,JA)*VECYG(3,JA)) PSCAGY = (VECZG(2,JA)*VECXG(3,JA)) - & (VECXG(2,JA)*VECZG(3,JA)) PSCAGZ = (VECXG(2,JA)*VECYG(3,JA)) - & (VECYG(2,JA)*VECXG(3,JA)) VOLUG(JA) = (VECXG(1,JA)* PSCAGX) + & (VECYG(1,JA)* PSCAGY) + & (VECZG(1,JA)* PSCAGZ) SURFAGX(KA) = 0.5D0* PSCAGX SURFAGY(KA) = 0.5D0* PSCAGY SURFAGZ(KA) = 0.5D0* PSCAGZ IF ( VOLUG(JA).GT.0) THEN SURFAGX(KA) = - SURFAGX(KA) SURFAGY(KA) = - SURFAGY(KA) SURFAGZ(KA) = - SURFAGZ(KA) ENDIF VOLUG(JA) = 1.D0/6.D0*ABS(VOLUG(JA)) ENDIF IF (KA.EQ.2) THEN PSCAGX = (VECYG(3,JA)*VECZG(1,JA)) - & (VECZG(3,JA)*VECYG(1,JA)) PSCAGY = (VECZG(3,JA)*VECXG(1,JA)) - & (VECXG(3,JA)*VECZG(1,JA)) PSCAGZ = (VECXG(3,JA)*VECYG(1,JA)) - & (VECYG(3,JA)*VECXG(1,JA)) SURFAGX(KA) = 0.5D0* PSCAGX SURFAGY(KA) = 0.5D0* PSCAGY SURFAGZ(KA) = 0.5D0* PSCAGZ PSCA = (VECXG(2,JA)* PSCAGX) + (VECYG(2,JA)* PSCAGY) + & (VECZG(2,JA)* PSCAGZ) IF ( PSCA.GT.0) THEN SURFAGX(KA) = - SURFAGX(KA) SURFAGY(KA) = - SURFAGY(KA) SURFAGZ(KA) = - SURFAGZ(KA) ENDIF ENDIF IF (KA.EQ.3) THEN PSCAGX = (VECYG(1,JA)*VECZG(2,JA)) - & (VECZG(1,JA)*VECYG(2,JA)) PSCAGY = (VECZG(1,JA)*VECXG(2,JA)) - & (VECXG(1,JA)*VECZG(2,JA)) PSCAGZ = (VECXG(1,JA)*VECYG(2,JA)) - & (VECYG(1,JA)*VECXG(2,JA)) SURFAGX(KA) = 0.5D0* PSCAGX SURFAGY(KA) = 0.5D0* PSCAGY SURFAGZ(KA) = 0.5D0* PSCAGZ PSCA = (VECXG(3,JA)* PSCAGX) + (VECYG(3,JA)* PSCAGY) + & (VECZG(3,JA)* PSCAGZ) IF ( PSCA.GT.0) THEN SURFAGX(KA) = - SURFAGX(KA) SURFAGY(KA) = - SURFAGY(KA) SURFAGZ(KA) = - SURFAGZ(KA) ENDIF ENDIF c CALCUL DE MATRICE POUR LE NOEUD D INDICE NS1 IF (ICHTE.EQ.0) THEN PX(KA,JA) = SURFAGX(KA)/(3.D0*VOLUG(JA)) PY(KA,JA) = SURFAGY(KA)/(3.D0*VOLUG(JA)) PZ(KA,JA) = SURFAGZ(KA)/(3.D0*VOLUG(JA)) ELSE IF (MPOTEN.VPOCHA(/2) .EQ.6) THEN C TENSEUR ANISOTROPE K11G = MPOTEN.VPOCHA(NLCG,1) K22G = MPOTEN.VPOCHA(NLCG,2) K33G = MPOTEN.VPOCHA(NLCG,3) K21G = MPOTEN.VPOCHA(NLCG,4) K31G = MPOTEN.VPOCHA(NLCG,5) K32G = MPOTEN.VPOCHA(NLCG,6) ELSEIF (MPOTEN.VPOCHA(/2) .EQ.1) THEN C TENSEUR DIAGONAL K11G = MPOTEN.VPOCHA(NLCG,1) K22G = K11G K33G = K11G K21G = 0.0D0 K31G = 0.0D0 K32G = 0.0D0 ELSE WRITE(6,*) 'TENSEUR NON PREVU' STOP ENDIF PSCAGX = (K11G*SURFAGX(KA)) + (K21G*SURFAGY(KA)) + & (K31G*SURFAGZ(KA)) PSCAGY = (K21G*SURFAGX(KA)) + (K22G*SURFAGY(KA)) + & (K32G*SURFAGZ(KA)) PSCAGZ = (K31G*SURFAGX(KA)) + (K32G*SURFAGY(KA)) + & (K33G*SURFAGZ(KA)) PX(KA,JA) = PSCAGX/(3.D0*VOLUG(JA)) PY(KA,JA) = PSCAGY/(3.D0*VOLUG(JA)) PZ(KA,JA) = PSCAGZ/(3.D0*VOLUG(JA)) ENDIF ENDDO ENDDO XPRO = 1.D0/NBNO DO JA = 1,NBNO C MARQ = 0 DO I5 = 1,INDLI.ID(NLS(JA)) INDAUX = IND2.NUME(I5,NLS(JA)) IF (INDAUX.EQ.NLOCFG(2,JA)) THEN IAFF = I5 IG2 = IAFF GOTO 62 ENDIF ENDDO 62 CONTINUE TRG2 = SCMB.MAT(IAFF,NLS(JA)) c WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG2 DO I5 = 1,INDLI.ID(NLS(JA)) INDAUX = IND2.NUME(I5,NLS(JA)) IF (INDAUX.EQ.NLOCFG(3,JA)) THEN IAFF = I5 IG3 = IAFF GOTO 629 ENDIF ENDDO 629 CONTINUE TRG3 = SCMB.MAT(IAFF,NLS(JA)) c WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG3 C MARQ = 0 DO I5 = 1,INDLI.ID(NLS(JA)) INDAUX = IND2.NUME(I5,NLS(JA)) IF (INDAUX.EQ.NGCF) THEN IAFF = I5 IG = IAFF GOTO 63 ENDIF ENDDO 63 CONTINUE TRG = SCMB.MAT(IAFF,NLS(JA)) TRGAUX(JA) = TRG IGNS(JA) = IAFF c WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG VAL = MPOCHP.VPOCHA(NLCG,1) C VALD = MPOCHP.VPOCHA(NLCD,1) c ICI AUX = (XPRO*( & (PX(1,JA) * (VAL - TRG)) & + (PX(2,JA) * (VAL - TRG2)) & + (PX(3,JA) * (VAL - TRG3)))) AUY = (XPRO*( & (PY(1,JA) * (VAL - TRG)) & + (PY(2,JA) * (VAL - TRG2)) & + (PY(3,JA) * (VAL - TRG3)))) AUZ = (XPRO*( & (PZ(1,JA) * (VAL - TRG)) & + (PZ(2,JA) * (VAL - TRG2)) & + (PZ(3,JA) * (VAL - TRG3)))) MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) + & (AUX*SCN1X) + (AUY*SCN1Y) + & (AUZ*SCN1Z) c IF (NLCF.EQ.791) THEN c WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1) c WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2) c WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3) c WRITE(6,*) 'PX= ',PX(1,JA),PX(2,JA),PX(3,JA) c WRITE(6,*) 'PY= ',PY(1,JA),PY(2,JA),PY(3,JA) c WRITE(6,*) 'PZ= ',PZ(1,JA),PZ(2,JA),PZ(3,JA) c WRITE(6,*) 'TR ',TRG,TRG2,TRG3,'VAL',VAL c WRITE(6,*) 'NLS= ',NLS(JA), c & 'IG= ',IG,'IG2= ',IG2,'IG3= ',IG3 c WRITE(6,*) 'NGCF= ',NGCF c ENDIF C ITROUVE = 0 INDIC = IPO3.POINT33(NLS(JA)) SEGACT INDIC *MOD MATR1 = IPO2.POINT(NLS(JA)) SEGACT MATR1 *MOD c NLS1 = NLS(JA) c DO IAUX = 1,INDLI.ID(NLS1) c DO IAUX2 = 1,TAB.ID(NLS1) c WRITE(6,*) 'NLS1= ',NLS1,'IAUX= ',IAUX ,'IAUX2= ', c & IAUX2,'VAUX',MATR1.MAT2(IAUX,IAUX2) c & ,'IND3= ',INDIC.NU(IAUX,IAUX2) c ENDDO c ENDDO DO ICOUR = 1,TAB.ID(NLS(JA)) IA = ICOUR J1 = INDIC.NU(IG,IA) DO IAUX2 = 2,NCON J2 = MELEME.NUM(IAUX2,NLCF) IF (J1.EQ.J2) THEN IAUX = IAUX2 C ITROUVE = 1 GOTO 5119 ENDIF ENDDO C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR NCON = NCON + 1 IAUX = NCON 5119 CONTINUE CX = MATR1.MAT2(IG,IA) CY = MATR1.MAT2(IG,IA) CZ = MATR1.MAT2(IG,IA) IF (J1.EQ.NGCG) THEN CX = CX - 1.D0 CY = CY - 1.D0 CZ = CZ - 1.D0 INDGA = IAUX ENDIF IF (J1.EQ.NGCD) THEN INDDR = IAUX ENDIF c MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - c & (XPRO*PX(1,JA)*CX) c MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) - c & (XPRO*PY(1,JA)*CY) c MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) - c & (XPRO*PZ(1,JA)*CZ) AUX = (XPRO*PX(1,JA)*CX*SCN1X) + & (XPRO*PY(1,JA)*CY*SCN1Y) + & (XPRO*PZ(1,JA)*CZ*SCN1Z) MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - & AUX INDAUX = INDIC.NU(IG,IA) MELEME.NUM(IAUX,NLCF) = INDAUX ENDDO C ITROUVE = 0 DO ICOUR = 1,TAB.ID(NLS(JA)) IA = ICOUR J1 = INDIC.NU(IG2,IA) DO IAUX2 = 2,NCON J2 = MELEME.NUM(IAUX2,NLCF) IF (J1.EQ.J2) THEN IAUX = IAUX2 C ITROUVE = 1 GOTO 511 ENDIF ENDDO C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR NCON = NCON + 1 IAUX = NCON 511 CONTINUE CX = MATR1.MAT2(IG2,IA) CY = MATR1.MAT2(IG2,IA) CZ = MATR1.MAT2(IG2,IA) IF (J1.EQ.NGCG) THEN CX = CX - 1.D0 CY = CY - 1.D0 CZ = CZ - 1.D0 ENDIF c MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - c & (XPRO*PX(2,JA)*CX) c MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) - c & (XPRO*PY(2,JA)*CY) c MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) - c & (XPRO*PZ(2,JA)*CZ) AUX = (XPRO*PX(2,JA)*CX*SCN1X) + & (XPRO*PY(2,JA)*CY*SCN1Y) + & (XPRO*PZ(2,JA)*CZ*SCN1Z) MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - & AUX INDAUX = INDIC.NU(IG2,IA) MELEME.NUM(IAUX,NLCF) = INDAUX ENDDO C ITROUVE = 0 DO ICOUR = 1,TAB.ID(NLS(JA)) IA = ICOUR J1 = INDIC.NU(IG3,IA) DO IAUX2 = 2,NCON J2 = MELEME.NUM(IAUX2,NLCF) IF (J1.EQ.J2) THEN IAUX = IAUX2 C ITROUVE = 1 GOTO 5118 ENDIF ENDDO C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR NCON = NCON + 1 IAUX = NCON 5118 CONTINUE CX = MATR1.MAT2(IG3,IA) CY = MATR1.MAT2(IG3,IA) CZ = MATR1.MAT2(IG3,IA) IF (J1.EQ.NGCG) THEN CX = CX - 1.D0 CY = CY - 1.D0 CZ = CZ - 1.D0 ENDIF c MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - c & (XPRO*PX(3,JA)*CX) c MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) - c & (XPRO*PY(3,JA)*CY) c MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) - c & (XPRO*PZ(3,JA)*CZ) AUX = (XPRO*PX(3,JA)*CX*SCN1X) + & (XPRO*PY(3,JA)*CY*SCN1Y) + & (XPRO*PZ(3,JA)*CZ*SCN1Z) MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - & AUX INDAUX = INDIC.NU(IG3,IA) MELEME.NUM(IAUX,NLCF) = INDAUX ENDDO c WRITE(6,*) 'NLCF= ',NLCF,'JA= ',JA,'VOLUG(JA) = ',VOLUG(JA), c & 'XPRO= ',XPRO c FIN DE LA BOUCLE SUR LES NOEUDS SEGDES MATR1 SEGDES INDIC ENDDO c MPOGRA.VPOCHA(NLCF,1) = 0.0D0 c ISUP = NCON c DO J= ISUP+1,NBNN c DO J= 2,NBNN c MELVA1.VELCHE(J,NLCF) = 0.0D0 c MELEME.NUM(J,NLCF) = MELEME.NUM(ISUP,NLCF) c ENDDO ISUP = NCON DO J= ISUP+1,NBNN MELVA1.VELCHE(J,NLCF) = 0.0D0 MELEME.NUM(J,NLCF) = MELEME.NUM(ISUP,NLCF) ENDDO CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C ON RAJOUTE LE CONVECTIF C CALCUL PLUS PRECIS c IF (NLCD.NE.NLCG) THEN c MELVA1.VELCHE(INDDR,NLCF) = - MELVA1.VELCHE(INDGA,NLCF) c ELSE c DO J= 1,ISUP c ICENT = MELEME.NUM(J,NLCF) c ICEN = MLECEN.LECT(ICENT) c DIFF = MELVA1.VELCHE(J,NLCF) + MELVA1.VELCHE(INDGA,NLCF) c DIFF = ABS(DIFF) c XX = ABS(MELVA1.VELCHE(INDGA,NLCF)) c IF ((ICEN.EQ.0).AND.(DIFF.LT.(1e-5*XX))) THEN c MELVA1.VELCHE(J,NLCF) = - MELVA1.VELCHE(INDGA,NLCF) c ENDIF c ENDDO c ENDIF IF (ICHCO.GT.0) THEN C BOUCLE POUR CALCUER INDGA,INDDR INDFR = 0 DO J= 1,ISUP IF (MELEME.NUM(J,NLCF).EQ.NGCG) INDGA = J IF (MELEME.NUM(J,NLCF).EQ.NGCD) INDDR = J IF (MELEME.NUM(J,NLCF).EQ.NGCF) INDFR = J ENDDO c WRITE(6,*) 'NGCF= ',NGCF UN = MPOVCO.VPOCHA(NLCF,1) c WRITE(6,*) 'UN= ',UN C OPTION CENTRE IF (IOP.EQ.2) THEN IF (NLCD.NE.NLCG) THEN VAL = 0.5D0*(MPOCHP.VPOCHA(NLCG,1) + & MPOCHP.VPOCHA(NLCD,1))*UN MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL ELSE C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION XPRO = 1.D0/NBNO TRACE = 0.0D0 DO JA = 1,NBNO TRACE =TRACE + (XPRO*TRGAUX(JA)) ENDDO c WRITE(6,*) 'NLCF= ',NLCF,'TRACE= ',TRACE VAL = TRACE*UN MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL ENDIF C ON COMPLETE MELVA1 POUR LE CONVECTIF IF (NLCD.NE.NLCG) THEN MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - & (0.5D0*UN) MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) - & (0.5D0*UN) C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDENCES DES TRACES c POUR LES CONDITIONS MIXTES OU DE NEUMAN ELSE NLFCL = MLENCL.LECT(NGCF) C ON RAJOUTE CECI POUR L OPTION GRADGEO IF (NLFCL.NE.0) THEN MELVA1.VELCHE(NCON+1,NLCF) = - UN MELEME.NUM(NCON+1,NLCF) = NGCF ENDIF IF (NLFCL.EQ.0) THEN XPRO = 1.D0/NBNO DO JA = 1,NBNO INDIC = IPO3.POINT33(NLS(JA)) SEGACT INDIC *MOD MATR1 = IPO2.POINT(NLS(JA)) SEGACT MATR1 *MOD DO ICOUR = 1,TAB.ID(NLS(JA)) IA = ICOUR J1 = INDIC.NU(IGNS(JA),IA) DO IAUX2 = 2,NCON J2 = MELEME.NUM(IAUX2,NLCF) IF (J1.EQ.J2) THEN IAUX = IAUX2 GOTO 5169 ENDIF ENDDO 5169 CONTINUE CX = MATR1.MAT2(IGNS(JA),IA) MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - & (UN*CX*XPRO) ENDDO SEGDES INDIC *MOD SEGDES MATR1 *MOD ENDDO ENDIF ENDIF C OPTION UPWIND ELSEIF (IOP.EQ.1) THEN IF (NLCD.NE.NLCG) THEN IF (UN.GT.0.0D0) THEN VAL = (MPOCHP.VPOCHA(NLCG,1)*UN) ELSE VAL = (MPOCHP.VPOCHA(NLCD,1)*UN) ENDIF c WRITE(6,*) 'VAL= ',VAL MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL ELSE C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION c ANCIENNE VERRUE IF (UN.GT.0.0D0) THEN VAL = (MPOCHP.VPOCHA(NLCG,1)*UN) ELSE XPRO = 1.D0/NBNO TRACE = 0.0D0 DO JA = 1,NBNO TRACE =TRACE + (XPRO*TRGAUX(JA)) ENDDO c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF c WRITE(6,*) 'UN= ',UN,'TRACE=',TRACE VAL = TRACE*UN ENDIF MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL ENDIF C ON COMPLETE MELVA1 POUR LE CONVECTIF IF (NLCD.NE.NLCG) THEN c WRITE(6,*) 'UN= ',UN IF (UN.GT.0.0D0) THEN MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - & (UN) ELSE MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) - & (UN) ENDIF C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDENCES DES TRACES c POUR LES CONDITIONS MIXTES OU DE NEUMAN ELSE c ANCIENNE VERRUE IF (UN.GT.0.0D0) THEN MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - & (UN) ELSE NLFCL = MLENCL.LECT(NGCF) C ON RAJOUTE CECI POUR L OPTION GRADGEO C IL Y A PROBABLEMENT REDONDANCE IF (NLFCL.NE.0) THEN MELVA1.VELCHE(NCON+1,NLCF) = - UN MELEME.NUM(NCON+1,NLCF) = NGCF ENDIF * ENDIF IF (NLFCL.EQ.0) THEN XPRO = 1.D0/NBNO DO JA = 1,NBNO INDIC = IPO3.POINT33(NLS(JA)) SEGACT INDIC *MOD MATR1 = IPO2.POINT(NLS(JA)) SEGACT MATR1 *MOD DO ICOUR = 1,TAB.ID(NLS(JA)) IA = ICOUR J1 = INDIC.NU(IGNS(JA),IA) DO IAUX2 = 2,NCON J2 = MELEME.NUM(IAUX2,NLCF) IF (J1.EQ.J2) THEN IAUX = IAUX2 GOTO 5129 ENDIF ENDDO 5129 CONTINUE CX = MATR1.MAT2(IGNS(JA),IA) MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - & (UN*CX*XPRO) ENDDO SEGDES INDIC *MOD SEGDES MATR1 *MOD ENDDO ENDIF ENDIF ENDIF C OPTION UPWICENT ELSEIF (IOP.EQ.3) THEN IF (NLCD.NE.NLCG) THEN c CALCUL DE THETA c DANS MELVA1 : INDDR EST LA NUM de NLCD C COEFD = MELVA1(INDDR,NLCF) COEFDD = MELVA1.VELCHE(INDDR,NLCF) THETA = 0.5D0 EXPR1 = - (THETA*UN) - COEFDD EXPR2 = - ((1.D0 - THETA)*UN) + COEFDD c WRITE(6,*) 'EXPR1',EXPR1 c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, c & 'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN c WRITE(6,*) 'EXPR2',EXPR2 c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, c & 'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN IF (EXPR1.GT.0.0D0) THEN IF (ABS(UN) .GT. 1e-20) THEN THETA = - (0.5D0*COEFDD / UN) ENDIF ENDIF IF (EXPR2.LT.0.0D0) THEN IF (ABS(UN) .GT. 1e-20) THEN THETA = 1.D0 - (0.5D0*COEFDD / UN) ENDIF ENDIF THETA = MIN(1.D0,THETA) THETA = MAX(0.D0,THETA) EXPR1 = (THETA*UN) + COEFDD EXPR2 = ((1.D0 - THETA)*UN) - COEFDD c IF (EXPR1.LT.0.0D0) THEN c WRITE(6,*) 'EXPR1',EXPR1 c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, c & 'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN c ENDIF c IF (EXPR2.GT.0.0D0) THEN c WRITE(6,*) 'EXPR2',EXPR2 c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, c & 'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN c ENDIF c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, c & 'EXPR1= ',EXPR1,'EXPR2= ',EXPR2 VAL = (MPOCHP.VPOCHA(NLCG,1)*UN*THETA) + & (MPOCHP.VPOCHA(NLCD,1)*(1.D0-THETA)*UN) MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL ELSE C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION IF (INDFR.NE.0) THEN COEFDD = MELVA1.VELCHE(INDFR,NLCF) ELSE COEFDD = 0.0D0 ENDIF C C'EST LE THETA LE MIEUX ADPATE POUR LES CONDITIONS AUX LIMITES THETA = 0.0D0 EXPR1 = - (THETA*UN) - COEFDD EXPR2 = - ((1.D0 - THETA)*UN) + COEFDD IF (EXPR1.GT.0.0D0) THEN IF (ABS(UN) .GT. 1e-20) THEN THETA = - (0.5D0*COEFDD / UN) ENDIF ENDIF IF (EXPR2.LT.0.0D0) THEN IF (ABS(UN) .GT. 1e-20) THEN THETA = 1.D0 - (0.5D0*COEFDD / UN) ENDIF ENDIF THETA = MIN(1.D0,THETA) THETA = MAX(0.D0,THETA) C ANCIENNE VERRUE c THETA = 0.0D0 XPRO = 1.D0/NBNO TRACE = 0.0D0 DO JA = 1,NBNO TRACE =TRACE + (XPRO*TRGAUX(JA)) ENDDO VAL = TRACE*UN VAL = (MPOCHP.VPOCHA(NLCG,1)*UN*THETA) + & (TRACE*(1.D0-THETA)*UN) MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, c & 'GR1= ',MPOGRA.VPOCHA(NLCF,1),'VAL = ',VAL, c & 'TRACE = ',TRACE,'UN= ',UN c MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, c & 'GR1= ',MPOGRA.VPOCHA(NLCF,1) ENDIF C ON COMPLETE MELVA1 POUR LE CONVECTIF IF (NLCD.NE.NLCG) THEN MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - & (THETA*UN) MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) - & ((1.D0-THETA)*UN) C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDANCES DES TRACES c POUR LES CONDITIONS MIXTES OU DE NEUMAN ELSE MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - & (THETA*UN) NLFCL = MLENCL.LECT(NGCF) C ON EST ICI : CORRIGER C ON RAJOUTE CECI POUR L OPTION GRADGEO c IF (NLFCL.NE.0) THEN c MELVA1.VELCHE(NCON+1,NLCF) = - UN c MELEME.NUM(NCON+1,NLCF) = NGCF c ENDIF c IF (NLFCL.EQ.0) THEN XPRO = 1.D0/NBNO DO JA = 1,NBNO INDIC = IPO3.POINT33(NLS(JA)) SEGACT INDIC *MOD MATR1 = IPO2.POINT(NLS(JA)) SEGACT MATR1 *MOD DO ICOUR = 1,TAB.ID(NLS(JA)) IA = ICOUR J1 = INDIC.NU(IGNS(JA),IA) DO IAUX2 = 2,NCON J2 = MELEME.NUM(IAUX2,NLCF) IF (J1.EQ.J2) THEN IAUX = IAUX2 GOTO 5159 ENDIF ENDDO 5159 CONTINUE CX = MATR1.MAT2(IGNS(JA),IA) MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - & ((1.D0-THETA)*UN*CX*XPRO) ENDDO SEGDES INDIC *MOD SEGDES MATR1 *MOD ENDDO c ENDIF ENDIF ENDIF ENDIF c DO J= 1,NBNN c PSCA = (MPOGRA.VPOCHA(NLCF,1)*SCN1X) + c & (MPOGRA.VPOCHA(NLCF,2)*SCN1Y) + c & (MPOGRA.VPOCHA(NLCF,3)*SCN1Z) c MPOGRA.VPOCHA(NLCF,1) = PSCA*SCN1X c MPOGRA.VPOCHA(NLCF,2) = PSCA*SCN1Y c MPOGRA.VPOCHA(NLCF,3) = PSCA*SCN1Z c ENDDO c IF (NLCF.EQ.791) THEN c WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1) c ENDIF c WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2) c WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3) c DO J= 1,NBNN c PSCA = (MELVA1.VELCHE(J,NLCF)*SCN1X) + c & (MELVA2.VELCHE(J,NLCF)*SCN1Y) + c & (MELVA3.VELCHE(J,NLCF)*SCN1Z) c MELVA1.VELCHE(J,NLCF) = PSCA*SCN1X c MELVA2.VELCHE(J,NLCF) = PSCA*SCN1Y c MELVA3.VELCHE(J,NLCF) = PSCA*SCN1Z c ENDDO c DO J= 1,NBNN c WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA1= ', c & MELVA1.VELCHE(J,NLCF) c WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA2= ', c & MELVA2.VELCHE(J,NLCF) c WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA3= ', c & MELVA3.VELCHE(J,NLCF) c WRITE(6,*) 'MELEME=',MELEME.NUM(J,NLCF) c ENDDO C IF (1.EQ.0) THEN C WRITE(6,*) 'NGCF= ',NGCF C WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1) C WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2) C WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3) C VALD = MPOCHP.VPOCHA(NLCD,1) C WRITE(6,*) 'NLCF= ',NLCF,'VALD= ',VALD C WRITE(6,*) 'NLCF= ',NLCF,'VAL= ',VAL C WRITE(6,*) 'NLCF= ',NLCF,'TRG= ',TRG C WRITE(6,*) 'NLCF= ',NLCF,'TRG2= ',TRG2 C WRITE(6,*) 'NLCF= ',NLCF,'TRG3= ',TRG3 C WRITE(6,*) 'NLCF= ',NLCF,'TRG= ',TRGAUX C WRite(6,*) 'AG1=',AG1 C WRite(6,*) 'AG2=',AG2 C WRite(6,*) 'AD1=',AD1 C WRite(6,*) 'AD2=',AD2 C ENDIF NAUX2 = MAX(NAUX2,NCON) NMOY = NMOY + NCON c WRITE(6,*) 'INDGA= ',INDGA c WRITE(6,*) 'INDDR= ',INDDR c DO J= 1,NBNN c WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA1= ', c & MELVA1.VELCHE(J,NLCF) c ENDDO ENDDO NMOY = NMOY/(NFAC*1.D0) IF (1.EQ.0) THEN DO NLCF = 1, NFAC, 1 MPOGRA.VPOCHA(NLCF,1) = 0.D0 SCNX=MPONOR.VPOCHA(NLCF,1) SCNY=MPONOR.VPOCHA(NLCF,2) SCNZ=MPONOR.VPOCHA(NLCF,3) SCN1X = SCNX SCN1Y = SCNY SCN1Z = SCNZ NGCF=MELEFL.NUM(2,NLCF) DO IVOI=2,MELEME.NUM(/1) ICENT = MELEME.NUM(IVOI,NLCF) ICEN = MLECEN.LECT(ICENT) VAL = 0.0D0 IF (ICEN.NE.0) THEN c WRITE(6,*) 'INTERIEUR' VAL = MPOCHP.VPOCHA(ICEN,1) ELSE ICENL = MLENCL.LECT(ICENT) IF (ICENL.GT.0) THEN c WRITE(6,*) 'DIRICHLET' VAL = MPOVCL.VPOCHA(ICENL,1) ENDIF ENDIF COEF1 = MELVA1.VELCHE(IVOI,NLCF) MPOGRA.VPOCHA(NLCF,1)= MPOGRA.VPOCHA(NLCF,1) + & (COEF1 * VAL) c WRITE(6,*) 'NLCF= ',NLCF,'VAL= ',VAL c WRITE(6,*) 'IVOI= ',IVOI,'MELEME= ', MELEME.NUM(IVOI,NLCF), c & 'COEF1 = ',COEF1,'COEF2= ',COEF2,'COEF3= ',COEF3 ENDDO c DO J= 1,NBNN c PSCA = (MPOGRA.VPOCHA(NLCF,1)*SCN1X) + c & (MPOGRA.VPOCHA(NLCF,2)*SCN1Y) + c & (MPOGRA.VPOCHA(NLCF,3)*SCN1Z) c MPOGRA.VPOCHA(NLCF,1) = PSCA*SCN1X c MPOGRA.VPOCHA(NLCF,2) = PSCA*SCN1Y c MPOGRA.VPOCHA(NLCF,3) = PSCA*SCN1Z c ENDDO c NGCF=MELEFL.NUM(2,NLCF) c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF c WRITE(6,*) 'MPOGRA1= ', MPOGRA.VPOCHA(NLCF,1) c WRITE(6,*) 'MPOGRA2= ', MPOGRA.VPOCHA(NLCF,2) c WRITE(6,*) 'MPOGRA3= ', MPOGRA.VPOCHA(NLCF,3) ENDDO ENDIF IF (NBSO.EQ.2) THEN SEGDES IPT1 SEGDES IPT2 ELSEIF (NBSO.EQ.3) THEN SEGDES IPT1 SEGDES IPT2 SEGDES IPT3 ELSEIF (NBSO.EQ.4) THEN SEGDES IPT1 SEGDES IPT2 SEGDES IPT3 SEGDES IPT4 ENDIF IF (NBSOF.EQ.2) THEN SEGDES IPT5 SEGDES IPT6 ENDIF C ON REJUSTE LE CHAMELEM c WRITE(6,*) 'NAUX2= ',NAUX2 c WRITE(6,*) 'NMOY2= ',NMOY c IF (NAUX2.GT.NESSAI) THEN c WRITE(6,*) 'NAUX2 = ',NAUX2,'NESSAI = ',NESSAI c WRITE(6,*) 'NESSAI TROP PETIT' c STOP c ENDIF N1PTEL=NAUX2 N1EL=NBELEM N2PTEL=0 N2EL=0 SEGADJ MELVA1 NBNN = NAUX2 SEGADJ MELEME K7 = NFAC K8 = NAUX2 SEGINI ITAB c ON SUPPRIME LES ZEROS INTERIEURS IND = 2 DO NLCF = 1,NFAC IND = 2 NFIN = NAUX2 C ON CALCULE D4ABORD LE MAXIMUM DE LA LIGNE AUXMA = 0 DO J=2,NFIN ICEN = MLECEN.LECT(MELEME.NUM(J,NLCF)) IF (ICEN.NE.0) THEN AUXMA = MAX(ABS(MELVA1.VELCHE(J,NLCF)),AUXMA) ENDIF ENDDO ITAB.XMAX(NLCF) = AUXMA DO J=IND, NFIN ICEN = MLECEN.LECT(MELEME.NUM(J,NLCF)) IF (ABS(MELVA1.VELCHE(J,NLCF)).Le.(1e-14*AUXMA) & .AND.(ICEN.NE.0) ) THEN DO K=1,NFIN-J AUX = MELVA1.VELCHE(J+K,NLCF) ICEN = MLECEN.LECT(MELEME.NUM(J+K,NLCF)) C ON DECALE DE K CRAN IF (ABS(AUX).gt.(1e-14*AUXMA).OR.(ICEN.EQ.0)) THEN DO I=J,NFIN - K MELVA1.VELCHE(I,NLCF) = MELVA1.VELCHE(I+K,NLCF) MELEME.NUM(I,NLCF) = MELEME.NUM(I+K,NLCF) ENDDO C MISE A ZERO DES TERMES DU BOUT DO I=NFIN-K+1,NFIN MELVA1.VELCHE(I,NLCF) = 0.0 ENDDO GOTO 2000 ENDIF ENDDO 2000 CONTINUE ENDIF ENDDO ENDDO c IF (1.EQ.0) THEN c ON RECONSTUIT UN CHAMELEM EN SUPPRIMANT LES 0 c TABLEAU CALCULANT LE NOMBRE DE VOISIN NON NUL POUR CHAQUE FACE NMOY = 0 C INF = NAUX2 DO NLCF = 1,NFAC IREC = 3 DO J=NAUX2,1,-1 IF (ABS(MELVA1.VELCHE(J,NLCF)).gt. & (1e-14*ITAB.XMAX(NLCF))) THEN IREC = J GOTO 1111 ENDIF ENDDO 1111 CONTINUE ITAB.TABL(NLCF) = IREC NMOY = NMOY + IREC ENDDO NMOY = NMOY/(NFAC*1.D0) c WRITE(6,*) 'NEWMOY2= ',NMOY C TAB1(U) TABLEAU QUI CONTIENT LE NOMBRE DE FACE AYANT U VOISIN C NMAX = 0 DO ICOUR = 1,NAUX2 ITAB.TABL1(ICOUR) = 0 ENDDO DO NLCF = 1,NFAC ICOUR = ITAB.TABL(NLCF) ITAB.TABL1(ICOUR)=ITAB.TABL1(ICOUR) + 1 ENDDO C ON COMPTE LE NOMVRE DE SOUS DOMAINE NTSOUS = 0 DO ICOUR = 1,NAUX2 IF (ITAB.TABL1(ICOUR) .NE.0) NTSOUS = NTSOUS +1 ENDDO C IPOS INDICE DE LA PREMIERE FACE AYANT I VOISIN C ICOUR INDICE COURANT INITIALISE A IPOS ITAB.IPOS(1) = 1 DO I = 2,NAUX2 ITAB.IPOS(I) = ITAB.IPOS(I-1) + ITAB.TABL1(I-1) ITAB.ICOURANT(I) = ITAB.IPOS(I) ENDDO c TABL2 TABLEAU QUI RANGE DANS L4ODRES DES SOUS DOMAINES LES FACES NLCF DO I =1,NFAC IHELP = ITAB.TABL(I) IAUX = ITAB.ICOURANT(IHELP) ITAB.TABL2(IAUX) = I IAUX2 = ITAB.TABL(I) ITAB.ICOURANT(IAUX2) = ITAB.ICOURANT(IAUX2) + 1 ENDDO C**** Initialisation du MCHELM C N1=NTSOUS N2=1 N3=6 L1=8 SEGINI MCHEL2 MCHEL2.TITCHE='Gradient' MCHEL2.IFOCHE=IFOUR C ISOUS=0 NBSOUS=0 NBREF=0 DO I1 = 1, NAUX2, 1 NBELEM=ITAB.TABL1(I1) IF(NBELEM .GT. 0)THEN ISOUS=ISOUS+1 NBNN=I1 SEGINI IPT8 C ITYPEL=32 -> 'POLY' ITYPEL=32 MCHEL2.IMACHE(ISOUS)=IPT8 SEGINI MCHAM2 MCHEL2.ICHAML(ISOUS)=MCHAM2 MCHAM2.NOMCHE(1)='SCAL' MCHAM2.TYPCHE(1)='REAL*8 ' N1PTEL=NBNN N1EL=NBELEM N2PTEL=0 N2EL=0 SEGINI MELVA2 MCHAM2.IELVAL(1)=MELVA2 DO I2=1,NBELEM,1 DO I3=1,NBNN,1 IAUX = ITAB.IPOS(I1) IP= ITAB.TABL2(I2 + IAUX - 1) c WRITE(6,*) 'IP= ',IP,'I1= ',I1,'I2=',I2,'I3=',I3 IPT8.NUM(I3,I2)= MELEME.NUM(I3,IP) MELVA2.VELCHE(I3,I2)=MELVA1.VELCHE(I3,IP) ENDDO ENDDO C SEGDES MCHAM2 SEGDES IPT8 SEGDES MELVA2 ENDIF ENDDO SEGDES MCHEL2 SEGDES ITAB c ENDIF c SEGDES MCHAML c SEGDES MELEME c SEGDES MELVA1 c SEGDES MCHELM C 3009 SEGSUP MCHAML SEGSUP MELEME SEGSUP MELVA1 SEGSUP MCHELM C SUPRESSION DES SEGMENTS K3 = NSOMM DO I = 1,K3 MATR1 = IPO2.POINT(I) SEGACT MATR1 SEGSUP MATR1 ENDDO SEGSUP IPO3 SEGSUP IPO2 SEGSUP INDLI SEGSUP TAB SEGSUP IND2 c SEGSUP IND c SEGSUP IND22 c SEGSUP VAL1 c SEGSUP VAL2 SEGSUP SCMB C 3009 SEGSUP NBCOT RETURN END