From a42eaedf741e3962f191a7025e2bcea208c837de Mon Sep 17 00:00:00 2001
From: Gaelle DELAUTIER <gaelle.delautier@meteo.fr>
Date: Tue, 15 May 2018 14:24:39 +0200
Subject: [PATCH] J.Colin 15/5/2018 : Viscosity introduction
---
src/MNH/modd_viscosity.f90 | 48 ++++++
src/MNH/modn_viscosity.f90 | 37 ++++
src/MNH/read_desfmn.f90 | 16 ++
src/MNH/read_exsegn.f90 | 15 +-
src/MNH/viscosity.f90 | 339 +++++++++++++++++++++++++++++++++++++
5 files changed, 452 insertions(+), 3 deletions(-)
create mode 100644 src/MNH/modd_viscosity.f90
create mode 100644 src/MNH/modn_viscosity.f90
create mode 100644 src/MNH/viscosity.f90
diff --git a/src/MNH/modd_viscosity.f90 b/src/MNH/modd_viscosity.f90
new file mode 100644
index 000000000..4fb0abcf0
--- /dev/null
+++ b/src/MNH/modd_viscosity.f90
@@ -0,0 +1,48 @@
+!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
+!MNH_LIC for details. version 1.
+! #######################
+ MODULE MODD_VISC
+! #######################
+!
+!!**** *MODD_VISCOSITY* - declaration of viscosity forces constants
+!!
+!! PURPOSE
+!! -------
+! The purpose of this declarative module is to declare the
+! viscosity forces constants
+!!
+!!** IMPLICIT ARGUMENTS
+!! ------------------
+!! NONE
+!!
+!!
+!! AUTHOR
+!! ------
+!1 Jeanne Colin * Meteo-France *
+!!
+!! MODIFICATIONS
+!! -------------
+!! Original 13/04/11
+!----------------------------------------------------------------------------
+!
+!* 0. DECLARATIONS
+! ------------
+!
+IMPLICIT NONE
+!
+LOGICAL :: LVISC ! Logical switch to activate viscosity
+
+LOGICAL :: LVISC_UVW ! Logical switch to activate viscosity for the
+! momentum
+LOGICAL :: LVISC_TH ! Logical switch to activate viscosity for the
+! potential temperature
+LOGICAL :: LVISC_SV ! Logical switch to activate viscosity for the
+! scalar tracer
+LOGICAL :: LVISC_R ! Logical switch to activate viscosity for the
+! moisture
+REAL, SAVE :: XMU_v ! Molecular (cinematic) viscosity
+REAL, SAVE :: XPRANDTL ! Prandtl number
+
+END MODULE MODD_VISC
diff --git a/src/MNH/modn_viscosity.f90 b/src/MNH/modn_viscosity.f90
new file mode 100644
index 000000000..d2a9dd8f3
--- /dev/null
+++ b/src/MNH/modn_viscosity.f90
@@ -0,0 +1,37 @@
+! ###################
+ MODULE MODN_VISC
+! ###################
+!
+!!**** *MODN_VISCOSITY* - declaration of namelist NAM_VISC
+!!
+!! PURPOSE
+!! -------
+! The purpose of this module is to specify the namelist NAM_VISC
+! which concern the parameters of the viscosity forces for all models
+!
+!!
+!!** IMPLICIT ARGUMENTS
+!! ------------------
+!!
+!! REFERENCE
+!! ---------
+!!
+!! AUTHOR
+!! ------
+!! J. Colin * Meteo-France *
+!!
+!! MODIFICATIONS
+!! -------------
+!! Original April 2011
+!-------------------------------------------------------------------------------
+!
+!* 0. DECLARATIONS
+! ------------
+!
+USE MODD_VISC
+!
+IMPLICIT NONE
+!
+NAMELIST/NAM_VISC/LVISC,LVISC_UVW,LVISC_TH,LVISC_SV,LVISC_R,XMU_v,XPRANDTL
+!
+END MODULE MODN_VISC
diff --git a/src/MNH/read_desfmn.f90 b/src/MNH/read_desfmn.f90
index 98abbf158..1dc1e41e0 100644
--- a/src/MNH/read_desfmn.f90
+++ b/src/MNH/read_desfmn.f90
@@ -250,6 +250,8 @@ USE MODN_CH_ORILAM
USE MODN_DUST
USE MODN_SALT
USE MODN_PASPOL
+USE MODN_VISC
+USE MODN_DRAG_n
#ifdef MNH_FOREFIRE
USE MODN_FOREFIRE
#endif
@@ -426,6 +428,12 @@ IF (GFOUND) THEN
READ(UNIT=ILUDES,NML=NAM_CH_SOLVERn)
CALL UPDATE_NAM_CH_SOLVERn
END IF
+CALL POSNAM(ILUDES,'NAM_DRAGN',GFOUND)
+CALL INIT_NAM_DRAGn
+IF (GFOUND) THEN
+ READ(UNIT=ILUDES,NML=NAM_DRAGn)
+ CALL UPDATE_NAM_DRAGn
+END IF
CALL POSNAM(ILUDES,'NAM_SERIESN',GFOUND,ILUOUT)
CALL INIT_NAM_SERIESn
@@ -575,6 +583,8 @@ IF (KMI == 1) THEN
LTEMPDEPOS_AER(:) = LDEPOS_AER(:)
CALL POSNAM(ILUDES,'NAM_LATZ_EDFLX',GFOUND)
IF (GFOUND) READ(UNIT=ILUDES,NML=NAM_LATZ_EDFLX)
+ CALL POSNAM(ILUDES,'NAM_VISC',GFOUND,ILUOUT)
+ IF (GFOUND) READ(UNIT=ILUDES,NML=NAM_VISC)
END IF
!
!-------------------------------------------------------------------------------
@@ -669,6 +679,9 @@ IF (NVERB >= 10) THEN
!
WRITE(UNIT=ILUOUT,FMT="('********** TURBn *******************')")
WRITE(UNIT=ILUOUT,NML=NAM_TURBn)
+!
+ WRITE(UNIT=ILUOUT,FMT="('********** DRAGn *******************')")
+ WRITE(UNIT=ILUOUT,NML=NAM_DRAGn)
!
WRITE(UNIT=ILUOUT,FMT="('********** NUDGINGn ****************')")
WRITE(UNIT=ILUOUT,NML=NAM_NUDGINGn)
@@ -760,6 +773,9 @@ IF (NVERB >= 10) THEN
!
WRITE(UNIT=ILUOUT,FMT="('************ PASSIVE POLLUTANT ***************')")
WRITE(UNIT=ILUOUT,NML=NAM_PASPOL)
+!
+ WRITE(UNIT=ILUOUT,FMT="('************ VISCOSITY ***************')")
+ WRITE(UNIT=ILUOUT,NML=NAM_VISC)
!
#ifdef MNH_FOREFIRE
WRITE(UNIT=ILUOUT,FMT="('************ FOREFIRE ***************')")
diff --git a/src/MNH/read_exsegn.f90 b/src/MNH/read_exsegn.f90
index 89b59328e..ea4aa6e26 100644
--- a/src/MNH/read_exsegn.f90
+++ b/src/MNH/read_exsegn.f90
@@ -322,6 +322,7 @@ USE MODN_PARAM_ICE
USE MODN_LUNIT_n
USE MODN_NUDGING_n
USE MODN_TURB_n
+USE MODN_DRAG_n
USE MODN_BLANK
USE MODN_CH_MNHC_n
USE MODN_CH_SOLVER_n
@@ -373,6 +374,9 @@ USE MODN_CONDSAMP
USE MODN_BLOWSNOW
USE MODN_BLOWSNOW_n
USE MODN_2D_FRC
+USE MODN_VISC
+USE MODD_VISC
+USE MODD_DRAG_n
!
IMPLICIT NONE
!
@@ -463,6 +467,7 @@ CALL INIT_NAM_PARAM_MFSHALLN
CALL INIT_NAM_LBCN
CALL INIT_NAM_NUDGINGN
CALL INIT_NAM_TURBN
+CALL INIT_NAM_DRAGN
CALL INIT_NAM_CH_MNHCN
CALL INIT_NAM_CH_SOLVERN
CALL INIT_NAM_SERIESN
@@ -495,6 +500,8 @@ CALL POSNAM(ILUSEG,'NAM_NUDGINGN',GFOUND,ILUOUT)
IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_NUDGINGn)
CALL POSNAM(ILUSEG,'NAM_TURBN',GFOUND,ILUOUT)
IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_TURBn)
+CALL POSNAM(ILUSEG,'NAM_DRAGN',GFOUND,ILUOUT)
+IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_DRAGn)
CALL POSNAM(ILUSEG,'NAM_CH_MNHCN',GFOUND,ILUOUT)
IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_CH_MNHCn)
CALL POSNAM(ILUSEG,'NAM_CH_SOLVERN',GFOUND,ILUOUT)
@@ -647,6 +654,8 @@ IF (KMI == 1) THEN
IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_LATZ_EDFLX)
CALL POSNAM(ILUSEG,'NAM_BLOWSNOW',GFOUND,ILUOUT)
IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_BLOWSNOW)
+ CALL POSNAM(ILUSEG,'NAM_VISC',GFOUND,ILUOUT)
+ IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_VISC)
END IF
!
!-------------------------------------------------------------------------------
@@ -669,7 +678,7 @@ CALL TEST_NAM_VAR(ILUOUT,'CRAD',CRAD,'NONE','FIXE','ECMW',&
#endif
'TOPA')
CALL TEST_NAM_VAR(ILUOUT,'CCLOUD',CCLOUD,'NONE','REVE','KESS', &
- & 'ICE2','ICE3','ICE4','C2R2','C3R5','KHKO','LIMA')
+ & 'ICE3','ICE4','C2R2','C3R5','KHKO','LIMA')
CALL TEST_NAM_VAR(ILUOUT,'CDCONV',CDCONV,'NONE','KAFR')
CALL TEST_NAM_VAR(ILUOUT,'CSCONV',CSCONV,'NONE','KAFR','EDKF')
CALL TEST_NAM_VAR(ILUOUT,'CELEC',CELEC,'NONE','ELE3','ELE4')
@@ -690,7 +699,7 @@ CALL TEST_NAM_VAR(ILUOUT,'CLBCY(1)',CLBCY(1),'CYCL','WALL','OPEN')
CALL TEST_NAM_VAR(ILUOUT,'CLBCY(2)',CLBCY(2),'CYCL','WALL','OPEN')
!
CALL TEST_NAM_VAR(ILUOUT,'CTURBDIM',CTURBDIM,'1DIM','3DIM')
-CALL TEST_NAM_VAR(ILUOUT,'CTURBLEN',CTURBLEN,'DELT','BL89','DEAR','BLKR')
+CALL TEST_NAM_VAR(ILUOUT,'CTURBLEN',CTURBLEN,'DELT','BL89','RM17','DEAR','BLKR')
CALL TEST_NAM_VAR(ILUOUT,'CTOM',CTOM,'NONE','TM06')
CALL TEST_NAM_VAR(ILUOUT,'CSUBG_AUCV',CSUBG_AUCV,'NONE','CLFR','SIGM')
!
@@ -718,7 +727,7 @@ CALL TEST_NAM_VAR(ILUOUT,'CMF_CLOUD',CMF_CLOUD,'NONE','STAT','DIRE')
! The test on the CSOLVER name is made elsewhere
!
CALL TEST_NAM_VAR(ILUOUT,'CPRISTINE_ICE',CPRISTINE_ICE,'PLAT','COLU','BURO')
-CALL TEST_NAM_VAR(ILUOUT,'CSEDIM',CSEDIM,'SPLI','STAT')
+CALL TEST_NAM_VAR(ILUOUT,'CSEDIM',CSEDIM,'SPLI','STAT','NONE')
IF( CCLOUD == 'C3R5' ) THEN
CALL TEST_NAM_VAR(ILUOUT,'CPRISTINE_ICE_C1R3',CPRISTINE_ICE_C1R3, &
'PLAT','COLU','BURO')
diff --git a/src/MNH/viscosity.f90 b/src/MNH/viscosity.f90
new file mode 100644
index 000000000..ae03806c0
--- /dev/null
+++ b/src/MNH/viscosity.f90
@@ -0,0 +1,339 @@
+!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
+!MNH_LIC for details. version 1.
+!
+! #####################
+ MODULE MODI_VISC
+! #####################
+!
+INTERFACE
+!
+ SUBROUTINE VISCOSITY(HLBCX, HLBCY, KRR, KSV, PNU, PPRANDTL, &
+ OVISC_UVW, OVISC_TH, OVISC_SV, OVISC_R, &
+ ODRAG, &
+ PUT, PVT, PWT, PTHT, PRT, PSVT, &
+ PRHODJ, PDXX, PDYY, PDZZ, PDZX, PDZY, &
+ PRUS, PRVS, PRWS, PRTHS, PRRS, PRSVS,PDRAG )
+!
+ IMPLICIT NONE
+!
+!* 0.1 Declarations of arguments:
+!
+! X and Y lateral boundary conditions
+ CHARACTER(LEN=4),DIMENSION(2),INTENT(IN):: HLBCX, HLBCY
+!
+ INTEGER, INTENT(IN) :: KRR ! number of moist variables
+ INTEGER, INTENT(IN) :: KSV ! number of scalar variables
+!
+ REAL, INTENT(IN) :: PNU ! viscosity coefficient
+ REAL, INTENT(IN) :: PPRANDTL ! Parandtl number
+!
+!
+! logical switches
+ LOGICAL, INTENT(IN) :: OVISC_UVW ! momentum
+ LOGICAL, INTENT(IN) :: OVISC_TH ! theta
+ LOGICAL, INTENT(IN) :: OVISC_SV ! scalar tracer
+ LOGICAL, INTENT(IN) :: OVISC_R ! moisture
+ LOGICAL, INTENT(IN) :: ODRAG ! noslip/freeslip
+!
+!
+! input variables at time t
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PVT
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PWT
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PTHT
+ REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRT
+ REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PSVT
+!
+!
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! rho_ref * Jacobian
+!
+ REAL, DIMENSION(:,:), INTENT(IN) :: PDRAG ! Array -1/1 defining where the no-slipcondition is applied
+! metric coefficients
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDYY
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDZZ
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDZX
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDZY
+!
+! output source terms
+ REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRUS, PRVS, PRWS
+ REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRTHS
+ REAL, DIMENSION(:,:,:,:), INTENT(INOUT) :: PRRS, PRSVS
+!
+ END SUBROUTINE VISCOSITY
+!
+END INTERFACE
+!
+END MODULE MODI_VISC
+!
+!-------------------------------------------------------------------------------
+!
+SUBROUTINE VISCOSITY(HLBCX, HLBCY, KRR, KSV, PNU, PPRANDTL, &
+ OVISC_UVW, OVISC_TH, OVISC_SV, OVISC_R, &
+ ODRAG, &
+ PUT, PVT, PWT, PTHT, PRT, PSVT, &
+ PRHODJ, PDXX, PDYY, PDZZ, PDZX, PDZY, &
+ PRUS, PRVS, PRWS, PRTHS, PRRS, PRSVS,PDRAG )
+!
+! IMPLICIT ARGUMENTS
+! ------------------
+! Module MODD_PARAMETERS: JPHEXT, JPVEXT
+! Module MODD_CONF: L2D
+!
+!! AUTHOR
+!! ------
+!! Jeanne Colin * CNRM *
+!!
+!! MODIFICATIONS
+!! -------------
+!! 01/18 (C.Lac) Add budgets
+!-------------------------------------------------------------------------------
+!
+!* 0. DECLARATIONS
+! ------------
+!
+ USE MODI_LAP_M
+ USE MODI_SHUMAN
+ USE MODD_PARAMETERS
+ USE MODD_CONF
+ USE MODD_VISC
+ USE MODD_DRAG_n
+ USE MODD_BUDGET
+ USE MODE_ll
+ USE MODD_ARGSLIST_ll, ONLY : LIST_ll
+ USE MODE_FM
+ USE MODI_BUDGET
+!
+!-------------------------------------------------------------------------------
+!
+ IMPLICIT NONE
+!
+!* 0.1 Declarations of dummy arguments :
+!
+!
+! X and Y lateral boundary conditions
+ CHARACTER(LEN=4),DIMENSION(2),INTENT(IN):: HLBCX, HLBCY
+!
+ INTEGER, INTENT(IN) :: KRR ! number of moist variables
+ INTEGER, INTENT(IN) :: KSV ! number of scalar variables
+!
+ REAL, INTENT(IN) :: PPRANDTL ! Parandtl number
+ REAL, INTENT(IN) :: PNU ! viscous diffusion rate
+!
+! logical switches
+ LOGICAL, INTENT(IN) :: OVISC_UVW ! momentum
+ LOGICAL, INTENT(IN) :: OVISC_TH ! theta
+ LOGICAL, INTENT(IN) :: OVISC_SV ! scalar tracer
+ LOGICAL, INTENT(IN) :: OVISC_R ! moisture
+ LOGICAL, INTENT(IN) :: ODRAG ! moisture
+!
+! input variables at time t
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PVT
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PWT
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PTHT
+ REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRT
+ REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PSVT
+!
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! rho_ref * Jacobian
+!
+!
+REAL, DIMENSION(:,:), INTENT(IN) :: PDRAG ! Array -1/1 defining where the no-slip condition is applied
+
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDYY
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDZZ
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDZX
+ REAL, DIMENSION(:,:,:), INTENT(IN) :: PDZY
+!
+! output source terms
+ REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRUS, PRVS, PRWS
+ REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRTHS
+ REAL, DIMENSION(:,:,:,:), INTENT(INOUT) :: PRRS, PRSVS
+!
+!
+!* 0.2 Declarations of local variables
+!
+ INTEGER :: IK ! counter
+ INTEGER :: IKB, IKE
+!
+ REAL, DIMENSION(SIZE(PWT,1),SIZE(PWT,2),SIZE(PWT,3)) :: ZTMP ! temp storage
+ REAL, DIMENSION(SIZE(PWT,1),SIZE(PWT,2),SIZE(PWT,3)) :: ZLAPu ! temp storage
+ REAL, DIMENSION(SIZE(PWT,1),SIZE(PWT,2),SIZE(PWT,3)) :: ZLAPv ! temp storage
+ REAL, DIMENSION(SIZE(PWT,1),SIZE(PWT,2),SIZE(PWT,3)) :: ZY1 ! temp storage
+ REAL, DIMENSION(SIZE(PWT,1),SIZE(PWT,2),SIZE(PWT,3)) :: ZY2 ! temp storage
+!
+!
+INTEGER :: IIU,IJU,IKU ! I,J,K array sizes
+!
+INTEGER :: JI,JJ,JK ! I loop index
+INTEGER :: IINFO_ll
+TYPE(LIST_ll), POINTER :: TZFIELDS_ll ! list of fields to exchange
+!
+!
+!-------------------------------------------------------------------------------
+!
+IIU=SIZE(PWT,1)
+IJU=SIZE(PWT,2)
+IKU=SIZE(PWT,3)
+
+!* 1. Viscous forcing for potential temperature
+! -----------------------------------------
+!
+!
+IF (OVISC_TH) THEN
+!
+!
+ PRTHS = PRTHS + PNU/PPRANDTL * &
+ LAP_M(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,PTHT)
+!
+!
+ END IF
+!
+IF (LBUDGET_TH) CALL BUDGET (PRTHS,4,'VISC_BU_RU')
+!
+!-------------------------------------------------------------------------------
+!
+!* 2. Viscous forcing for moisture
+! ----------------------------
+!
+ IF (OVISC_R .AND. (SIZE(PRT,1) > 0)) THEN
+!
+!
+ DO IK = 1, KRR
+ PRRS(:,:,:,IK) = PRRS(:,:,:,IK) + PNU/PPRANDTL * &
+ LAP_M(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,PRT(:,:,:,IK))
+ END DO
+!
+!
+ END IF
+!
+IF (LBUDGET_RV) CALL BUDGET (PRRS(:,:,:,1),6,'VISC_BU_RRV')
+IF (LBUDGET_RC) CALL BUDGET (PRRS(:,:,:,2),7,'VISC_BU_RRC')
+IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:,3),8,'VISC_BU_RRR')
+IF (LBUDGET_RI) CALL BUDGET (PRRS(:,:,:,4),9,'VISC_BU_RRI')
+IF (LBUDGET_RS) CALL BUDGET (PRRS(:,:,:,5),10,'VISC_BU_RRS')
+IF (LBUDGET_RG) CALL BUDGET (PRRS(:,:,:,6),11,'VISC_BU_RRG')
+IF (LBUDGET_RH) CALL BUDGET (PRRS(:,:,:,7),12,'VISC_BU_RRH')
+!
+!-------------------------------------------------------------------------------
+!
+!* 3. Viscous forcing for passive scalars
+! -----------------------------------
+!
+ IF (OVISC_SV .AND. (SIZE(PSVT,1) > 0)) THEN
+!
+!
+ DO IK = 1, KSV
+ PRSVS(:,:,:,IK) = PRSVS(:,:,:,IK) + PNU/PPRANDTL * &
+ LAP_M(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,PSVT(:,:,:,IK))
+ END DO
+!
+ END IF
+!
+IF (LBUDGET_SV) THEN
+ DO IK = 1, KSV
+ CALL BUDGET (PRSVS(:,:,:,IK), 12+IK, 'VISC_BU_RSV')
+ END DO
+END IF
+!
+
+!-------------------------------------------------------------------------------
+!
+!* 4. Viscous forcing for momentum
+! ----------------------------
+!
+IF (OVISC_UVW) THEN
+!
+!* 4.1 U - component
+!
+!
+ ZY1 = MXF(PUT)
+ IF (ODRAG) THEN
+ ZY1(:,:,1) = PDRAG * ZY1(:,:,2)
+ ENDIF
+!
+!
+ ZLAPu = LAP_M(HLBCX,HLBCY,PDXX,PDYY,PDZX, &
+ PDZY,PDZZ,PRHODJ,ZY1)
+!! Update halo to compute the source term
+NULLIFY(TZFIELDS_ll)
+CALL ADD3DFIELD_ll(TZFIELDS_ll,ZLAPu)
+CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll)
+CALL CLEANLIST_ll(TZFIELDS_ll)
+!
+PRUS = PRUS + MXM(PNU*ZLAPu)
+!
+!* 4.2 V - component
+! -------------
+!
+ IF (.NOT. L2D) THEN
+
+ ZY2 = MYF(PVT)
+ IF (ODRAG) THEN
+ ZY2(:,:,1) = PDRAG * ZY2(:,:,2)
+ ENDIF
+ ELSE
+!
+ ZLAPv = LAP_M(HLBCX,HLBCY,PDXX,PDYY,PDZX, &
+ PDZY,PDZZ,PRHODJ,ZY2)
+!! Update halo to compute the source term
+!
+NULLIFY(TZFIELDS_ll)
+CALL ADD3DFIELD_ll(TZFIELDS_ll,ZLAPv)
+CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll)
+CALL CLEANLIST_ll(TZFIELDS_ll)
+!
+PRVS = PRVS + MYM(PNU*ZLAPv)
+
+ENDIF
+
+!
+!* 4.3 W - component
+! -------------
+!
+ IKB = JPVEXT + 1
+ IKE = SIZE(PWT,3) - JPVEXT
+
+ ZTMP = PWT
+!
+IF (ODRAG) THEN
+ WHERE (PDRAG==-1)
+ ZTMP(:,:,IKB) = 0.
+ ENDWHERE
+ENDIF
+!
+ DO IK = 1,JPVEXT
+ ZTMP(:,:,IK) = ZTMP(:,:,IKB)
+ ZTMP(:,:,IKE+IK) = ZTMP(:,:,IKE)
+ END DO
+!
+ ZTMP = MZM(1,IKU,1, PNU * &
+ LAP_M(HLBCX,HLBCY,PDXX,PDYY,PDZX,PDZY,PDZZ,PRHODJ,ZTMP) )
+
+ DO IK = 1,JPVEXT
+ ZTMP(:,:,IK) = ZTMP(:,:,IKB)
+ ZTMP(:,:,IKE+IK) = ZTMP(:,:,IKE)
+ END DO
+PRWS = PRWS + ZTMP
+!
+!!! Debug provisoire dans le cas ou le noslip est applique jusqu'au bord de
+!sortie de flux en OPEN
+ IF ( LWEST_ll().AND.(ODRAG).AND.(MINVAL(PDRAG(IIU,:))== -1)) THEN
+ DO JK=1,IKU
+ WHERE(PDRAG(IIU,:)== -1)
+ PRUS(IIU,:,JK) = PRUS(IIU-1,:,JK)
+ PRVS(IIU,:,JK) = PRVS(IIU-1,:,JK)
+ PRWS(IIU,:,JK) = PRWS(IIU-1,:,JK)
+ ENDWHERE
+ END DO
+ ENDIF
+END IF
+!
+IF (LBUDGET_U) CALL BUDGET (PRUS,1,'VISC_BU_RU')
+IF (LBUDGET_V) CALL BUDGET (PRVS,2,'VISC_BU_RV')
+IF (LBUDGET_W) CALL BUDGET (PRWS,2,'VISC_BU_RW')
+!
+END SUBROUTINE VISCOSITY
--
GitLab