kepcl3.f

Go to the documentation of this file.
00001 C:\opentelemac\v7p0\sources\telemac3d\kepcl3.f
00002 !
00073                      SUBROUTINE KEPCL3
00074 !                    *****************
00075 !
00076      &(KBORF,EBORF,LIKBOF,LIEBOF,LIUBOF,
00077      & KBORL,EBORL,LIKBOL,LIEBOL,LIUBOL,RUGOL,
00078      & KBORS,EBORS,LIKBOS,LIEBOS,LIUBOS,
00079      & DISBOR,AK,U,V,H,Z,NBOR,NPOIN2,NPLAN,NPTFR,
00080      & DNUVIH,DNUVIV,KARMAN,CMU,LISRUF,LISRUL,
00081      & VIRT,KMIN,KMAX,EMIN,EMAX,KENT,KENTU,KSORT,KADH,KLOG,
00082      & UETCAR,UETCAL,FICTIF)
00083 !
00084 !***********************************************************************
00085 ! TELEMAC3D   V7P0                                   21/08/2010
00086 !***********************************************************************
00087 !
00088 !
00089 !
00090 !
00091 !
00092 !
00093 !
00094 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00095 !| AK             |-->| TURBULENT ENERGY
00096 !| CMU            |-->| CONSTANT FOR MODELE K-EPSILON MODEL
00097 !| DISBOR         |-->| DISTANCE TO BOUNDARY OF POINTS CLOSE TO BOUNDARY
00098 !| DNUVIH         |-->| COEFFICIENT FOR HORIZONTAL DIFFUSION OF VELOCITIES
00099 !| DNUVIV         |-->| COEFFICIENT FOR VERTICAL DIFFUSION OF VELOCITIES
00100 !| EBORF          |<->| EPSILON ON BOTTOM
00101 !| EBORL          |<->| EPSILON ON LATERAL SOLID BOUNDARIES
00102 !| EBORS          |<->| EPSILON AT SURFACE
00103 !| EMAX           |-->| MAXIMUM VALUE FOR EPSILON WHEN CLIPPING
00104 !| EMIN           |-->| MINIMUM VALUE FOR EPSILON WHEN CLIPPING
00105 !| H              |-->| WATER DEPTH AT TIME N
00106 !| KADH           |-->| CONVENTION FOR NO SLIP BOUNDARY CONDITION
00107 !| KARMAN         |-->| KARMAN CONSTANT
00108 !| KBORF          |<->| K ON BOTTOM
00109 !| KBORL          |<->| K ON LATERAL SOLID BOUNDARIES
00110 !| KBORS          |<->| K AT SURFACE
00111 !| KENT           |-->| CONVENTION FOR LIQUID INPUT WITH PRESCRIBED VALUE
00112 !| KENTU          |-->| CONVENTION FOR LIQUID INPUT WITH PRESCRIBED VELOCITY
00113 !| KLOG           |-->| CONVENTION FOR LOGARITHMIC WALL
00114 !| KMAX           |-->| MAXIMUM VALUE FOR K WHEN CLIPPING
00115 !| KMIN           |-->| MINIMUM VALUE FOR K WHEN CLIPPING
00116 !| KSORT          |-->| CONVENTION FOR FREE OUTPUT
00117 !| LIEBOF         |<->| TYPE OF BOUNDARY CONDITIONS ON EPSILON AT THE BOTTOM
00118 !| LIEBOL         |<->| TYPE OF BOUNDARY CONDITIONS ON EPSILON ON THE LATERAL WALLS
00119 !| LIEBOS         |<->| TYPE OF BOUNDARY CONDITIONS ON EPSILON AT THE SURFACE
00120 !| LIKBOF         |<->| TYPE OF BOUNDARY CONDITIONS ON K AT THE BOTTOM
00121 !| LIKBOL         |<->| TYPE OF BOUNDARY CONDITIONS ON K ON THE LATERAL WALLS
00122 !| LIKBOS         |<->| TYPE OF BOUNDARY CONDITIONS ON K AT THE SURFACE
00123 !| LISRUF         |-->| TURBULENCE MODEL FOR BOTTOM
00124 !|                |   | 1: SMOOTH  2: ROUGH  3: ROUGH (CHEZY)
00125 !| LISRUL         |-->| TURBULENCE MODEL FOR SOLID BOUNDARIES
00126 !|                |   | 1: SMOOTH  2: ROUGH  3: ROUGH (CHEZY)
00127 !| LIUBOF         |<->| TYPE OF BOUNDARY CONDITIONS ON U AT THE BOTTOM
00128 !| LIUBOL         |<->| TYPE OF BOUNDARY CONDITIONS ON U ON THE LATERAL WALLS
00129 !| LIUBOS         |<->| TYPE OF BOUNDARY CONDITIONS ON U AT THE SURFACE
00130 !| NBOR           |-->| GLOBAL NUMBER OF BOUNDARY POINTS
00131 !| NPLAN          |-->| NUMBER OF PLANES IN THE 3D MESH OF PRISMS
00132 !| NPOIN2         |-->| NUMBER OF POINTS IN 2D
00133 !| NPTFR          |-->| NUMBER OF BOUNDARY POINTS IN 2D
00134 !| RUGOL          |-->| NOT USED
00135 !| U              |-->| X COMPONENT OF VELOCITY AT TIME N
00136 !| UETCAL         |-->| (UETUTA*UTANG(IPTFR))**2: IN COMMENT
00137 !| UETCAR         |-->| USTAR**2
00138 !| V              |-->| Y COMPONENT OF VELOCITY AT TIME N
00139 !| VIRT           |-->| VIRTUAL ORIGIN FOR EPSILON (TELEMAC 3D): COMMENT
00140 !| Z              |-->| ELEVATION OF REAL 3D MESH POINTS
00141 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00142 !
00143       USE DECLARATIONS_TELEMAC3D, ONLY: IPBOT,AEBORF,BEBORF,SIGMAE,RUGOF
00144 !
00145       IMPLICIT NONE
00146 !
00147       INTEGER LNG,LU
00148       COMMON/INFO/LNG,LU
00149 !
00150 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00151 !
00152       INTEGER, INTENT(IN) :: NPTFR, NPLAN, NPOIN2,KENTU
00153       INTEGER, INTENT(IN) :: LISRUF,LISRUL,KENT,KSORT,KADH,KLOG
00154 !
00155       INTEGER, INTENT(INOUT) :: LIKBOF(NPOIN2), LIKBOS(NPOIN2)
00156       INTEGER, INTENT(INOUT) :: LIKBOL(NPTFR,NPLAN)
00157       INTEGER, INTENT(INOUT) :: LIEBOF(NPOIN2), LIEBOS(NPOIN2)
00158       INTEGER, INTENT(INOUT) :: LIEBOL(NPTFR,NPLAN)
00159       INTEGER, INTENT(INOUT) :: LIUBOF(NPOIN2), LIUBOS(NPOIN2)
00160       INTEGER, INTENT(INOUT) :: LIUBOL(NPTFR,NPLAN)
00161 !
00162       INTEGER, INTENT(IN) :: NBOR(NPTFR)
00163 !
00164       DOUBLE PRECISION, INTENT(IN) :: U(NPOIN2,NPLAN), V(NPOIN2,NPLAN)
00165       DOUBLE PRECISION, INTENT(IN) :: Z(NPOIN2,NPLAN), AK(NPOIN2,NPLAN)
00166       DOUBLE PRECISION, INTENT(IN) :: H(NPOIN2)      , UETCAR(NPOIN2)
00167 !
00168       DOUBLE PRECISION, INTENT(INOUT) :: KBORF(NPOIN2), KBORS(NPOIN2)
00169       DOUBLE PRECISION, INTENT(INOUT) :: KBORL(NPTFR,NPLAN)
00170       DOUBLE PRECISION, INTENT(INOUT) :: EBORF(NPOIN2), EBORS(NPOIN2)
00171       DOUBLE PRECISION, INTENT(INOUT) :: EBORL(NPTFR,NPLAN)
00172       DOUBLE PRECISION, INTENT(IN)    :: UETCAL(NPTFR,NPLAN)
00173 !
00174       DOUBLE PRECISION, INTENT(IN)    :: DISBOR(NPTFR)
00175       DOUBLE PRECISION, INTENT(IN)    :: RUGOL(NPTFR,NPLAN)
00176 !
00177       DOUBLE PRECISION, INTENT(IN) :: VIRT, DNUVIH, DNUVIV, KARMAN
00178       DOUBLE PRECISION, INTENT(IN) :: CMU
00179       DOUBLE PRECISION, INTENT(IN) :: KMIN, KMAX, EMIN, EMAX
00180       DOUBLE PRECISION, INTENT(IN) :: FICTIF
00181 !
00182 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
00183 !
00184       INTRINSIC EXP
00185 !
00186 !-----------------------------------------------------------------------
00187 !
00188       INTEGER IPTFR,IPLAN,IPOIN2,IP,IBOT
00189 !
00190       DOUBLE PRECISION ESURF,HAUT,Z0,SSQCMU,DIST,PROPNU,DISTFOND
00191 !
00192       DOUBLE PRECISION, PARAMETER :: NIVTURB = 0.005D0
00193       DOUBLE PRECISION, PARAMETER :: TESTREICH = 1.D-4
00194       INTEGER, PARAMETER :: MAXITEREICH = 30
00195 !
00196       INTRINSIC SQRT,MAX,LOG
00197 !
00198 !-----------------------------------------------------------------------
00199 !
00200       SSQCMU = 1.D0 /SQRT(CMU)
00201       PROPNU = (2*DNUVIH + DNUVIV) /3.D0
00202 !
00203 !=======================================================================
00204 !     BOTTOM
00205 !=======================================================================
00206 !
00207 !     THIS IS DONE IN LIMI3D
00208 !     AEBORF%TYPR='0'
00209 !     BEBORF%TYPR='0'
00210 !     AKBORF%TYPR='0'
00211 !     BKBORF%TYPR='0'
00212 !
00213       DO IPOIN2=1,NPOIN2
00214 !
00215         IF(IPBOT%I(IPOIN2).EQ.0) THEN
00216 !         NORMAL CASE
00217           DIST =(Z(IPOIN2,2)-Z(IPOIN2,1))/FICTIF
00218           IF(LIEBOF(IPOIN2).EQ.KENT) THEN
00219             EBORF(IPOIN2)=MAX(UETCAR(IPOIN2)*SQRT(UETCAR(IPOIN2))
00220      &                                             /(KARMAN*DIST),EMIN)
00221           ENDIF
00222           IF(LIKBOF(IPOIN2).EQ.KENT) THEN
00223             KBORF(IPOIN2) = MAX(SSQCMU*UETCAR(IPOIN2),KMIN)
00224           ENDIF
00225         ELSE
00226 !         RISK OF SMASHED PLANES OR TIDAL FLATS
00227           IPLAN=IPBOT%I(IPOIN2)+1
00228           IF(IPLAN.EQ.NPLAN) THEN
00229 !           CASE OF TIDAL FLATS
00230             IF(LIEBOF(IPOIN2).EQ.KENT) THEN
00231               EBORF(IPOIN2)=EMIN
00232             ENDIF
00233 !           IN THIS CASE KBORF COMPUTED ABOVE MAY YIELD
00234 !           ABNORMAL VALUES OF VISCOSITY
00235             IF(LIKBOF(IPOIN2).EQ.KENT) KBORF(IPOIN2)=KMIN
00236           ELSE
00237 !           CASE OF SMASHED PLANES : DIST COMPUTED ON FIRST FREE LAYER
00238             DIST =(Z(IPOIN2,IPLAN+1)-Z(IPOIN2,IPLAN))/FICTIF
00239             IF(LIEBOF(IPOIN2).EQ.KENT) THEN
00240               EBORF(IPOIN2)=MAX(UETCAR(IPOIN2)*
00241      &                         SQRT(UETCAR(IPOIN2))/(KARMAN*DIST),EMIN)
00242             ENDIF
00243             IF(LIKBOF(IPOIN2).EQ.KENT) THEN
00244               KBORF(IPOIN2) = MAX(SSQCMU*UETCAR(IPOIN2),KMIN)
00245             ENDIF
00246           ENDIF
00247         ENDIF
00248 !
00249       ENDDO
00250 !
00251 !=======================================================================
00252 !     FREE SURFACE
00253 !=======================================================================
00254 !
00255       DO IPOIN2=1,NPOIN2
00256 !
00257 !       DIRICHLET ON EPSILON
00258 !       ---------------------
00259 !
00260         IF(LIEBOS(IPOIN2).EQ.KENT) THEN
00261 !
00262 !         NEZU & NAKAGAWA: TURBULENCE IN OPEN CHANNEL FLOWS
00263 !
00264 !         SEE KEPICL: THIS VALUE IS NOT USED IF CONDITION
00265 !                     IS NOT SET TO KENT ON EPSILON
00266           HAUT = MAX(H(IPOIN2),1.D-7)
00267           ESURF = SQRT(AK(IPOIN2,NPLAN)**3) / (0.18D0*HAUT)
00268 !
00269 !         5.9
00270 !         ESURF = CMU**0.75D0 / KARMAN
00271 !    &          * AK(IPOIN2,NPLAN)**1.5D0 / (VIRT*HAUT)
00272           EBORS(IPOIN2) = MAX(ESURF,EMIN)
00273 !
00274         ENDIF
00275 !
00276       ENDDO
00277 !
00278 !=======================================================================
00279 !     LATERAL BOUNDARIES
00280 !=======================================================================
00281 !
00282       DO IPTFR=1,NPTFR
00283 !
00284         IPOIN2 = NBOR(IPTFR)
00285         DIST   = DISBOR(IPTFR) / FICTIF
00286         HAUT   = MAX(H(IPOIN2),1.D-7)
00287 !
00288         Z0=RUGOF%R(IPOIN2)/30.D0
00289 !
00290         DO IPLAN=1,NPLAN
00291 !
00292           IP=MAX(IPLAN,2)
00293           IBOT=MIN(IPBOT%I(IPOIN2)+1,NPLAN-1)
00294 !         DISTANCE TO BOTTOM (WILL BE 0 WITH TIDAL FLATS)
00295           DISTFOND = (Z(IPOIN2,IP)-Z(IPOIN2,IBOT))
00296 !
00297 !         DIRICHLET ON K
00298 !         ---------------
00299 !
00300           IF(LIKBOL(IPTFR,IPLAN).EQ.KENT) THEN
00301 !         ------------------------------------
00302 !
00303 !           ************************************
00304             IF(LIUBOL(IPTFR,IPLAN).EQ.KENT.OR.
00305      &         LIUBOL(IPTFR,IPLAN).EQ.KENTU.OR.
00306      &         LIUBOL(IPTFR,IPLAN).EQ.KSORT     ) THEN
00307 !           ************************************
00308 !
00309 !          THEORY BY VINCENT BOYER (SEE ALSO KEPINI)
00310 !
00311 !          KBORL(IPTFR,IPLAN) = MAX(NIVTURB*U(IPOIN2,IPLAN)**2,KMIN)
00312 !
00313 !          NO TURBULENCE
00314 !
00315 !           KBORL(IPTFR,IPLAN) = KMIN
00316 !           CV  HANS AND BURCHARD CL FOR K
00317             KBORL(IPTFR,IPLAN)=UETCAR(IPOIN2)
00318      &                        *(1.D0-DISTFOND/HAUT)/SQRT(CMU)
00319 !
00320 !            ****************************************
00321             ELSEIF(LIUBOL(IPTFR,IPLAN).EQ.KLOG .OR.
00322      &             LIUBOL(IPTFR,IPLAN).EQ.KADH) THEN
00323 !           ****************************************
00324 !
00325 !             WALL
00326 !
00327 !             KBORL(IPTFR,IPLAN)=MAX(SSQCMU*UETCAL(IPTFR,IPLAN),KMIN)
00328               KBORL(IPTFR,IPLAN)=KMIN
00329 !
00330 !           ****
00331             ELSE
00332 !           ****
00333 !
00334               IF (LNG.EQ.1) WRITE(LU,111) IPTFR,LIUBOL(IPTFR,IPLAN)
00335               IF (LNG.EQ.2) WRITE(LU,112) IPTFR,LIUBOL(IPTFR,IPLAN)
00336               CALL PLANTE(1)
00337               STOP
00338 !
00339 !           *****
00340             ENDIF
00341 !           *****
00342 !
00343           ENDIF
00344 !         -----
00345 !
00346 !         DIRICHLET ON EPSILON
00347 !         ---------------------
00348 !
00349           IF(LIEBOL(IPTFR,IPLAN).EQ.KENT) THEN
00350 !         ------------------------------------
00351 !
00352 !           ************************************
00353             IF(LIUBOL(IPTFR,IPLAN).EQ.KENT.OR.
00354      &         LIUBOL(IPTFR,IPLAN).EQ.KENTU.OR.
00355      &         LIUBOL(IPTFR,IPLAN).EQ.KSORT     ) THEN
00356 !           ************************************
00357 !
00358 !              COMING IN THE DOMAIN: TURBULENCE DUE TO THE
00359 !              BOTTOM AS IN KEPINI; COMPUTES EBORL ACCORDING
00360 !              TO KBORL AT THE BOTTOM
00361 !
00362 !              EBORL(IPTFR,IPLAN)=CMU**0.75*SQRT(KBORL(IPTFR,1)**3)
00363 !    &                           /KARMAN/MAX(DISTFOND,1.D-6)
00364 !              EBORL(IPTFR,IPLAN)= MAX(EBORL(IPTFR,IPLAN),EMIN)
00365 !               EBORL(IPTFR,IPLAN)=EMIN
00366 !             Hans et Burchard
00367 !             CV ...
00368               EBORL(IPTFR,IPLAN)=SQRT(UETCAR(IPOIN2))**3
00369      &                          *(1.D0-DISTFOND/HAUT)
00370      &                          /KARMAN/MAX(DISTFOND,Z0)
00371 !
00372 !           ****************************************
00373             ELSEIF(LIUBOL(IPTFR,IPLAN).EQ.KLOG .OR.
00374      &             LIUBOL(IPTFR,IPLAN).EQ.KADH) THEN
00375 !           ****************************************
00376 !
00377 !             WALL
00378 !
00379 !             EBORL(IPTFR,IPLAN) =
00380 !    &        MAX(UETCAL(IPTFR,IPLAN)*SQRT(UETCAL(IPTFR,IPLAN))/
00381 !    &        (KARMAN*DIST*FICTIFUET/FICTIFEPS),EMIN)
00382               EBORL(IPTFR,IPLAN)=EMIN
00383 !
00384 !           ****
00385             ELSE
00386 !           ****
00387 !
00388 !             OTHER
00389 !
00390               IF (LNG.EQ.1) WRITE(LU,121) IPTFR,LIUBOL(IPTFR,IPLAN)
00391               IF (LNG.EQ.2) WRITE(LU,122) IPTFR,LIUBOL(IPTFR,IPLAN)
00392               CALL PLANTE(1)
00393               STOP
00394 !
00395 !           *****
00396             ENDIF
00397 !           *****
00398 !
00399           ENDIF
00400 !          -----
00401 !
00402         ENDDO
00403       ENDDO
00404 !
00405 !-----------------------------------------------------------------------
00406 !
00407 111   FORMAT(' KEPCL3 : POINT DE BORD',I6,
00408      &       ' - CAS NON PREVU POUR KBOR : LIUBOR =',I6)
00409 112   FORMAT(' KEPCL3 : BOUNDARY NODE',I6,
00410      &       ' - UNEXPECTED CONDITION FOR KBOR : LIUBOR =',I6)
00411 121   FORMAT(' KEPCL3 : POINT DE BORD',I6,
00412      &       ' - CAS NON PREVU POUR EBOR : LIUBOR =',I6)
00413 122   FORMAT(' KEPCL3 : BOUNDARY NODE',I6,
00414      &       ' - UNEXPECTED CONDITION FOR EBOR : LIUBOR =',I6)
00415 !
00416 !-----------------------------------------------------------------------
00417 !
00418       RETURN
00419       END

Generated on Fri Aug 31 2013 18:12:58 by S.E.Bourban (HRW) using doxygen 1.7.0