source: trunk/NEMO/OPA_SRC/DOM/domvvl.F90 @ 1983

Last change on this file since 1983 was 1983, checked in by rblod, 11 years ago

Intialize fse3 for ln_zco, see ticket #686

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 8.6 KB
Line 
1MODULE domvvl
2   !!======================================================================
3   !!                       ***  MODULE domvvl   ***
4   !! Ocean :
5   !!======================================================================
6   !! History :  2.0  !  2006-06  (B. Levier, L. Marie)  original code
7   !!            3.1  !  2009-02  (G. Madec, M. Leclair, R. Benshila)  pure z* coordinate
8   !!----------------------------------------------------------------------
9#if defined key_vvl
10   !!----------------------------------------------------------------------
11   !!   'key_vvl'                              variable volume
12   !!----------------------------------------------------------------------
13   !!   dom_vvl     : defined coefficients to distribute ssh on each layers
14   !!----------------------------------------------------------------------
15   USE oce             ! ocean dynamics and tracers
16   USE dom_oce         ! ocean space and time domain
17   USE sbc_oce         ! surface boundary condition: ocean
18   USE phycst          ! physical constants
19   USE in_out_manager  ! I/O manager
20   USE lib_mpp         ! distributed memory computing library
21   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC   dom_vvl    ! called by domain.F90
27
28   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)     ::   ee_t, ee_u, ee_v, ee_f   !: ???
29   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   mut, muu, muv, muf       !: ???
30
31   REAL(wp), DIMENSION(jpk) ::   r2dt   ! vertical profile time-step, = 2 rdttra
32      !                                 ! except at nit000 (=rdttra) if neuler=0
33
34   !! * Substitutions
35#  include "domzgr_substitute.h90"
36#  include "vectopt_loop_substitute.h90"
37   !!----------------------------------------------------------------------
38   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)
39   !! $Id$
40   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
41   !!----------------------------------------------------------------------
42
43CONTAINS       
44
45   SUBROUTINE dom_vvl
46      !!----------------------------------------------------------------------
47      !!                ***  ROUTINE dom_vvl  ***
48      !!                   
49      !! ** Purpose :  compute coefficients muX at T-U-V-F points to spread
50      !!               ssh over the whole water column (scale factors)
51      !!----------------------------------------------------------------------
52      INTEGER  ::   ji, jj, jk
53      REAL(wp) ::   zcoefu, zcoefv, zcoeff
54      !!----------------------------------------------------------------------
55
56      IF(lwp)   THEN
57         WRITE(numout,*)
58         WRITE(numout,*) 'dom_vvl : Variable volume activated'
59         WRITE(numout,*) '~~~~~~~~  compute coef. used to spread ssh over each layers'
60      ENDIF
61
62      IF( lk_zco )   CALL ctl_stop( 'dom_vvl : key_zco is incompatible with variable volume option key_vvl')
63
64      IF( ln_zco) THEN
65         DO jk = 1, jpk
66            gdept(:,:,jk) = gdept_0(jk)
67            gdepw(:,:,jk) = gdepw_0(jk)
68            gdep3w(:,:,jk) = gdepw_0(jk)
69            e3t (:,:,jk) = e3t_0(jk)
70            e3u (:,:,jk) = e3t_0(jk)
71            e3v (:,:,jk) = e3t_0(jk)
72            e3f (:,:,jk) = e3t_0(jk)
73            e3w (:,:,jk) = e3w_0(jk)
74            e3uw(:,:,jk) = e3w_0(jk)
75            e3vw(:,:,jk) = e3w_0(jk)
76         END DO
77      ELSE
78         fsdept(:,:,:) = gdept (:,:,:)
79         fsdepw(:,:,:) = gdepw (:,:,:)
80         fsde3w(:,:,:) = gdep3w(:,:,:)
81         fse3t (:,:,:) = e3t   (:,:,:)
82         fse3u (:,:,:) = e3u   (:,:,:)
83         fse3v (:,:,:) = e3v   (:,:,:)
84         fse3f (:,:,:) = e3f   (:,:,:)
85         fse3w (:,:,:) = e3w   (:,:,:)
86         fse3uw(:,:,:) = e3uw  (:,:,:)
87         fse3vw(:,:,:) = e3vw  (:,:,:)
88      ENDIF
89
90      !                                 !==  mu computation  ==!
91      ee_t(:,:) = fse3t_0(:,:,1)                ! Lower bound : thickness of the first model level
92      ee_u(:,:) = fse3u_0(:,:,1)
93      ee_v(:,:) = fse3v_0(:,:,1)
94      ee_f(:,:) = fse3f_0(:,:,1)
95      DO jk = 2, jpkm1                          ! Sum of the masked vertical scale factors
96         ee_t(:,:) = ee_t(:,:) + fse3t_0(:,:,jk) * tmask(:,:,jk)
97         ee_u(:,:) = ee_u(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk)
98         ee_v(:,:) = ee_v(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk)
99         DO jj = 1, jpjm1                      ! f-point : fmask=shlat at coasts, use the product of umask
100            ee_f(:,jj) = ee_f(:,jj) + fse3f_0(:,jj,jk) *  umask(:,jj,jk) * umask(:,jj+1,jk)
101         END DO
102      END DO 
103      !                                         ! Compute and mask the inverse of the local depth at T, U, V and F points
104      ee_t(:,:) = 1. / ee_t(:,:) * tmask(:,:,1)
105      ee_u(:,:) = 1. / ee_u(:,:) * umask(:,:,1)
106      ee_v(:,:) = 1. / ee_v(:,:) * vmask(:,:,1)
107      DO jj = 1, jpjm1                               ! f-point case fmask cannot be used
108         ee_f(:,jj) = 1. / ee_f(:,jj) * umask(:,jj,1) * umask(:,jj+1,1)
109      END DO
110      CALL lbc_lnk( ee_f, 'F', 1. )                  ! lateral boundary condition on ee_f
111      !
112      DO jk = 1, jpk                            ! mu coefficients
113         mut(:,:,jk) = ee_t(:,:) * tmask(:,:,jk)     ! T-point at T levels
114         muu(:,:,jk) = ee_u(:,:) * umask(:,:,jk)     ! U-point at T levels
115         muv(:,:,jk) = ee_v(:,:) * vmask(:,:,jk)     ! V-point at T levels
116      END DO
117      DO jk = 1, jpk                                 ! F-point : fmask=shlat at coasts, use the product of umask
118         DO jj = 1, jpjm1
119               muf(:,jj,jk) = ee_f(:,jj) * umask(:,jj,jk) * umask(:,jj+1,jk)   ! at T levels
120         END DO
121         muf(:,jpj,jk) = 0.e0
122      END DO
123      CALL lbc_lnk( muf, 'F', 1. )                   ! lateral boundary condition
124
125
126      hu_0(:,:) = 0.e0                          ! Reference ocean depth at U- and V-points
127      hv_0(:,:) = 0.e0
128      DO jk = 1, jpk
129         hu_0(:,:) = hu_0(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk)
130         hv_0(:,:) = hv_0(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk)
131      END DO
132
133      DO jj = 1, jpjm1                          ! initialise before and now Sea Surface Height at u-, v-, f-points
134         DO ji = 1, jpim1   ! NO vector opt.
135            zcoefu = 0.5  * umask(ji,jj,1) / ( e1u(ji,jj) * e2u(ji,jj) )
136            zcoefv = 0.5  * vmask(ji,jj,1) / ( e1v(ji,jj) * e2v(ji,jj) )
137            zcoeff = 0.25 * umask(ji,jj,1) * umask(ji,jj+1,1)
138            sshu_b(ji,jj) = zcoefu * ( e1t(ji  ,jj) * e2t(ji  ,jj) * sshb(ji  ,jj)     &
139               &                     + e1t(ji+1,jj) * e2t(ji+1,jj) * sshb(ji+1,jj) )   
140            sshv_b(ji,jj) = zcoefv * ( e1t(ji,jj  ) * e2t(ji,jj  ) * sshb(ji,jj  )     &
141               &                     + e1t(ji,jj+1) * e2t(ji,jj+1) * sshb(ji,jj+1) )   
142            sshf_b(ji,jj) = zcoeff * ( sshb(ji  ,jj) + sshb(ji  ,jj+1)                 &
143               &                     + sshb(ji+1,jj) + sshb(ji+1,jj+1) )               
144            sshu_n(ji,jj) = zcoefu * ( e1t(ji  ,jj) * e2t(ji  ,jj) * sshn(ji  ,jj)     &
145               &                     + e1t(ji+1,jj) * e2t(ji+1,jj) * sshn(ji+1,jj) )   
146            sshv_n(ji,jj) = zcoefv * ( e1t(ji,jj  ) * e2t(ji,jj  ) * sshn(ji,jj  )     & 
147               &                     + e1t(ji,jj+1) * e2t(ji,jj+1) * sshn(ji,jj+1) )     
148            sshf_n(ji,jj) = zcoeff * ( sshn(ji  ,jj) + sshn(ji  ,jj+1)                 &
149               &                     + sshn(ji+1,jj) + sshn(ji+1,jj+1) )               
150         END DO
151      END DO
152      CALL lbc_lnk( sshu_b, 'U', 1. )   ;   CALL lbc_lnk( sshu_n, 'U', 1. )      ! lateral boundary conditions
153      CALL lbc_lnk( sshv_b, 'V', 1. )   ;   CALL lbc_lnk( sshv_n, 'V', 1. )
154      CALL lbc_lnk( sshf_b, 'F', 1. )   ;   CALL lbc_lnk( sshf_n, 'F', 1. )
155      !
156         DO jk = 1, jpkm1
157            fsdept(:,:,jk) = fsdept_n(:,:,jk)          ! now local depths stored in fsdep. arrays
158            fsdepw(:,:,jk) = fsdepw_n(:,:,jk)
159            fsde3w(:,:,jk) = fsde3w_n(:,:,jk)
160            !
161            fse3t (:,:,jk) = fse3t_n (:,:,jk)          ! vertical scale factors stored in fse3. arrays
162            fse3u (:,:,jk) = fse3u_n (:,:,jk)
163            fse3v (:,:,jk) = fse3v_n (:,:,jk)
164            fse3f (:,:,jk) = fse3f_n (:,:,jk)
165            fse3w (:,:,jk) = fse3w_n (:,:,jk)
166            fse3uw(:,:,jk) = fse3uw_n(:,:,jk)
167            fse3vw(:,:,jk) = fse3vw_n(:,:,jk)
168         END DO
169
170
171
172   END SUBROUTINE dom_vvl
173
174#else
175   !!----------------------------------------------------------------------
176   !!   Default option :                                      Empty routine
177   !!----------------------------------------------------------------------
178CONTAINS
179   SUBROUTINE dom_vvl
180   END SUBROUTINE dom_vvl
181#endif
182
183   !!======================================================================
184END MODULE domvvl
Note: See TracBrowser for help on using the repository browser.