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 LIMGRD C ***************** C ------------------------------------------------------------- & ( NDIM , NCELET , NCEL , NFAC , & IMRGRA , IMLIGP , IDIMTE , ITENSO , IWARNP , NFECRA , CLIMGP , & IFACEL , IPCVSE , IELVSE , & XYZCEN , PVAR , & DPDX , DPDY , DPDZ , & FACLIM , DENOM , DENUM ) C ------------------------------------------------------------- C*********************************************************************** C FONCTION : C ---------- c@foncb CFONC LIMITATION DU GRADIENT CELLULE CFONC CFONC MAX | IJ.GRADP | <= k MAX | P -P | avec k>=1 CFONC j e V I J I CFONC i c@fonce C----------------------------------------------------------------------- C ARGUMENTS c@argub CARGU .______________.____._____.______________________________________. CARGU ! NOM !TYPE!MODE ! ROLE ! CARGU !______________!____!_____!______________________________________! 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 ! IMRGRA ! E ! -> ! METHODE DE RECONSTRUCTION DU GRADIENT! CARGU ! ! ! ! 0 RECONSTRUCTION 97 ! CARGU ! ! ! ! 1 MOINDRE CARRE ! CARGU ! ! ! ! 2 MOINDRE CARRE SUPPORT ETENDU ! CARGU ! ! ! ! 3 MOINDRE CARRE SUPPORT ETENDU REDUI! CARGU ! IMLIGP ! E ! -> ! METHODE DE LIMITATION DU GRADIENT ! CARGU ! ! ! ! < 0 PAS DE LIMITATION ! CARGU ! ! ! ! = 0 A PARTIR DES GRADIENTS VOISINS ! CARGU ! ! ! ! = 1 A PARTIR DU GRADIENT MOYEN ! CARGU ! IDIMTE ! E ! -> ! DIMENSION DE LA VARIBLE (MAXIMUM 3) ! CARGU ! ! ! ! 0 : SCALAIRE (VAR11) OU ASSIMILE ! CARGU ! ! ! ! SCALAIRE ! CARGU ! ! ! ! 1 : VECTEUR (VAR11,VAR22,VAR33) ! CARGU ! ! ! ! 2 : TENSEUR D'ORDRE 2 (VARIJ) ! CARGU ! ITENSO ! E ! -> ! POUR L'EXPLICITATION DE LA ROTATION ! CARGU ! ! ! ! 0 : SCALAIRE (VAR11) ! CARGU ! ! ! ! 1 : COMPOSANTE DE VECTEUR OU DE ! CARGU ! ! ! ! TENSEUR (VAR11) IMPLCITE POUR LA ! CARGU ! ! ! ! TRANSLATION ! CARGU ! ! ! !11 : REPREND LE TRAITEMENT ITENSO=1 ET! CARGU ! ! ! ! COMPOSANTE DE VECTEUR OU DE ! CARGU ! ! ! ! TENSEUR (VAR11) ANNULEE POUR LA ! CARGU ! ! ! ! ROTATION ! CARGU ! ! ! ! 2 : VECTEUR (VAR11 ET VAR22 ET VAR33)! CARGU ! ! ! ! IMPLICITE POUR LA TRANSLATION ! CARGU ! IWARNP ! E ! -> ! NIVEAU D'IMPRESSION ! CARGU ! NFECRA ! E ! -> ! UNITE DU FICHIER SORTIE STD ! CARGU ! CLIMGP ! R ! -> ! COEF GRADIENT*DISTANCE/ECART ! CARGU ! IFACEL ! TE ! -> ! ELEMENTS VOISINS D'UNE FACE INTERNE ! CARGU ! (2, NFAC) ! ! ! ! CARGU ! IPCVSE ! TE ! -> ! POSITION DANS IELVSE DES VOISINS ! CARGU ! (NCEL ) ! ! ! ETENDUS DES CELLULES ! CARGU ! IELVSE (*) ! TE ! -> ! NUMERO DES VOISINS ETENDUS DES CELLUL! CARGU ! XYZCEN ! TR ! -> ! POINT ASSOCIES AUX VOLUMES DE CONTROL! CARGU ! (NDIM,NCELET ! ! ! ! CARGU ! PVAR (NCELET! TR ! -> ! VARIABLE ! CARGU ! DPD. (NCELET! TR ! <-> ! GRADIENT DE PVAR ! CARGU ! FACLIM(NCELET! TR ! - ! FACTEUR DE LIMITATION (IMLIGP = 1) ! CARGU ! DENOM (NCELET! TR ! - ! VAR MAX BASEE SUR LA VARIABLE ! CARGU ! DENUM (NCELET! TR ! - ! VAR MAX BASEE SUR LE GRADIENT ! 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 "vector.h" INCLUDE "parall.h" INCLUDE "period.h" C C*********************************************************************** C C ARGUMENTS C INTEGER NDIM INTEGER NCELET , NCEL , NFAC INTEGER IMRGRA , IMLIGP INTEGER IDIMTE , ITENSO INTEGER IWARNP , NFECRA DOUBLE PRECISION CLIMGP C INTEGER IFACEL(2,NFAC) INTEGER IPCVSE(*), IELVSE(*) DOUBLE PRECISION XYZCEN(NDIM,*) DOUBLE PRECISION PVAR(*) DOUBLE PRECISION DPDX (*),DPDY (*),DPDZ (*) DOUBLE PRECISION FACLIM(NCELET),DENOM(*),DENUM(*) C C VARIABLES LOCALES C INTEGER IFAC, IEL, II, JJ, ICLIP INTEGER IDIMDE, ITENDE, IPOS DOUBLE PRECISION RI, RJ, RMAX, RMIN, RIJ, R1, R2, R3, ADP, VECFAC DOUBLE PRECISION DPDXF, DPDYF, DPDZF C C*********************************************************************** C C======================================================================= C 0. PAS DE LIMITATION C======================================================================= C IF (IMLIGP.LT.0) THEN RETURN ENDIF C C======================================================================= C 1. INITIALISATION C======================================================================= C DO IEL = 1, NCELET DENOM(IEL) = 0.D0 DENUM(IEL) = 0.D0 ENDDO C C C======================================================================= C 2. CALCULS PRELIMINAIRES C C DENUM CONTIENT LA VARIATION MAX BASEE SUR LE GRADIENT C DENOM CONTIENT LA VARIATION MAX BASEE SUR LA VARIABLE C C======================================================================= C IF( IMLIGP .EQ. 0 ) THEN C C IF (IVECTI.EQ.1) THEN C !OCL NOVREC,VRL(16) DO IFAC = 1, NFAC II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) R1 = XYZCEN(1,II)-XYZCEN(1,JJ) R2 = XYZCEN(2,II)-XYZCEN(2,JJ) R3 = XYZCEN(3,II)-XYZCEN(3,JJ) RI = ABS( R1*DPDX(II)+R2*DPDY(II)+R3*DPDZ(II) ) RJ = ABS( R1*DPDX(JJ)+R2*DPDY(JJ)+R3*DPDZ(JJ) ) DENUM(II) = MAX( RI, DENUM(II) ) DENUM(JJ) = MAX( RJ, DENUM(JJ) ) ENDDO C ELSE C C VECTORISATION NON FORCEE DO IFAC = 1, NFAC II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) R1 = XYZCEN(1,II)-XYZCEN(1,JJ) R2 = XYZCEN(2,II)-XYZCEN(2,JJ) R3 = XYZCEN(3,II)-XYZCEN(3,JJ) RI = ABS( R1*DPDX(II)+R2*DPDY(II)+R3*DPDZ(II) ) RJ = ABS( R1*DPDX(JJ)+R2*DPDY(JJ)+R3*DPDZ(JJ) ) DENUM(II) = MAX( RI, DENUM(II) ) DENUM(JJ) = MAX( RJ, DENUM(JJ) ) ENDDO C ENDIF C C COMPLEMENT POUR LE VOISINAGE ETENDU C PAS DE VECTORISATION PARTICULIERE A IMPOSER A PRIORI C IF (IMRGRA.EQ.2.OR.IMRGRA.EQ.3) THEN C DO II = 1, NCEL DO IPOS = IPCVSE(II), IPCVSE(II+1)-1 JJ = IELVSE(IPOS) R1 = XYZCEN(1,II)-XYZCEN(1,JJ) R2 = XYZCEN(2,II)-XYZCEN(2,JJ) R3 = XYZCEN(3,II)-XYZCEN(3,JJ) RI = ABS( R1*DPDX(II)+R2*DPDY(II)+R3*DPDZ(II) ) DENUM(II) = MAX( RI, DENUM(II) ) ENDDO ENDDO C ENDIF C C IF (IVECTI.EQ.1) THEN C !OCL NOVREC,VRL(16) DO IFAC = 1, NFAC II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) ADP = ABS(PVAR(II)-PVAR(JJ)) DENOM(II) = MAX ( DENOM(II), ADP ) DENOM(JJ) = MAX ( DENOM(JJ), ADP ) ENDDO C ELSE C C VECTORISATION NON FORCEE DO IFAC = 1, NFAC II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) ADP = ABS(PVAR(II)-PVAR(JJ)) DENOM(II) = MAX ( DENOM(II), ADP ) DENOM(JJ) = MAX ( DENOM(JJ), ADP ) ENDDO C ENDIF C C COMPLEMENT POUR LE VOISINAGE ETENDU C PAS DE VECTORISATION PARTICULIERE A IMPOSER A PRIORI C IF (IMRGRA.EQ.2.OR.IMRGRA.EQ.3) THEN C DO II = 1, NCEL DO IPOS = IPCVSE(II), IPCVSE(II+1)-1 JJ = IELVSE(IPOS) ADP = ABS(PVAR(II)-PVAR(JJ)) DENOM(II) = MAX ( DENOM(II), ADP ) ENDDO ENDDO C ENDIF C C ELSEIF(IMLIGP.EQ.1) THEN C C IF (IVECTI.EQ.1) THEN C !OCL NOVREC,VRL(16) DO IFAC =1, NFAC II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) DPDXF = 0.5D0*( DPDX(II) +DPDX(JJ) ) DPDYF = 0.5D0*( DPDY(II) +DPDY(JJ) ) DPDZF = 0.5D0*( DPDZ(II) +DPDZ(JJ) ) RIJ = ABS( (XYZCEN(1,II)-XYZCEN(1,JJ))*DPDXF & +(XYZCEN(2,II)-XYZCEN(2,JJ))*DPDYF & +(XYZCEN(3,II)-XYZCEN(3,JJ))*DPDZF ) DENUM(II) = MAX( RIJ, DENUM(II) ) DENUM(JJ) = MAX( RIJ, DENUM(JJ) ) ENDDO C ELSE C C VECTORISATION NON FORCEE DO IFAC =1, NFAC II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) DPDXF = 0.5D0*( DPDX(II) +DPDX(JJ) ) DPDYF = 0.5D0*( DPDY(II) +DPDY(JJ) ) DPDZF = 0.5D0*( DPDZ(II) +DPDZ(JJ) ) RIJ = ABS( (XYZCEN(1,II)-XYZCEN(1,JJ))*DPDXF & +(XYZCEN(2,II)-XYZCEN(2,JJ))*DPDYF & +(XYZCEN(3,II)-XYZCEN(3,JJ))*DPDZF ) DENUM(II) = MAX( RIJ, DENUM(II) ) DENUM(JJ) = MAX( RIJ, DENUM(JJ) ) ENDDO C ENDIF C C COMPLEMENT POUR LE VOISINAGE ETENDU C PAS DE VECTORISATION PARTICULIERE A IMPOSER A PRIORI C IF (IMRGRA.EQ.2.OR.IMRGRA.EQ.3) THEN C DO II = 1, NCEL DO IPOS = IPCVSE(II), IPCVSE(II+1)-1 JJ = IELVSE(IPOS) DPDXF = 0.5D0*( DPDX(II) +DPDX(JJ) ) DPDYF = 0.5D0*( DPDY(II) +DPDY(JJ) ) DPDZF = 0.5D0*( DPDZ(II) +DPDZ(JJ) ) RIJ = ABS( (XYZCEN(1,II)-XYZCEN(1,JJ))*DPDXF & +(XYZCEN(2,II)-XYZCEN(2,JJ))*DPDYF & +(XYZCEN(3,II)-XYZCEN(3,JJ))*DPDZF ) DENUM(II) = MAX( RIJ, DENUM(II) ) ENDDO ENDDO C ENDIF C C IF (IVECTI.EQ.1) THEN C !OCL NOVREC,VRL(16) DO IFAC =1, NFAC II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) ADP = ABS(PVAR(II)-PVAR(JJ)) DENOM(II) = MAX ( DENOM(II), ADP ) DENOM(JJ) = MAX ( DENOM(JJ), ADP ) ENDDO C ELSE C C VECTORISATION NON FORCEE DO IFAC =1, NFAC II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) ADP = ABS(PVAR(II)-PVAR(JJ)) DENOM(II) = MAX ( DENOM(II), ADP ) DENOM(JJ) = MAX ( DENOM(JJ), ADP ) ENDDO C ENDIF C C COMPLEMENT POUR LE VOISINAGE ETENDU C PAS DE VECTORISATION PARTICULIERE A IMPOSER A PRIORI C IF (IMRGRA.EQ.2.OR.IMRGRA.EQ.3) THEN C DO II = 1, NCEL DO IPOS = IPCVSE(II), IPCVSE(II+1)-1 JJ = IELVSE(IPOS) ADP = ABS(PVAR(II)-PVAR(JJ)) DENOM(II) = MAX ( DENOM(II), ADP ) ENDDO ENDDO C ENDIF C C ENDIF C C C======================================================================= C 3. LIMITATION DU GRADIENT C ON LIMITE LE RAPPORT ENTRE DENOM ET DENUM C A UN FACTEUR MULTIPLICATIF C======================================================================= C ICLIP = 0 RMAX = 0.D0 RMIN = 1.D0 C IF( IMLIGP.EQ.0 ) THEN C DO IEL = 1, NCEL IF( DENUM(IEL) .GT. CLIMGP*DENOM(IEL) ) THEN RI = CLIMGP*DENOM(IEL)/DENUM(IEL) DPDX(IEL) = RI*DPDX(IEL) DPDY(IEL) = RI*DPDY(IEL) DPDZ(IEL) = RI*DPDZ(IEL) RMAX = MAX( RMAX, RI ) RMIN = MIN( RMIN, RI ) ICLIP = ICLIP +1 ENDIF ENDDO IF (IRANGP.GE.0) THEN CALL PARMAX (RMAX) C =========== CALL PARMIN (RMIN) C =========== ENDIF C ELSEIF(IMLIGP.EQ.1) THEN C DO IEL = 1, NCELET FACLIM(IEL) = 1.D+12 ENDDO C C ECHANGES PARALLELES ET PERIODIQUES POUR DENOM ET DENUM C IF(IRANGP.GE.0) THEN CALL PARCOM (DENOM) C =========== CALL PARCOM (DENUM) C =========== ENDIF C IF(IPERIO.EQ.1) THEN IDIMDE = 0 ITENDE = 0 CALL PERCOM C =========== & ( IDIMDE, ITENDE , & DENOM , DENOM , DENOM , & DENOM , DENOM , DENOM , & DENOM , DENOM , DENOM ) CALL PERCOM C =========== & ( IDIMDE, ITENDE , & DENUM , DENUM , DENUM , & DENUM , DENUM , DENUM , & DENUM , DENUM , DENUM ) ENDIF C C COMPLEMENT POUR LE VOISINAGE ETENDU C IF (IMRGRA.EQ.2.OR.IMRGRA.EQ.3) THEN C IF(IRANGP.GE.0) THEN CALL PARCVE (DENOM) C =========== CALL PARCVE (DENUM) C =========== ENDIF C IF(IPERIO.EQ.1) THEN CALL PERCVE (DENOM) C =========== CALL PERCVE (DENUM) C =========== ENDIF C ENDIF C C IF (IVECTI.EQ.1) THEN C !OCL NOVREC,VRL(16) DO IFAC = 1, NFAC C II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) C RI = 1.D0 IF( DENUM(II) .GT. CLIMGP*DENOM(II) ) THEN RI = CLIMGP*DENOM(II)/DENUM(II) ENDIF RJ = 1.D0 IF( DENUM(JJ) .GT. CLIMGP*DENOM(JJ) ) THEN RJ = CLIMGP*DENOM(JJ)/DENUM(JJ) ENDIF VECFAC = MIN(RI,RJ) C FACLIM(II) = MIN( FACLIM(II), VECFAC ) FACLIM(JJ) = MIN( FACLIM(JJ), VECFAC ) C ENDDO C ELSE C C VECTORISATION NON FORCEE DO IFAC = 1, NFAC C II = IFACEL(1,IFAC) JJ = IFACEL(2,IFAC) C RI = 1.D0 IF( DENUM(II) .GT. CLIMGP*DENOM(II) ) THEN RI = CLIMGP*DENOM(II)/DENUM(II) ENDIF RJ = 1.D0 IF( DENUM(JJ) .GT. CLIMGP*DENOM(JJ) ) THEN RJ = CLIMGP*DENOM(JJ)/DENUM(JJ) ENDIF VECFAC = MIN(RI,RJ) C FACLIM(II) = MIN( FACLIM(II), VECFAC ) FACLIM(JJ) = MIN( FACLIM(JJ), VECFAC ) C ENDDO C ENDIF C C COMPLEMENT POUR LE VOISINAGE ETENDU C PAS DE VECTORISATION PARTICULIERE A IMPOSER A PRIORI C ON NE CALCULE RIEN SUR LE HALO (ON NE S'EN SERT PAS) C NOTER QUE LE VOISINAGE ETENDU PEUT ETRE NON SYMETRIQUE C IF (IMRGRA.EQ.2.OR.IMRGRA.EQ.3) THEN C DO II = 1, NCEL VECFAC = 1.D0 DO IPOS = IPCVSE(II), IPCVSE(II+1)-1 JJ = IELVSE(IPOS) RJ = 1.D0 IF( DENUM(JJ) .GT. CLIMGP*DENOM(JJ) ) THEN RJ = CLIMGP*DENOM(JJ)/DENUM(JJ) ENDIF VECFAC = MIN(VECFAC,RJ) ENDDO FACLIM(II) = MIN( FACLIM(II), VECFAC ) ENDDO C ENDIF C C DO IEL = 1, NCEL DPDX(IEL) = FACLIM(IEL)*DPDX(IEL) DPDY(IEL) = FACLIM(IEL)*DPDY(IEL) DPDZ(IEL) = FACLIM(IEL)*DPDZ(IEL) IF( FACLIM(IEL).LT.0.99D0 ) THEN ICLIP = ICLIP+1 RMAX = MAX( RMAX, FACLIM(IEL) ) RMIN = MIN( RMIN, FACLIM(IEL) ) ENDIF ENDDO C IF (IRANGP.GE.0) THEN CALL PARMAX (RMAX) C =========== CALL PARMIN (RMIN) C =========== ENDIF C ENDIF C IF (IRANGP.GE.0) THEN CALL PARCPT (ICLIP) C =========== ENDIF C IF(IWARNP.GE.2) THEN WRITE(NFECRA,1000) ICLIP,RMIN,RMAX ENDIF C C C C C TRAITEMENT DU PARALLELISME C IF(IRANGP.GE.0) THEN CALL PARCOM (DPDX) C =========== CALL PARCOM (DPDY) C =========== CALL PARCOM (DPDZ) C =========== ENDIF C C TRAITEMENT DE LA PERIODICITE C IF(IPERIO.EQ.1) THEN CALL PERCOM C =========== & ( IDIMTE , ITENSO , & DPDX , DPDX , DPDX , & DPDY , DPDY , DPDY , & DPDZ , DPDZ , DPDZ ) ENDIF C C C-------- C FORMATS C-------- C 1000 FORMAT(1X,'LIMITATION DU GRADIENT EN ',I10,' CELLULES',/, & ' FACTEUR MINIMUM = ',E14.5,' ; FACTEUR MAXIMUM = ',E14.5) C C---- C FIN C---- C RETURN END c@z