Newer
Older
!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.
!-----------------------------------------------------------------
!--------------- special set of characters for RCS information
!-----------------------------------------------------------------
! $Source$ $Revision$
! MASDEV4_7 lbc 2006/06/06 15:29:52
!-----------------------------------------------------------------
!####################
MODULE MODI_RAD_BOUND
!####################
!
INTERFACE
!
SUBROUTINE RAD_BOUND (HLBCX,HLBCY,HTURB,PRIMKMAX, &
PTSTEP,PDXHAT,PDYHAT,PZHAT, &
PLBXUM,PLBYVM,PLBXUS,PLBYVS, &
PCPHASE,PCPHASE_PBL,PRHODJ, &
!
CHARACTER(LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX,HLBCY ! X and Y-direc. LBC type
CHARACTER(LEN=4), INTENT(IN) :: HTURB ! Turbulence scheme
!
!
REAL, INTENT(IN) :: PRIMKMAX ! Rayleigh damping amplitude
REAL, INTENT(IN) :: PTSTEP ! time step dt
REAL, DIMENSION(:), INTENT(IN) :: PDXHAT ! X-direc. meshlength
REAL, DIMENSION(:), INTENT(IN) :: PDYHAT ! Y-direc. meshlength
REAL, DIMENSION(:), INTENT(IN) :: PZHAT ! height level without orography
!
REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT,PVT ! at t
!
! Lateral Boundary fields at time t
REAL, DIMENSION(:,:,:), INTENT(IN) :: PLBXUM,PLBYVM
! temporal derivative of the Lateral Boundary fields
REAL, DIMENSION(:,:,:), INTENT(IN) :: PLBXUS,PLBYVS
!
REAL, INTENT(IN) :: PCPHASE ! prescribed phase velocity
REAL, INTENT(IN) :: PCPHASE_PBL ! prescribed PBL phase velocity
!
REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! Jacobian * dry density
! of the reference state
REAL, DIMENSION(:,:,:), INTENT(IN) :: PTKET ! TKE at t
!
REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRUS,PRVS ! Horizontal and Vertical
REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRWS ! momentum tendencies
!
END SUBROUTINE RAD_BOUND
!
END INTERFACE
!
END MODULE MODI_RAD_BOUND
!
! #################################################################
SUBROUTINE RAD_BOUND (HLBCX,HLBCY,HTURB,PRIMKMAX, &
PTSTEP,PDXHAT,PDYHAT,PZHAT, &
PLBXUM,PLBYVM,PLBXUS,PLBYVS, &
PCPHASE,PCPHASE_PBL,PRHODJ, &
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
! #################################################################
!
!!**** *RAD_BOUND* - routine computing the velocity components normal to
!! Lateral boundaries, at the boundary localization
!! for time t+dt.
!!
!! PURPOSE
!! -------
! This computation of wind normal component to the Lateral Boundaries,
! is required to deduce the pressure gradient at the boundaries, used
! as boundary conditions to the elliptic pressure equation to
! be solved in the following step.
!
!!** METHOD
!! ------
!! The computation of the adequate Normal velocity depends on the Lateral
!! Boundary Condition type (variables CLBCX and CLBCY).
!!
!! For the 'OPEN' case, a radiation Sommerfeld equation's type
!! is used:
!! dUn dUn | [ dUn dUn | ]
!! ---- = -----| - C [----- - -----| ]
!! dt dt | LS [ dx dx | LS ]
!!
!!
!! to account for Large Scale (LS) forcings, as proposed by
!! Carpenter (1982), where C is a "magic" phase velocity.
!!
!! A semi-implicit scheme is adopted as in Orlanski (1976).
!!
!! For the 'DAVI' case, the normal velocity is taken as the one of the
!! LS field.
!! Un = Un |
!! | LS
!!
!! EXTERNAL
!! --------
!! GET_INDICE_ll : get physical sub-domain bounds
!! LWEAST_ll,LEAST_ll,LNORTH_ll,LSOUTH_ll : position functions
!!
!! IMPLICIT ARGUMENTS
!! ------------------
!!
!! Module MODD_CONF : contains configuration variable
!! CCONF : Configuration of models
!!
!! Module MODD_PARAMETERS : contains parameters commun to all models
!! JPHEXT : Horizontal EXTernal points number (JPHEXT=1 for this version)
!! JPVEXT : Vertical EXTernal points number (JPVEXT=1 for this version)
!!
!! REFERENCE
!! ---------
!! Book1 and book2 of documentation (routine RAD_BOUND)
!! Orlanski (1976), Carpenter (1982)
!!
!! AUTHOR
!! ------
!! J.-P. Lafore J. Stein * Meteo France *
!!
!! MODIFICATIONS
!! -------------
!! Original 17/10/94
!! Modification 13/02/95 (Lafore) to account for the OPEN case and
!! for the LS fields introduction
!! Modification 16/02/95 (Mallet) bug in computating source terms
!! (depends on START or RESTA configuration)
!! Modification 06/03/95 " change the cphase computation
!! (depends on outflow or inflow)
!! 16/03/95 (J.Stein) remove R from the historical variables
!! Modification 31/05/95 (Lafore) MASTER_DEV2.1 preparation after the
!! LBC tests performed by I. Mallet
!! Modification 03/09/96 (Lafore) computing of LS forcing at central time t
!! previously it was done at t_dt
!! removing R from the LS sources
!! Modification 22/11/96 (Lafore) tests based on LSTEADYLS replaced by the
!! size of LS sources
!! Modification 12/12/96 (Lafore) add a relaxation toward LS fields in the
!! Carpenter equation (semi-implicit treatment)
!! Modification 06/01/97 (Masson) bug at J=2 frontier for relaxation in
!! Carpenter equation
!! Modification 20/10/97 (Lafore, Stein) introduction of the DAVI type of lbc
!! Modification 12/11/97 ( Stein) introduction of LB FIELDS
!! Modification 05/05/04 (Escobar) INOUT argument for PRUS,PRVS,PRWS
!! Modification 05/06 (C.Lac) Remove the DAVI type of lbc
!! Modification 11/09 (C.Lac) Impose a zero prescribed phase velocity
!! in the PBL
!! Juan 25/02/2010: BUG add ZTKEX = 0.0
!! Modification 08/10 (V.Masson) Bug correction and add cphase_profile
!! Escobar 9/11/2010 : cphas_profile : array bound problem if NO Turb => PTKET optional
!! Lac.C. 2011 : Adaptation to FIT temporal scheme
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
!! Modification 06/13 (C.Lac) Introduction of cphase_pbl
!!
!-------------------------------------------------------------------------------
!
!* 0. DECLARATIONS
! ------------
!
USE MODD_CONF
USE MODD_PARAMETERS
USE MODD_CTURB
!
USE MODI_CPHASE_PROFILE
!
USE MODE_ll
!
IMPLICIT NONE
!
!
!* 0.1 declarations of arguments
!
!
!
CHARACTER(LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX,HLBCY ! X and Y-direc. LBC type
CHARACTER(LEN=4), INTENT(IN) :: HTURB ! Turbulence scheme
!
!
REAL, INTENT(IN) :: PRIMKMAX ! Rayleigh damping amplitude
REAL, INTENT(IN) :: PTSTEP ! time step dt
REAL, DIMENSION(:), INTENT(IN) :: PDXHAT ! X-direc. meshlength
REAL, DIMENSION(:), INTENT(IN) :: PDYHAT ! Y-direc. meshlength
REAL, DIMENSION(:), INTENT(IN) :: PZHAT ! height level without orography
!
REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT,PVT ! at t
!
! Lateral Boundary fields at time t
REAL, DIMENSION(:,:,:), INTENT(IN) :: PLBXUM,PLBYVM
! temporal derivative of the Lateral Boundary fields
REAL, DIMENSION(:,:,:), INTENT(IN) :: PLBXUS,PLBYVS
!
REAL, INTENT(IN) :: PCPHASE ! prescribed phase velocity
REAL, INTENT(IN) :: PCPHASE_PBL ! prescribed PBL phase velocity
!
REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! Jacobian * dry density
! of the reference state
REAL, DIMENSION(:,:,:), INTENT(IN) :: PTKET ! TKE at t
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
!
REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRUS,PRVS ! Horizontal and Vertical
REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRWS ! momentum tendancies
!
!* 0.2 declarations of local variables
!
INTEGER :: IIB ! indice I Beginning in x direction
INTEGER :: IJB ! indice J Beginning in y direction
INTEGER :: IIE ! indice I End in x direction
INTEGER :: IJE ! indice J End in y direction
INTEGER :: IKB ! indice K Beginning in z direction
INTEGER :: IKE ! indice K End in z direction
INTEGER :: ILBX,ILBY ! number of points of the RIM arrays
!
REAL :: ZINVTSTEP ! Inverse of the applicable timestep
REAL :: ZKTSTEP ! Rayleigh damping by the timestep
!
REAL, DIMENSION(SIZE(PUT,2),SIZE(PUT,3)) :: ZLBGU ! LB x-gradient
REAL, DIMENSION(SIZE(PUT,2),SIZE(PUT,3)) :: ZLBEU ! LB temporal Evolution
REAL, DIMENSION(SIZE(PUT,2),SIZE(PUT,3)) :: ZLBXU ! LB field at t or t+1
REAL, DIMENSION(SIZE(PUT,2),SIZE(PUT,3)) :: ZCPHASX! Normalized Phase velocity
! ! for U field at X-boundaries
REAL, DIMENSION(SIZE(PUT,2),SIZE(PUT,3)) :: ZPHASX! Phase velocity
! ! for U field at X-boundaries
REAL, DIMENSION(SIZE(PUT,1),SIZE(PUT,3)) :: ZLBGV ! LB y-gradient
REAL, DIMENSION(SIZE(PUT,1),SIZE(PUT,3)) :: ZLBEV ! LB temporal Evolution
REAL, DIMENSION(SIZE(PUT,1),SIZE(PUT,3)) :: ZLBYV ! LB field at t or t+1
REAL, DIMENSION(SIZE(PUT,1),SIZE(PUT,3)) :: ZCPHASY! Normalized Phase velocity
! ! for V field at Y-boundaries
REAL, DIMENSION(SIZE(PUT,1),SIZE(PUT,3)) :: ZPHASY ! Phase velocity
! ! for V field at Y-boundaries
REAL :: ZTSTEP ! effective time step
!
!-------------------------------------------------------------------------------
!
!* 1. PROLOGUE
! --------
!
!* 1.1 Compute dimensions of arrays and other indices
!
CALL GET_INDICE_ll (IIB,IJB,IIE,IJE)
IKB = 1 + JPVEXT
IKE = SIZE(PUT,3) - JPVEXT
!
!
!* 1.2 Compute the inverse of the applicable timestep
!
!
! ZALPHA2 = O : explicit ; ZALPHA2 = 1 : implicit ; ZALPHA2 = 0.5 SI
ZALPHA2 = 1.
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
!
!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
!
!* 2. LBC FILLING IN THE X DIRECTION (LEFT WEST SIDE):
! ------------------------------
! ====> It only concernes U component
! -----------
!
IF (LWEST_ll( )) THEN
!
SELECT CASE ( HLBCX(1) )
!
!* 2.1 WALL CASE:
! ==========
!
CASE ('WALL')
!
PRUS (IIB,:,:) = 0.
!
!* 2.2 OPEN CASE:
! =========
!
CASE ('OPEN')
!
IF (HTURB /= "NONE" ) THEN
CALL CPHASE_PROFILE(PZHAT,PCPHASE,PCPHASE_PBL,ZPHASX,PTKET(IIB,:,:))
ELSE
CALL CPHASE_PROFILE(PZHAT,PCPHASE,PCPHASE_PBL,ZPHASX)
END IF
ZCPHASX(:,:) = MAX ( 0., MIN ( 1., &
(-PUT(IIB,:,:) + ZPHASX(:,:) ) * ZTSTEP / PDXHAT(IIB) ) )
! notice that ZCPHASX=0. when ZPHASX < PUT(IIB,:,:)
!
!
IF ( SIZE(PLBXUS,1) == 0 ) THEN
ZLBEU (:,:) = 0.
ZLBGU (:,:) = PLBXUM(2,:,:) - PLBXUM(1,:,:)
ZLBXU(:,:) = PLBXUM(1,:,:)
ELSE
ZLBEU (:,:) = PLBXUS(1,:,:)
ZLBGU (:,:) = PLBXUM(2,:,:) - PLBXUM(1,:,:) + &
ZTSTEP * (PLBXUS(2,:,:) - PLBXUS(1,:,:))
ZLBXU(:,:) = PLBXUM(1,:,:) + ZTSTEP * PLBXUS(1,:,:)
END IF
!
! ============================================================
!
! PRUS (IIB,:,:) =(PRHODJ(IIB-1,:,:) + PRHODJ(IIB,:,:)) * 0.5 * &
! ( (1. - ZCPHASX(:,:) - ZKTSTEP) * PUM(IIB ,:,:) &
! + 2. * ZCPHASX(:,:) * PUT(IIB+1,:,:) &
! +2.* ( ZLBEU (:,:) * ZTSTEP &
! - ZLBGU (:,:) * ZCPHASX(:,:) &
! + ZKTSTEP*( ZLBXU(:,:) ) ) &
! ) * ZINVTSTEP / (1.+ ZCPHASX(:,:) +ZKTSTEP)
!
PRUS (IIB,:,:) =(PRHODJ(IIB-1,:,:) + PRHODJ(IIB,:,:)) * 0.5 * &
ZINVTSTEP / (1.+ ZKTSTEP * ZALPHA2 ) * &
( (1. - ZCPHASX(:,:) - ZKTSTEP * (1. - ZALPHA2)) * PUT(IIB,:,:) &
+ ZCPHASX(:,:) * PUT(IIB+1 ,:,:) &
+ ( ZLBEU (:,:) * ZTSTEP &
- ZLBGU (:,:) * ZCPHASX(:,:) &
+ ZKTSTEP*ZLBXU(:,:) ) )
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
!
!
END SELECT
!
END IF
!-------------------------------------------------------------------------------
!
!* 3. LBC FILLING IN THE X DIRECTION (RIGHT EAST SIDE):
! ------------------------------
! ====> It only concernes U component
! -----------
!
IF (LEAST_ll( )) THEN
!
!
SELECT CASE ( HLBCX(2) )
!
!* 3.1 WALL CASE:
! =========
!
CASE ('WALL')
!
PRUS (IIE+1,:,:) = 0.
!
!* 3.2 OPEN CASE:
! =========
!
CASE ('OPEN')
!
IF (HTURB /= "NONE" ) THEN
CALL CPHASE_PROFILE(PZHAT,PCPHASE,PCPHASE_PBL,ZPHASX,PTKET(IIE,:,:))
ELSE
CALL CPHASE_PROFILE(PZHAT,PCPHASE,PCPHASE_PBL,ZPHASX)
END IF
!
ZCPHASX(:,:) = MAX ( 0., MIN ( 1., &
( PUT(IIE+1,:,:) + ZPHASX(:,:) ) * ZTSTEP/PDXHAT(IIE) ) )
!
!
ILBX=SIZE(PLBXUM,1)
IF (SIZE(PLBXUS,1) == 0 ) THEN
ZLBEU (:,:) = 0.
ZLBGU (:,:) = PLBXUM(ILBX,:,:) - PLBXUM(ILBX-1,:,:)
ZLBXU(:,:) = PLBXUM(ILBX,:,:)
ELSE
ZLBEU (:,:) = PLBXUS(ILBX,:,:)
ZLBGU (:,:) = PLBXUM(ILBX,:,:) - PLBXUM(ILBX-1,:,:) + &
ZTSTEP * (PLBXUS(ILBX,:,:) - PLBXUS(ILBX-1,:,:))
ZLBXU(:,:) = PLBXUM(ILBX,:,:) + ZTSTEP * PLBXUS(ILBX,:,:)
END IF
!
! ============================================================
!
! PRUS (IIE+1,:,:) =(PRHODJ(IIE+1,:,:) + PRHODJ(IIE,:,:)) * 0.5 * &
! ( (1. - ZCPHASX(:,:) - ZKTSTEP) * PUM(IIE+1,:,:) &
! + 2. * ZCPHASX(:,:) * PUT(IIE ,:,:) &
! +2.* ( ZLBEU (:,:) * ZTSTEP &
! + ZLBGU (:,:) * ZCPHASX(:,:) &
! + ZKTSTEP*ZLBXU(:,:) ) &
! ) * ZINVTSTEP / (1.+ZCPHASX(:,:) +ZKTSTEP)
!
PRUS (IIE+1,:,:) =(PRHODJ(IIE+1,:,:) + PRHODJ(IIE,:,:)) * 0.5 * &
ZINVTSTEP / (1.+ ZKTSTEP * ZALPHA2 ) * &
( (1. - ZCPHASX(:,:) - ZKTSTEP * (1. - ZALPHA2) ) * PUT(IIE+1,:,:) &
+ ZCPHASX(:,:) * PUT(IIE ,:,:) &
+ ( ZLBEU (:,:) * ZTSTEP &
+ ZLBGU (:,:) * ZCPHASX(:,:) &
+ ZKTSTEP*ZLBXU(:,:) ) )
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
!
!
!
END SELECT
!
END IF
!
!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
!
!* 4. LBC FILLING IN THE Y DIRECTION (BOTTOM SOUTH SIDE):
! ------------------------------
! ====> It only concernes V component
! -----------
!
IF (LSOUTH_ll( )) THEN
!
SELECT CASE ( HLBCY(1) )
!
!* 4.1 WALL CASE:
! =========
!
CASE ('WALL')
!
PRVS (:,IJB,:) = 0.
!
!* 4.2 OPEN CASE:
! =========
!
CASE ('OPEN')
!
IF (HTURB /= "NONE" ) THEN
CALL CPHASE_PROFILE(PZHAT,PCPHASE,PCPHASE_PBL,ZPHASY,PTKET(:,IJB,:))
ELSE
CALL CPHASE_PROFILE(PZHAT,PCPHASE,PCPHASE_PBL,ZPHASY)
END IF
!
ZCPHASY(:,:) = MAX ( 0., MIN ( 1., &
(-PVT(:,IJB,:) + ZPHASY(:,:) ) * ZTSTEP/ PDYHAT(IJB) ) )
!
IF ( SIZE(PLBYVS,1) == 0 ) THEN
ZLBEV (:,:) = 0.
ZLBGV (:,:) = PLBYVM(:,2,:) - PLBYVM(:,1,:)
ZLBYV(:,:) = PLBYVM(:,1,:)
ELSE
ZLBEV (:,:) = PLBYVS(:,1,:)
ZLBGV (:,:) = PLBYVM(:,2,:) - PLBYVM(:,1,:) + &
ZTSTEP * (PLBYVS(:,2,:) - PLBYVS(:,1,:))
ZLBYV(:,:) = PLBYVM(:,1,:) + ZTSTEP * PLBYVS(:,1,:)
END IF
!
! ============================================================
!
! PRVS (:,IJB,:) =(PRHODJ(:,IJB-1,:) + PRHODJ(:,IJB,:)) * 0.5 * &
! ( (1. - ZCPHASY(:,:) - ZKTSTEP) * PVM(:,IJB ,:) &
! + 2. * ZCPHASY(:,:) * PVT(:,IJB+1,:) &
! +2.* ( ZLBEV (:,:) * ZTSTEP &
! - ZLBGV (:,:) * ZCPHASY(:,:) &
! + ZKTSTEP*ZLBYV(:,:) ) &
! ) * ZINVTSTEP / (1.+ ZCPHASY(:,:) +ZKTSTEP)
PRVS (:,IJB,:) =(PRHODJ(:,IJB-1,:) + PRHODJ(:,IJB,:)) * 0.5 * &
ZINVTSTEP / (1.+ ZKTSTEP * ZALPHA2 ) * &
( (1. - ZCPHASY(:,:) - ZKTSTEP * (1. - ZALPHA2) ) * PVT(:,IJB,:)&
+ ZCPHASY(:,:) * PVT(:,IJB+1,:) &
+ ( ZLBEV (:,:) * ZTSTEP &
- ZLBGV (:,:) * ZCPHASY(:,:) &
+ ZKTSTEP*ZLBYV(:,:) ) )
!
!
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
!
!
END SELECT
!
END IF
!-------------------------------------------------------------------------------
!
!* 5. LBC FILLING IN THE Y DIRECTION (TOP NORTH SIDE):
! ------------------------------
! ====> It only concernes V component
! -----------
!
IF (LNORTH_ll( )) THEN
!
!
SELECT CASE ( HLBCY(2) )
!
!* 5.1 WALL CASE:
! =========
!
CASE ('WALL')
!
PRVS (:,IJE+1,:) = 0.
!
!* 5.2 OPEN CASE:
! =========
!
CASE ('OPEN')
!
IF (HTURB /= "NONE" ) THEN
CALL CPHASE_PROFILE(PZHAT,PCPHASE,PCPHASE_PBL,ZPHASY,PTKET(:,IJE,:))
ELSE
CALL CPHASE_PROFILE(PZHAT,PCPHASE,PCPHASE_PBL,ZPHASY)
END IF
!
ZCPHASY(:,:) = MAX ( 0., MIN ( 1., &
( PVT(:,IJE+1,:) + ZPHASY(:,:) ) * ZTSTEP/PDYHAT(IJE) ) )
!
ILBY=SIZE(PLBYVM,2)
IF ( SIZE(PLBYVS,1) == 0 ) THEN
ZLBEV (:,:) = 0.
ZLBGV (:,:) = PLBYVM(:,ILBY,:) - PLBYVM(:,ILBY-1,:)
ZLBYV(:,:) = PLBYVM(:,ILBY,:)
ELSE
ZLBEV (:,:) = PLBYVS(:,ILBY,:)
ZLBGV (:,:) = PLBYVM(:,ILBY,:) - PLBYVM(:,ILBY-1,:) + &
ZTSTEP * (PLBYVS(:,ILBY,:) - PLBYVS(:,ILBY-1,:))
ZLBYV(:,:) = PLBYVM(:,ILBY,:) + ZTSTEP * PLBYVS(:,ILBY,:)
END IF
!
! ============================================================
!
! PRVS (:,IJE+1,:) =(PRHODJ(:,IJE+1,:) + PRHODJ(:,IJE,:)) * 0.5 * &
! ( (1. - ZCPHASY(:,:) - ZKTSTEP) * PVM(:,IJE+1,:) &
! + 2. * ZCPHASY(:,:) * PVT(:,IJE ,:) &
! +2.* ( ZLBEV (:,:) * ZTSTEP &
! + ZLBGV (:,:) * ZCPHASY(:,:) &
! + ZKTSTEP* ZLBYV(:,:) ) &
! ) * ZINVTSTEP / (1.+ ZCPHASY(:,:) +ZKTSTEP)
!
PRVS (:,IJE+1,:) =(PRHODJ(:,IJE+1,:) + PRHODJ(:,IJE,:)) * 0.5 * &
ZINVTSTEP / (1.+ ZKTSTEP * ZALPHA2 ) * &
( (1. - ZCPHASY(:,:) - ZKTSTEP * (1. - ZALPHA2) ) * PVT(:,IJE+1,:)&
+ ZCPHASY(:,:) * PVT(:,IJE,:) &
+ ( ZLBEV (:,:) * ZTSTEP &
!
END SELECT
!
END IF
!-------------------------------------------------------------------------------
!
!* 6. UPPER BOUNDARY (FLAT SURFACE):
! ------------------------------
!
PRWS (:,:,IKE+1) = 0.
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE RAD_BOUND