New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
domvvl.F90 in branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/DOM – NEMO

source: branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/DOM/domvvl.F90 @ 2236

Last change on this file since 2236 was 2236, checked in by cetlod, 14 years ago

First guess of NEMO_v3.3
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 8.8 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  (NEMOGCM/License_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                   ! temporary scalars
54      REAL(wp) ::   zv_t_ij, zv_t_ip1j, zv_t_ijp1, zv_t_ip1jp1   !     -        -
55      REAL(wp), DIMENSION(jpi,jpj) ::  zs_t, zs_u_1, zs_v_1      !     -     2D workspace
56      !!----------------------------------------------------------------------
57
58      IF(lwp)   THEN
59         WRITE(numout,*)
60         WRITE(numout,*) 'dom_vvl : Variable volume activated'
61         WRITE(numout,*) '~~~~~~~~  compute coef. used to spread ssh over each layers'
62      ENDIF
63
64      IF( lk_zco )   CALL ctl_stop( 'dom_vvl : key_zco is incompatible with variable volume option key_vvl')
65
66      fsdept(:,:,:) = gdept (:,:,:)
67      fsdepw(:,:,:) = gdepw (:,:,:)
68      fsde3w(:,:,:) = gdep3w(:,:,:)
69      fse3t (:,:,:) = e3t   (:,:,:)
70      fse3u (:,:,:) = e3u   (:,:,:)
71      fse3v (:,:,:) = e3v   (:,:,:)
72      fse3f (:,:,:) = e3f   (:,:,:)
73      fse3w (:,:,:) = e3w   (:,:,:)
74      fse3uw(:,:,:) = e3uw  (:,:,:)
75      fse3vw(:,:,:) = e3vw  (:,:,:)
76
77      !                                 !==  mu computation  ==!
78      ee_t(:,:) = fse3t_0(:,:,1)                ! Lower bound : thickness of the first model level
79      ee_u(:,:) = fse3u_0(:,:,1)
80      ee_v(:,:) = fse3v_0(:,:,1)
81      ee_f(:,:) = fse3f_0(:,:,1)
82      DO jk = 2, jpkm1                          ! Sum of the masked vertical scale factors
83         ee_t(:,:) = ee_t(:,:) + fse3t_0(:,:,jk) * tmask(:,:,jk)
84         ee_u(:,:) = ee_u(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk)
85         ee_v(:,:) = ee_v(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk)
86         DO jj = 1, jpjm1                      ! f-point : fmask=shlat at coasts, use the product of umask
87            ee_f(:,jj) = ee_f(:,jj) + fse3f_0(:,jj,jk) *  umask(:,jj,jk) * umask(:,jj+1,jk)
88         END DO
89      END DO 
90      !                                         ! Compute and mask the inverse of the local depth at T, U, V and F points
91      ee_t(:,:) = 1. / ee_t(:,:) * tmask(:,:,1)
92      ee_u(:,:) = 1. / ee_u(:,:) * umask(:,:,1)
93      ee_v(:,:) = 1. / ee_v(:,:) * vmask(:,:,1)
94      DO jj = 1, jpjm1                               ! f-point case fmask cannot be used
95         ee_f(:,jj) = 1. / ee_f(:,jj) * umask(:,jj,1) * umask(:,jj+1,1)
96      END DO
97      CALL lbc_lnk( ee_f, 'F', 1. )                  ! lateral boundary condition on ee_f
98      !
99      DO jk = 1, jpk                            ! mu coefficients
100         mut(:,:,jk) = ee_t(:,:) * tmask(:,:,jk)     ! T-point at T levels
101         muu(:,:,jk) = ee_u(:,:) * umask(:,:,jk)     ! U-point at T levels
102         muv(:,:,jk) = ee_v(:,:) * vmask(:,:,jk)     ! V-point at T levels
103      END DO
104      DO jk = 1, jpk                                 ! F-point : fmask=shlat at coasts, use the product of umask
105         DO jj = 1, jpjm1
106               muf(:,jj,jk) = ee_f(:,jj) * umask(:,jj,jk) * umask(:,jj+1,jk)   ! at T levels
107         END DO
108         muf(:,jpj,jk) = 0.e0
109      END DO
110      CALL lbc_lnk( muf, 'F', 1. )                   ! lateral boundary condition
111
112
113      hu_0(:,:) = 0.e0                          ! Reference ocean depth at U- and V-points
114      hv_0(:,:) = 0.e0
115      DO jk = 1, jpk
116         hu_0(:,:) = hu_0(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk)
117         hv_0(:,:) = hv_0(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk)
118      END DO
119     
120      ! surface at t-points and inverse surface at (u/v)-points used in surface averaging computations
121      ! for ssh and scale factors
122      zs_t  (:,:) =       e1t(:,:) * e2t(:,:)
123      zs_u_1(:,:) = 0.5 / e1u(:,:) * e2u(:,:)
124      zs_v_1(:,:) = 0.5 / e1v(:,:) * e2v(:,:)
125
126      DO jj = 1, jpjm1                          ! initialise before and now Sea Surface Height at u-, v-, f-points
127         DO ji = 1, jpim1   ! NO vector opt.
128            zcoefu = umask(ji,jj,1) * zs_u_1(ji,jj)
129            zcoefv = vmask(ji,jj,1) * zs_v_1(ji,jj)
130            zcoeff = 0.5 * umask(ji,jj,1) * umask(ji,jj+1,1) / ( e1f(ji,jj) * e2f(ji,jj) )
131            ! before fields
132            zv_t_ij       = zs_t(ji  ,jj  ) * sshb(ji  ,jj  )
133            zv_t_ip1j     = zs_t(ji+1,jj  ) * sshb(ji+1,jj  )
134            zv_t_ijp1     = zs_t(ji  ,jj+1) * sshb(ji  ,jj+1)
135            sshu_b(ji,jj) = zcoefu * ( zv_t_ij + zv_t_ip1j )
136            sshv_b(ji,jj) = zcoefv * ( zv_t_ij + zv_t_ijp1 )
137            ! now fields
138            zv_t_ij       = zs_t(ji  ,jj  ) * sshn(ji  ,jj  )
139            zv_t_ip1j     = zs_t(ji+1,jj  ) * sshn(ji+1,jj  )
140            zv_t_ijp1     = zs_t(ji  ,jj+1) * sshn(ji  ,jj+1)
141            zv_t_ip1jp1   = zs_t(ji  ,jj+1) * sshn(ji  ,jj+1)
142            sshu_n(ji,jj) = zcoefu * ( zv_t_ij + zv_t_ip1j )
143            sshv_n(ji,jj) = zcoefv * ( zv_t_ij + zv_t_ijp1 )
144            sshf_n(ji,jj) = zcoeff * ( zv_t_ij + zv_t_ip1j + zv_t_ijp1 + zv_t_ip1jp1 )
145         END DO
146      END DO
147      CALL lbc_lnk( sshu_n, 'U', 1. )   ;   CALL lbc_lnk( sshu_b, 'U', 1. )      ! lateral boundary conditions
148      CALL lbc_lnk( sshv_n, 'V', 1. )   ;   CALL lbc_lnk( sshv_b, 'V', 1. )
149      CALL lbc_lnk( sshf_n, 'F', 1. )
150
151                                                ! initialise before scale factors at (u/v)-points
152      ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points
153      DO jk = 1, jpkm1
154         DO jj = 1, jpjm1
155            DO ji = 1, jpim1
156               zv_t_ij           = zs_t(ji  ,jj  ) * fse3t_b(ji  ,jj  ,jk)
157               zv_t_ip1j         = zs_t(ji+1,jj  ) * fse3t_b(ji+1,jj  ,jk)
158               zv_t_ijp1         = zs_t(ji  ,jj+1) * fse3t_b(ji  ,jj+1,jk)
159               fse3u_b(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) )
160               fse3v_b(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) )
161            END DO
162         END DO
163      END DO
164      CALL lbc_lnk( fse3u_b(:,:,:), 'U', 1. )               ! lateral boundary conditions
165      CALL lbc_lnk( fse3v_b(:,:,:), 'V', 1. )
166      ! Add initial scale factor to scale factor anomaly
167      fse3u_b(:,:,:) = fse3u_b(:,:,:) + fse3u_0(:,:,:)
168      fse3v_b(:,:,:) = fse3v_b(:,:,:) + fse3v_0(:,:,:)
169      !
170   END SUBROUTINE dom_vvl
171
172#else
173   !!----------------------------------------------------------------------
174   !!   Default option :                                      Empty routine
175   !!----------------------------------------------------------------------
176CONTAINS
177   SUBROUTINE dom_vvl
178   END SUBROUTINE dom_vvl
179#endif
180
181   !!======================================================================
182END MODULE domvvl
Note: See TracBrowser for help on using the repository browser.