c@a c@versb C----------------------------------------------------------------------- C CVERS Code_Saturne version 1.3 C ------------------------ C C This file is part of the Code_Saturne Kernel, element of the C Code_Saturne CFD tool. C C Copyright (C) 1998-2007 EDF S.A., France C C contact: saturne-support@edf.fr C C The Code_Saturne Kernel is free software; you can redistribute it C and/or modify it under the terms of the GNU General Public License C as published by the Free Software Foundation; either version 2 of C the License, or (at your option) any later version. C C The Code_Saturne Kernel is distributed in the hope that it will be C useful, but WITHOUT ANY WARRANTY; without even the implied warranty C of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the C GNU General Public License for more details. C C You should have received a copy of the GNU General Public License C along with the Code_Saturne Kernel; if not, write to the C Free Software Foundation, Inc., C 51 Franklin St, Fifth Floor, C Boston, MA 02110-1301 USA C C----------------------------------------------------------------------- c@verse SUBROUTINE RIJECH C ***************** C ------------------------------------------------------------- & ( IDBIA0 , IDBRA0 , & NDIM , NCELET , NCEL , NFAC , NFABOR , NFML , NPRFML , & NNOD , LNDFAC , LNDFBR , NCELBR , & NVAR , NSCAL , NPHAS , & NIDEVE , NRDEVE , NITUSE , NRTUSE , & IPHAS , IVAR , ISOU , IPP , & IFACEL , IFABOR , IFMFBR , IFMCEL , IPRFML , & IPNFAC , NODFAC , IPNFBR , NODFBR , & IDEVEL , ITUSER , IA , & XYZCEN , SURFAC , SURFBO , CDGFAC , CDGFBO , XYZNOD , VOLUME , & RTP , RTPA , PROPCE , PROPFA , PROPFB , & COEFA , COEFB , PRODUC , SMBR , & COEFAX , COEFBX , & PRODUK , W2 , W3 , W4 , EPSK , W6 , & RDEVEL , RTUSER , RA ) C ------------------------------------------------------------- C*********************************************************************** C FONCTION : C ---------- c@foncb CFONC CFONC TERMES D'ECHO DE PAROI CFONC POUR Rij CFONC VAR = R11 R22 R33 R12 R13 R23 CFONC ISOU = 1 2 3 4 5 6 CFONC c@fonce C----------------------------------------------------------------------- c@argub CARGU ARGUMENTS CARGU .______________.____._____.______________________________________. CARGU ! NOM !TYPE!MODE ! ROLE ! CARGU !______________!____!_____!______________________________________! CARGU ! IDBIA0 ! E ! -> ! NUMERO DE LA 1ERE CASE LIBRE DANS IA ! CARGU ! IDBRA0 ! E ! -> ! NUMERO DE LA 1ERE CASE LIBRE DANS RA ! CARGU ! NDIM ! E ! -> ! DIMENSION DE L'ESPACE ! CARGU ! NCELET ! E ! -> ! NOMBRE D'ELEMENTS HALO COMPRIS ! CARGU ! NCEL ! E ! -> ! NOMBRE D'ELEMENTS ACTIFS ! CARGU ! NFAC ! E ! -> ! NOMBRE DE FACES INTERNES ! CARGU ! NFABOR ! E ! -> ! NOMBRE DE FACES DE BORD ! CARGU ! NFML ! E ! -> ! NOMBRE DE FAMILLES D ENTITES ! CARGU ! NPRFML ! E ! -> ! NOMBRE DE PROPRIETESE DES FAMILLES ! CARGU ! NNOD ! E ! -> ! NOMBRE DE SOMMETS ! CARGU ! LNDFAC ! E ! -> ! LONGUEUR DU TABLEAU NODFAC (OPTIONNEL! CARGU ! LNDFBR ! E ! -> ! LONGUEUR DU TABLEAU NODFBR (OPTIONNEL! CARGU ! NCELBR ! E ! -> ! NOMBRE D'ELEMENTS AYANT AU MOINS UNE ! CARGU ! ! ! ! FACE DE BORD ! CARGU ! NVAR ! E ! -> ! NOMBRE TOTAL DE VARIABLES ! CARGU ! NSCAL ! E ! -> ! NOMBRE TOTAL DE SCALAIRES ! CARGU ! NPHAS ! E ! -> ! NOMBRE DE PHASES ! CARGU ! NIDEVE NRDEVE! E ! -> ! LONGUEUR DE IDEVEL RDEVEL ! CARGU ! NITUSE NRTUSE! E ! -> ! LONGUEUR DE ITUSER RTUSER ! CARGU ! IPHAS ! E ! -> ! NUMERO DE PHASE ! CARGU ! IVAR ! E ! -> ! NUMERO DE VARIABLE ! CARGU ! ISOU ! E ! -> ! NUMERO DE PASSAGE ! CARGU ! IPP ! E ! -> ! NUMERO DE VARIABLE POUR SORTIES POST ! CARGU ! IFACEL ! TE ! -> ! ELEMENTS VOISINS D'UNE FACE INTERNE ! CARGU ! (2, NFAC) ! ! ! ! CARGU ! IFABOR ! TE ! -> ! ELEMENT VOISIN D'UNE FACE DE BORD ! CARGU ! (NFABOR) ! ! ! ! CARGU ! IFMFBR ! TE ! -> ! NUMERO DE FAMILLE D'UNE FACE DE BORD ! CARGU ! (NFABOR) ! ! ! ! CARGU ! IFMCEL ! TE ! -> ! NUMERO DE FAMILLE D'UNE CELLULE ! CARGU ! (NCELET) ! ! ! ! CARGU ! IPRFML ! TE ! -> ! PROPRIETES D'UNE FAMILLE ! CARGU ! NFML ,NPRFML! ! ! ! CARGU ! IPNFAC ! TE ! -> ! POSITION DU PREMIER NOEUD DE CHAQUE ! CARGU ! (LNDFAC) ! ! ! FACE INTERNE DANS NODFAC (OPTIONNEL)! CARGU ! NODFAC ! TE ! -> ! CONNECTIVITE FACES INTERNES/NOEUDS ! CARGU ! (NFAC+1) ! ! ! (OPTIONNEL) ! CARGU ! IPNFBR ! TE ! -> ! POSITION DU PREMIER NOEUD DE CHAQUE ! CARGU ! (LNDFBR) ! ! ! FACE DE BORD DANS NODFBR (OPTIONNEL)! CARGU ! NODFBR ! TE ! -> ! CONNECTIVITE FACES DE BORD/NOEUDS ! CARGU ! (NFABOR+1) ! ! ! (OPTIONNEL) ! CARGU ! IDEVEL(NIDEVE! TE ! <-> ! TAB ENTIER COMPLEMENTAIRE DEVELOPEMT ! CARGU ! ITUSER(NITUSE! TE ! <-> ! TAB ENTIER COMPLEMENTAIRE UTILISATEUR! CARGU ! IA(*) ! TR ! - ! MACRO TABLEAU ENTIER ! CARGU ! XYZCEN ! TR ! -> ! POINT ASSOCIES AUX VOLUMES DE CONTROL! CARGU ! (NDIM,NCELET ! ! ! ! CARGU ! SURFAC ! TR ! -> ! VECTEUR SURFACE DES FACES INTERNES ! CARGU ! (NDIM,NFAC) ! ! ! ! CARGU ! SURFBO ! TR ! -> ! VECTEUR SURFACE DES FACES DE BORD ! CARGU ! (NDIM,NFABOR)! ! ! ! CARGU ! CDGFAC ! TR ! -> ! CENTRE DE GRAVITE DES FACES INTERNES ! CARGU ! (NDIM,NFAC) ! ! ! ! CARGU ! CDGFBO ! TR ! -> ! CENTRE DE GRAVITE DES FACES DE BORD ! CARGU ! (NDIM,NFABOR)! ! ! ! CARGU ! XYZNOD ! TR ! -> ! COORDONNES DES NOEUDS (OPTIONNEL) ! CARGU ! (NDIM,NNOD) ! ! ! ! CARGU ! VOLUME ! TR ! -> ! VOLUME D'UN DES NCELET ELEMENTS ! CARGU ! (NCELET ! ! ! ! CARGU ! RTP, RTPA ! TR ! -> ! VARIABLES DE CALCUL AU CENTRE DES ! CARGU ! (NCELET,*) ! ! ! CELLULES (INSTANT COURANT OU PREC)! CARGU ! PROPCE ! TR ! -> ! PROPRIETES PHYSIQUES AU CENTRE DES ! CARGU ! (NCELET,*) ! ! ! CELLULES ! CARGU ! PROPFA ! TR ! -> ! PROPRIETES PHYSIQUES AU CENTRE DES ! CARGU ! (NFAC,*) ! ! ! FACES INTERNES ! CARGU ! PROPFB ! TR ! -> ! PROPRIETES PHYSIQUES AU CENTRE DES ! CARGU ! (NFABOR,*) ! ! ! FACES DE BORD ! CARGU ! COEFA, COEFB ! TR ! -> ! CONDITIONS AUX LIMITES AUX ! CARGU ! (NFABOR,*) ! ! ! FACES DE BORD ! CARGU ! PRODUC ! TR ! -> ! PRODUCTION ! CARGU ! (6,NCELET) ! ! ! ! CARGU ! SMBR(NCELET ! TR ! <-> ! TABLEAU DE TRAVAIL POUR SEC MEM ! CARGU ! COEFAX,COEFBX! TR ! - ! TABLEAU DE TRAVAIL POUR LES COND. ! CARGU ! (NFABOR) ! ! ! AUX LIMITES DE LA DIST. PAROI ! CARGU ! PRODUK(NCELET! TR ! - ! TABLEAU DE TRAVAIL PRODUCTION ! CARGU ! EPSK (NCELET! TR ! - ! TABLEAU DE TRAVAIL EPSILON/K ! CARGU ! W2...6(NCELET! TR ! - ! TABLEAU DE TRAVAIL PRODUCTION ! CARGU ! RDEVEL(NRDEVE! TR ! <-> ! TAB REEL COMPLEMENTAIRE DEVELOPEMT ! CARGU ! RTUSER(NRTUSE! TR ! <-> ! TAB REEL COMPLEMENTAIRE UTILISATEUR ! CARGU ! RA(*) ! TR ! - ! MACRO TABLEAU REEL ! CARGU !______________!____!_____!______________________________________! c@argue C c@commb CCOMM COMMONS CCOMM .______________.____._____.______________________________________. CCOMM ! NOM !TYPE!MODE ! ROLE ! CCOMM !______________!____!_____!______________________________________! CCOMM !______________!____!_____!______________________________________! c@comme C C TYPE : E (ENTIER), R (REEL), A (ALPHANUMERIQUE), T (TABLEAU) C L (LOGIQUE) .. ET TYPES COMPOSES (EX : TR TABLEAU REEL) C MODE : -> DONNEE, <- RESULTAT, <-> DONNEE MODIFIEE, C - TABLEAU DE TRAVAIL C----------------------------------------------------------------------- C*********************************************************************** C IMPLICIT NONE C C*********************************************************************** C DONNEES EN COMMON C*********************************************************************** C INCLUDE "paramx.h" INCLUDE "numvar.h" INCLUDE "entsor.h" INCLUDE "optcal.h" INCLUDE "cstphy.h" INCLUDE "cstnum.h" INCLUDE "pointe.h" INCLUDE "period.h" INCLUDE "parall.h" C C*********************************************************************** C C ARGUMENTS C INTEGER IDBIA0 , IDBRA0 INTEGER NDIM , NCELET , NCEL , NFAC , NFABOR INTEGER NFML , NPRFML INTEGER NNOD , LNDFAC , LNDFBR , NCELBR INTEGER NVAR , NSCAL , NPHAS INTEGER NIDEVE , NRDEVE , NITUSE , NRTUSE INTEGER IPHAS , IVAR , ISOU , IPP C INTEGER IFACEL(2,NFAC) , IFABOR(NFABOR) INTEGER IFMFBR(NFABOR) , IFMCEL(NCELET) INTEGER IPRFML(NFML,NPRFML) INTEGER IPNFAC(NFAC+1), NODFAC(LNDFAC) INTEGER IPNFBR(NFABOR+1), NODFBR(LNDFBR) INTEGER IDEVEL(NIDEVE), ITUSER(NITUSE), IA(*) C DOUBLE PRECISION XYZCEN(NDIM,NCELET) DOUBLE PRECISION SURFAC(NDIM,NFAC), SURFBO(NDIM,NFABOR) DOUBLE PRECISION CDGFAC(NDIM,NFAC), CDGFBO(NDIM,NFABOR) DOUBLE PRECISION XYZNOD(NDIM,NNOD), VOLUME(NCELET) DOUBLE PRECISION RTP(NCELET,*), RTPA(NCELET,*) DOUBLE PRECISION PROPCE(NCELET,*) DOUBLE PRECISION PROPFA(NFAC,*), PROPFB(NFABOR,*) DOUBLE PRECISION COEFA(NFABOR,*), COEFB(NFABOR,*) DOUBLE PRECISION PRODUC(6,NCELET) DOUBLE PRECISION SMBR(NCELET) DOUBLE PRECISION COEFAX(NFABOR), COEFBX(NFABOR) DOUBLE PRECISION PRODUK(NCELET),W2(NCELET),W3(NCELET) DOUBLE PRECISION W4(NCELET),EPSK(NCELET),W6(NCELET) DOUBLE PRECISION RDEVEL(NRDEVE), RTUSER(NRTUSE), RA(*) C C VARIABLES LOCALES C INTEGER IDEBIA, IDEBRA INTEGER IFACPT, IEL , II , JJ , KK , MM INTEGER IRKM , IRKI , IRKJ , ISKM , ISKI , ISKJ INTEGER IR11IP, IR22IP, IR33IP, IR12IP, IR13IP, IR23IP INTEGER IEIPH , IPCROM, IPCROO INTEGER IFAC , IDIMTE, ITENSO INTEGER INC , ICCOCG, IPHYDP, IVAR0 , IITYPH C DOUBLE PRECISION CMU075, DISTXN, D2S3 , TRRIJ , XK DOUBLE PRECISION UNSSUR, VNK , VNM , VNI , VNJ DOUBLE PRECISION DELTKI, DELTKJ, DELTKM, DELTIJ DOUBLE PRECISION AA , BB , SURFBN, XNORME C C C*********************************************************************** C C======================================================================= C 1. INITIALISATION C======================================================================= C IDEBIA = IDBIA0 IDEBRA = IDBRA0 C IR11IP = IR11(IPHAS) IR22IP = IR22(IPHAS) IR33IP = IR33(IPHAS) IR12IP = IR12(IPHAS) IR13IP = IR13(IPHAS) IR23IP = IR23(IPHAS) IEIPH = IEP (IPHAS) C IPCROM = IPPROC(IROM (IPHAS)) IPCROO = IPCROM IF(ISTO2T(IPHAS).GT.0.AND.IROEXT(IPHAS).GT.0) THEN IPCROO = IPPROC(IROMA(IPHAS)) ENDIF C DELTIJ = 1.0D0 IF(ISOU.GT.3) THEN DELTIJ = 0.0D0 ENDIF C CMU075 = CMU**0.75D0 D2S3 = 2.D0/3.D0 C C======================================================================= C 2. CALCUL AUX CELLULES DES NORMALES EN PAROI CORRESPONDANTES C======================================================================= C C On stocke les composantes dans les tableaux de travail W2, W3 et W4 C IF(ABS(ICDPAR).EQ.2) THEN C C On connait la face de paroi correspondante C DO IEL = 1, NCEL IFACPT = IA(IIFAPA(IPHAS)-1+IEL) SURFBN = RA(ISRFBN-1+IFACPT) UNSSUR = 1.D0/SURFBN W2(IEL)= SURFBO(1,IFACPT)*UNSSUR W3(IEL)= SURFBO(2,IFACPT)*UNSSUR W4(IEL)= SURFBO(3,IFACPT)*UNSSUR ENDDO C ELSEIF(ABS(ICDPAR).EQ.1) THEN C C La normale est definie comme - gradient de la distance C a la paroi C C La distance a la paroi vaut 0 en paroi C par definition et obeit a un flux nul ailleurs C IITYPH = IITYPF+(IPHAS-1)*NFABOR DO IFAC = 1, NFABOR IF(IA(IITYPH+IFAC-1).EQ.IPAROI) THEN COEFAX(IFAC) = 0.0D0 COEFBX(IFAC) = 0.0D0 ELSE COEFAX(IFAC) = 0.0D0 COEFBX(IFAC) = 1.0D0 ENDIF ENDDO C C Calcul du gradient C IF(IRANGP.GE.0) CALL PARCOM (RA(IDIPAR)) C =========== IF(IPERIO.EQ.1) THEN IDIMTE = 0 ITENSO = 0 CALL PERCOM C =========== & ( IDIMTE , ITENSO , & RA(IDIPAR), RA(IDIPAR), RA(IDIPAR), & RA(IDIPAR), RA(IDIPAR), RA(IDIPAR), & RA(IDIPAR), RA(IDIPAR), RA(IDIPAR)) ENDIF C INC = 1 ICCOCG = 1 IPHYDP = 0 IVAR0 = 0 C CALL GRDCEL C =========== & ( IDEBIA , IDEBRA , & NDIM , NCELET , NCEL , NFAC , NFABOR , NFML , NPRFML , & NNOD , LNDFAC , LNDFBR , NCELBR , NPHAS , & NIDEVE , NRDEVE , NITUSE , NRTUSE , & IVAR0 , IMRGRA , INC , ICCOCG , NSWRGY , IMLIGY , IPHYDP , & IWARNY , NFECRA , & EPSRGY , CLIMGY , EXTRAY , & IFACEL , IFABOR , IFMFBR , IFMCEL , IPRFML , & IPNFAC , NODFAC , IPNFBR , NODFBR , & IDEVEL , ITUSER , IA , & XYZCEN , SURFAC , SURFBO , CDGFAC , CDGFBO , XYZNOD , VOLUME , & RA , RA , RA , & RA(IDIPAR) , COEFAX , COEFBX , & W2 , W3 , W4 , C ------ ------ ------ & PRODUK , EPSK , W6 , & RDEVEL , RTUSER , RA ) C C C Normalisation (attention, le gradient peut etre nul, parfois) C DO IEL = 1 ,NCEL XNORME = MAX(SQRT(W2(IEL)**2+W3(IEL)**2+W4(IEL)**2),EPZERO) W2(IEL) = -W2(IEL)/XNORME W3(IEL) = -W3(IEL)/XNORME W4(IEL) = -W4(IEL)/XNORME ENDDO ENDIF C C======================================================================= C 2. CALCUL DE VARIABLES DE TRAVAIL C======================================================================= C C C C ---> Production et k C DO IEL = 1 , NCEL PRODUK(IEL) = 0.5D0 * (PRODUC(1,IEL) + PRODUC(2,IEL) + & PRODUC(3,IEL)) XK = 0.5D0 * (RTPA(IEL,IR11IP) + RTPA(IEL,IR22IP) + & RTPA(IEL,IR33IP)) EPSK(IEL) = RTPA(IEL,IEIPH)/XK ENDDO C C C C ---> Indices des tensions C IF ((ISOU.EQ.1).OR.(ISOU.EQ.4).OR.(ISOU.EQ.5)) THEN II = 1 ELSEIF ((ISOU.EQ.2).OR.(ISOU.EQ.6)) THEN II = 2 ELSEIF (ISOU.EQ.3) THEN II = 3 ENDIF C IF ((ISOU.EQ.3).OR.(ISOU.EQ.5).OR.(ISOU.EQ.6)) THEN JJ = 3 ELSEIF ((ISOU.EQ.2).OR.(ISOU.EQ.4)) THEN JJ = 2 ELSEIF (ISOU.EQ.1) THEN JJ = 1 ENDIF C C ---> Boucle pour construction des termes sources C DO IEL = 1, NCEL W6(IEL) = 0.D0 ENDDO C DO KK = 1, 3 C C ---> Sommes sur m C DO MM = 1, 3 C C --> Delta km C IF(KK.EQ.MM) THEN DELTKM = 1.0D0 ELSE DELTKM = 0.0D0 ENDIF C C --> R km C IF ((KK*MM).EQ.1) THEN IRKM = IR11IP ISKM = 1 ELSEIF ((KK*MM).EQ.4) THEN IRKM = IR22IP ISKM = 2 ELSEIF ((KK*MM).EQ.9) THEN IRKM = IR33IP ISKM = 3 ELSEIF ((KK*MM).EQ.2) THEN IRKM = IR12IP ISKM = 4 ELSEIF ((KK*MM).EQ.3) THEN IRKM = IR13IP ISKM = 5 ELSEIF ((KK*MM).EQ.6) THEN IRKM = IR23IP ISKM = 6 ENDIF C C --> Termes en R km et Phi km C DO IEL = 1, NCEL C IF (KK.EQ.1) THEN VNK = W2(IEL) ELSEIF(KK.EQ.2) THEN VNK = W3(IEL) ELSEIF(KK.EQ.3) THEN VNK = W4(IEL) ENDIF C IF (MM.EQ.1) THEN VNM = W2(IEL) ELSEIF(MM.EQ.2) THEN VNM = W3(IEL) ELSEIF(MM.EQ.3) THEN VNM = W4(IEL) ENDIF C W6(IEL) = W6(IEL) + VNK*VNM*DELTIJ*( & CRIJP1*RTPA(IEL,IRKM)*EPSK(IEL) & -CRIJP2 & *CRIJ2*(PRODUC(ISKM,IEL)-D2S3*PRODUK(IEL)*DELTKM) ) ENDDO C ENDDO C C ---> Fin des sommes sur m C C C --> R ki C IF ((KK*II).EQ.1) THEN IRKI = IR11IP ISKI = 1 ELSEIF ((KK*II).EQ.4) THEN IRKI = IR22IP ISKI = 2 ELSEIF ((KK*II).EQ.9) THEN IRKI = IR33IP ISKI = 3 ELSEIF ((KK*II).EQ.2) THEN IRKI = IR12IP ISKI = 4 ELSEIF ((KK*II).EQ.3) THEN IRKI = IR13IP ISKI = 5 ELSEIF ((KK*II).EQ.6) THEN IRKI = IR23IP ISKI = 6 ENDIF C C --> R kj C IF ((KK*JJ).EQ.1) THEN IRKJ = IR11IP ISKJ = 1 ELSEIF ((KK*JJ).EQ.4) THEN IRKJ = IR22IP ISKJ = 2 ELSEIF ((KK*JJ).EQ.9) THEN IRKJ = IR33IP ISKJ = 3 ELSEIF ((KK*JJ).EQ.2) THEN IRKJ = IR12IP ISKJ = 4 ELSEIF ((KK*JJ).EQ.3) THEN IRKJ = IR13IP ISKJ = 5 ELSEIF ((KK*JJ).EQ.6) THEN IRKJ = IR23IP ISKJ = 6 ENDIF C C --> Delta ki C IF (KK.EQ.II) THEN DELTKI = 1.D0 ELSE DELTKI = 0.D0 ENDIF C C --> Delta kj C IF (KK.EQ.JJ) THEN DELTKJ = 1.D0 ELSE DELTKJ = 0.D0 ENDIF C DO IEL = 1, NCEL C IF (KK.EQ.1) THEN VNK = W2(IEL) ELSEIF(KK.EQ.2) THEN VNK = W3(IEL) ELSEIF(KK.EQ.3) THEN VNK = W4(IEL) ENDIF C IF (II.EQ.1) THEN VNI = W2(IEL) ELSEIF(II.EQ.2) THEN VNI = W3(IEL) ELSEIF(II.EQ.3) THEN VNI = W4(IEL) ENDIF C IF (JJ.EQ.1) THEN VNJ = W2(IEL) ELSEIF(JJ.EQ.2) THEN VNJ = W3(IEL) ELSEIF(JJ.EQ.3) THEN VNJ = W4(IEL) ENDIF C W6(IEL) = W6(IEL) + 1.5D0*VNK*( & -CRIJP1*(RTPA(IEL,IRKI)*VNJ+RTPA(IEL,IRKJ)*VNI)*EPSK(IEL) & +CRIJP2 & *CRIJ2*((PRODUC(ISKI,IEL)-D2S3*PRODUK(IEL)*DELTKI)*VNJ & +(PRODUC(ISKJ,IEL)-D2S3*PRODUK(IEL)*DELTKJ)*VNI) ) C ENDDO C ENDDO C C C ---> Distance a la paroi et fonction d'amortissement : W3 C Pour chaque mode de calcul : meme code, test C en dehors de la boucle C IF(ABS(ICDPAR).EQ.2) THEN DO IEL = 1 , NCEL IFACPT = IA(IIFAPA(IPHAS)-1+IEL) DISTXN = & (CDGFBO(1,IFACPT)-XYZCEN(1,IEL))**2 & +(CDGFBO(2,IFACPT)-XYZCEN(2,IEL))**2 & +(CDGFBO(3,IFACPT)-XYZCEN(3,IEL))**2 DISTXN = SQRT(DISTXN) TRRIJ = 0.5D0 * (RTPA(IEL,IR11IP) + RTPA(IEL,IR22IP) + & RTPA(IEL,IR33IP)) AA = 1.D0 BB = CMU075*TRRIJ**1.5D0/(XKAPPA*RTPA(IEL,IEIPH)*DISTXN) W3(IEL) = MIN(AA, BB) ENDDO ELSE DO IEL = 1 , NCEL DISTXN = MAX(RA(IDIPAR+IEL-1),EPZERO) TRRIJ = 0.5D0 * (RTPA(IEL,IR11IP) + RTPA(IEL,IR22IP) + & RTPA(IEL,IR33IP)) AA = 1.D0 BB = CMU075*TRRIJ**1.5D0/(XKAPPA*RTPA(IEL,IEIPH)*DISTXN) W3(IEL) = MIN(AA, BB) ENDDO ENDIF C C C ---> Increment du terme source C DO IEL = 1, NCEL SMBR(IEL) = SMBR(IEL) & + PROPCE(IEL,IPCROO)*VOLUME(IEL)*W6(IEL)*W3(IEL) ENDDO C C RETURN C END c@z