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.
Changeset 2915 – NEMO

Changeset 2915


Ignore:
Timestamp:
2011-10-13T17:25:00+02:00 (12 years ago)
Author:
acc
Message:

Branch NOCS_NEPTUNE final style changes and documention. See ticket #843

Location:
branches/2011/dev_r2787_NOCS_NEPTUNE
Files:
1 added
8 edited

Legend:

Unmodified
Added
Removed
  • branches/2011/dev_r2787_NOCS_NEPTUNE/DOC/TexFiles/Biblio/Biblio.bib

    r2541 r2915  
    12751275  url = {http://dx.doi.org/10.1016/j.ocemod.2009.12.003}, 
    12761276  issn = {1463-5003}, 
     1277} 
     1278 
     1279@ARTICLE{HollowayOM86, 
     1280  author = {Greg Holloway}, 
     1281  title = {A Shelf Wave/Topographic Pump Drives Mean Coastal Circulation (part I)}, 
     1282  journal = OM, 
     1283  year = {1986}, 
     1284  volume = {68},   
     1285} 
     1286 
     1287@ARTICLE{HollowayJPO92, 
     1288  author = {Greg Holloway}, 
     1289  title = {Representing Topographic Stress for Large-Scale Ocean Models}, 
     1290  journal = JPO, 
     1291  year = {1992}, 
     1292  volume = {22},   
     1293  pages = {1033--1046}, 
     1294} 
     1295 
     1296@ARTICLE{HollowayJPO94, 
     1297  author = {Michael Eby and Greg Holloway}, 
     1298  title = {Sensitivity of a Large-Scale Ocean Model to a Parameterization of Topographic Stress}, 
     1299  journal = JPO, 
     1300  year = {1994}, 
     1301  volume = {24},   
     1302  pages = {2577--2587}, 
     1303} 
     1304 
     1305@ARTICLE{HollowayJGR09, 
     1306  author = {Greg Holloway and Zeliang Wang}, 
     1307  title = {Representing eddy stress in an Arctic Ocean model}, 
     1308  journal = JGR, 
     1309  year = {2009}, 
     1310  doi = {10.1029/2008JC005169},   
     1311} 
     1312 
     1313@ARTICLE{HollowayOM08, 
     1314  author = {Mathew Maltrud and Greg Holloway}, 
     1315  title = {Implementing biharmonic neptune in a global eddying ocean model}, 
     1316  journal = OM, 
     1317  year = {2008}, 
     1318  volume = {21},   
     1319  pages = {22--34}, 
    12771320} 
    12781321 
  • branches/2011/dev_r2787_NOCS_NEPTUNE/DOC/TexFiles/Chapters/Chap_DYN.tex

    r2541 r2915  
    11621162 
    11631163% ================================================================ 
     1164% Neptune effect  
     1165% ================================================================ 
     1166\section  [Neptune effect (\textit{dynnept})] 
     1167                {Neptune effect (\mdl{dynnept})} 
     1168\label{DYN_nept} 
     1169 
     1170The "Neptune effect" (thus named in \citep{HollowayOM86}) is a 
     1171parameterisation of the potentially large effect of topographic form stress 
     1172(caused by eddies) in driving the ocean circulation. Originally developed for 
     1173low-resolution models, in which it was applied via a Laplacian (second-order) 
     1174diffusion-like term in the momentum equation, it can also be applied in eddy 
     1175permitting or resolving models, in which a more scale-selective bilaplacian 
     1176(fourth-order) implementation is preferred. This mechanism has a 
     1177significant effect on boundary currents (including undercurrents), and the 
     1178upwelling of deep water near continental shelves. 
     1179 
     1180The theoretical basis for the method can be found in  
     1181\citep{HollowayJPO92}, including the explanation of why form stress is not 
     1182necessarily a drag force, but may actually drive the flow.  
     1183\citep{HollowayJPO94} demonstrate the effects of the parameterisation in 
     1184the GFDL-MOM model, at a horizontal resolution of about 1.8 degrees.  
     1185\citep{HollowayOM08} demonstrate the biharmonic version of the 
     1186parameterisation in a global run of the POP model, with an average horizontal 
     1187grid spacing of about 32km. 
     1188 
     1189The NEMO implementation is a simplified form of that supplied by 
     1190Greg Holloway, the testing of which was described in \citep{HollowayJGR09}. 
     1191The major simplification is that a time invariant Neptune velocity 
     1192field is assumed.  This is computed only once, during start-up, and 
     1193made available to the rest of the code via a module.  Vertical 
     1194diffusive terms are also ignored, and the model topography itself 
     1195is used, rather than a separate topographic dataset as in 
     1196\citep{HollowayOM08}.  This implementation is only in the iso-level 
     1197formulation, as is the case anyway for the bilaplacian operator. 
     1198 
     1199The velocity field is derived from a transport stream function given by: 
     1200 
     1201\begin{equation} \label{Eq_dynnept_sf} 
     1202\psi = -fL^2H 
     1203\end{equation} 
     1204 
     1205where $L$ is a latitude-dependant length scale given by: 
     1206 
     1207\begin{equation} \label{Eq_dynnept_ls} 
     1208L = l_1 + (l_2 -l_1)\left ( {1 + \cos 2\phi \over 2 } \right ) 
     1209\end{equation} 
     1210 
     1211where $\phi$ is latitude and $l_1$ and $l_2$ are polar and equatorial length scales respectively. 
     1212Neptune velocity components, $u^*$, $v^*$ are derived from the stremfunction as: 
     1213 
     1214\begin{equation} \label{Eq_dynnept_vel} 
     1215u^* = -{1\over H} {\partial \psi \over \partial y}\ \ \  ,\ \ \ v^* = {1\over H} {\partial \psi \over \partial x} 
     1216\end{equation} 
     1217 
     1218\smallskip 
     1219%----------------------------------------------namdom---------------------------------------------------- 
     1220\namdisplay{namdyn_nept} 
     1221%-------------------------------------------------------------------------------------------------------- 
     1222\smallskip 
     1223 
     1224The Neptune effect is enabled when \np{ln\_neptsimp}=true (default=false). 
     1225\np{ln\_smooth\_neptvel} controls whether a scale-selective smoothing is applied 
     1226to the Neptune effect flow field (default=false) (this smoothing method is as 
     1227used by Holloway).  \np{rn\_tslse} and \np{rn\_tslsp} are the equatorial and 
     1228polar values respectively of the length-scale parameter $L$ used in determining 
     1229the Neptune stream function \eqref{Eq_dynnept_sf} and \eqref{Eq_dynnept_ls}. 
     1230Values at intermediate latitudes are given by a cosine fit, mimicking the 
     1231variation of the deformation radius with latitude.  The default values of 12km 
     1232and 3km are those given in \citep{HollowayJPO94}, appropriate for a coarse 
     1233resolution model. The finer resolution study of \citep{HollowayOM08} increased 
     1234the values of L by a factor of $\sqrt 2$ to 17km and 4.2km, thus doubling the 
     1235stream function for a given topography. 
     1236 
     1237The simple formulation for ($u^*$, $v^*$) can give unacceptably large velocities 
     1238in shallow water, and \citep{HollowayOM08} add an offset to the depth in the 
     1239denominator to control this problem. In this implementation we offer instead (at 
     1240the suggestion of G. Madec) the option of ramping down the Neptune flow field to 
     1241zero over a finite depth range. The switch \np{ln\_neptramp} activates this 
     1242option (default=false), in which case velocities at depths greater than 
     1243\np{rn\_htrmax} are unaltered, but ramp down linearly with depth to zero at a 
     1244depth of \np{rn\_htrmin} (and shallower). 
     1245 
     1246% ================================================================ 
  • branches/2011/dev_r2787_NOCS_NEPTUNE/NEMOGCM/CONFIG/GYRE/EXP00/namelist

    r2825 r2915  
    864864/ 
    865865!----------------------------------------------------------------------- 
    866 &nam_dynnept  !   Neptune effect (simplified: lateral and vertical diffusions removed) 
    867 !----------------------------------------------------------------------- 
    868 ! Suggested lengthscale values are those of Eby & Holloway (1994) for a coarse model 
    869    ln_neptsimp    = .false.  ! yes/no use simplified neptune 
    870  
    871    ln_smoothtopo  = .false.  ! yes/no smooth tsu, tsv 
    872    rn_tslse       =  1.2e4   ! value of lengthscale L at the equator 
    873    rn_tslsp       =  3.0e3   ! value of lengthscale L at the pole 
    874 ! Specify whether to ramp down the Neptune velocity in shallow 
    875 ! water, and if so the depth range controlling such ramping down 
    876    ln_neptramp    = .false.  ! ramp down Neptune velocity in shallow water 
    877    rn_htrmin      =  100.0   ! min. depth of transition range 
    878    rn_htrmax      =  200.0   ! max. depth of transition range 
    879 / 
     866&namdyn_nept  !   Neptune effect (simplified: lateral and vertical diffusions removed) 
     867!----------------------------------------------------------------------- 
     868   ! Suggested lengthscale values are those of Eby & Holloway (1994) for a coarse model 
     869   ln_neptsimp       = .false.  ! yes/no use simplified neptune 
     870 
     871   ln_smooth_neptvel = .false.  ! yes/no smooth zunep, zvnep 
     872   rn_tslse          =  1.2e4   ! value of lengthscale L at the equator 
     873   rn_tslsp          =  3.0e3   ! value of lengthscale L at the pole 
     874   ! Specify whether to ramp down the Neptune velocity in shallow 
     875   ! water, and if so the depth range controlling such ramping down 
     876   ln_neptramp       = .false.  ! ramp down Neptune velocity in shallow water 
     877   rn_htrmin         =  100.0   ! min. depth of transition range 
     878   rn_htrmax         =  200.0   ! max. depth of transition range 
     879/ 
  • branches/2011/dev_r2787_NOCS_NEPTUNE/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist

    r2825 r2915  
    864864/ 
    865865!----------------------------------------------------------------------- 
    866 &nam_dynnept  !   Neptune effect (simplified: lateral and vertical diffusions removed) 
    867 !----------------------------------------------------------------------- 
    868 ! Suggested lengthscale values are those of Eby & Holloway (1994) for a coarse model 
    869    ln_neptsimp    = .false.  ! yes/no use simplified neptune 
    870  
    871    ln_smoothtopo  = .false.  ! yes/no smooth tsu, tsv 
    872    rn_tslse       =  1.2e4   ! value of lengthscale L at the equator 
    873    rn_tslsp       =  3.0e3   ! value of lengthscale L at the pole 
    874 ! Specify whether to ramp down the Neptune velocity in shallow 
    875 ! water, and if so the depth range controlling such ramping down 
    876    ln_neptramp    = .false.  ! ramp down Neptune velocity in shallow water 
    877    rn_htrmin      =  100.0   ! min. depth of transition range 
    878    rn_htrmax      =  200.0   ! max. depth of transition range 
    879 / 
     866&namdyn_nept  !   Neptune effect (simplified: lateral and vertical diffusions removed) 
     867!----------------------------------------------------------------------- 
     868   ! Suggested lengthscale values are those of Eby & Holloway (1994) for a coarse model 
     869   ln_neptsimp       = .true.   ! yes/no use simplified neptune 
     870 
     871   ln_smooth_neptvel = .false.  ! yes/no smooth zunep, zvnep 
     872   rn_tslse          =  1.2e4   ! value of lengthscale L at the equator 
     873   rn_tslsp          =  3.0e3   ! value of lengthscale L at the pole 
     874   ! Specify whether to ramp down the Neptune velocity in shallow 
     875   ! water, and if so the depth range controlling such ramping down 
     876   ln_neptramp       = .true.   ! ramp down Neptune velocity in shallow water 
     877   rn_htrmin         =  100.0   ! min. depth of transition range 
     878   rn_htrmax         =  200.0   ! max. depth of transition range 
     879/ 
  • branches/2011/dev_r2787_NOCS_NEPTUNE/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist

    r2795 r2915  
    879879/ 
    880880!----------------------------------------------------------------------- 
    881 &nam_dynnept  !   Neptune effect (simplified: lateral and vertical diffusions removed) 
    882 !----------------------------------------------------------------------- 
    883    ln_neptsimp          = .false.  ! yes/no use simplified neptune 
    884  
    885    ln_smoothtopo        = .false.  ! yes/no smooth tsu, tsv 
    886    rn_tslse             =  1.2e4   ! value of L at the equator 
    887    rn_tslsp             =  3.0e3   ! value of L at the pole 
    888 / 
     881&namdyn_nept  !   Neptune effect (simplified: lateral and vertical diffusions removed) 
     882!----------------------------------------------------------------------- 
     883   ! Suggested lengthscale values are those of Eby & Holloway (1994) for a coarse model 
     884   ln_neptsimp       = .false.  ! yes/no use simplified neptune 
     885 
     886   ln_smooth_neptvel = .false.  ! yes/no smooth zunep, zvnep 
     887   rn_tslse          =  1.2e4   ! value of lengthscale L at the equator 
     888   rn_tslsp          =  3.0e3   ! value of lengthscale L at the pole 
     889   ! Specify whether to ramp down the Neptune velocity in shallow 
     890   ! water, and if so the depth range controlling such ramping down 
     891   ln_neptramp       = .false.  ! ramp down Neptune velocity in shallow water 
     892   rn_htrmin         =  100.0   ! min. depth of transition range 
     893   rn_htrmax         =  200.0   ! max. depth of transition range 
     894/ 
  • branches/2011/dev_r2787_NOCS_NEPTUNE/NEMOGCM/CONFIG/POMME/EXP00/namelist

    r2825 r2915  
    869869/ 
    870870!----------------------------------------------------------------------- 
    871 &nam_dynnept  !   Neptune effect (simplified: lateral and vertical diffusions removed) 
    872 !----------------------------------------------------------------------- 
    873 ! Suggested lengthscale values are those of Eby & Holloway (1994) for a coarse model 
    874    ln_neptsimp    = .false.  ! yes/no use simplified neptune 
    875  
    876    ln_smoothtopo  = .false.  ! yes/no smooth tsu, tsv 
    877    rn_tslse       =  1.2e4   ! value of lengthscale L at the equator 
    878    rn_tslsp       =  3.0e3   ! value of lengthscale L at the pole 
    879 ! Specify whether to ramp down the Neptune velocity in shallow 
    880 ! water, and if so the depth range controlling such ramping down 
    881    ln_neptramp    = .false.  ! ramp down Neptune velocity in shallow water 
    882    rn_htrmin      =  100.0   ! min. depth of transition range 
    883    rn_htrmax      =  200.0   ! max. depth of transition range 
    884 / 
     871&namdyn_nept  !   Neptune effect (simplified: lateral and vertical diffusions removed) 
     872!----------------------------------------------------------------------- 
     873   ! Suggested lengthscale values are those of Eby & Holloway (1994) for a coarse model 
     874   ln_neptsimp       = .false.  ! yes/no use simplified neptune 
     875 
     876   ln_smooth_neptvel = .false.  ! yes/no smooth zunep, zvnep 
     877   rn_tslse          =  1.2e4   ! value of lengthscale L at the equator 
     878   rn_tslsp          =  3.0e3   ! value of lengthscale L at the pole 
     879   ! Specify whether to ramp down the Neptune velocity in shallow 
     880   ! water, and if so the depth range controlling such ramping down 
     881   ln_neptramp       = .false.  ! ramp down Neptune velocity in shallow water 
     882   rn_htrmin         =  100.0   ! min. depth of transition range 
     883   rn_htrmax         =  200.0   ! max. depth of transition range 
     884/ 
  • branches/2011/dev_r2787_NOCS_NEPTUNE/NEMOGCM/NEMO/OPA_SRC/DYN/dynnept.F90

    r2825 r2915  
    1515   !!                                               to zero added in shallow depths added 
    1616   !!---------------------------------------------------------------------- 
     17   !! dynnept_alloc        : 
    1718   !! dyn_nept_init        : 
    18    !! div_cur_nept_init    : 
    19    !! dyn_cor_topo         : 
    20    !! dyn_topo_neptunevel  : 
    21    !! smooth_topo2         : 
     19   !! dyn_nept_div_cur_init: 
     20   !! dyn_nept_cor         : 
     21   !! dyn_nept_vel         : 
     22   !! dyn_nept_smooth_vel  : 
    2223   !!---------------------------------------------------------------------- 
    2324   USE oce              ! ocean dynamics and tracers 
     
    3536   !! * Routine accessibility 
    3637   PUBLIC dyn_nept_init      ! routine called by nemogcm.F90 
    37    PUBLIC dyn_cor_topo       ! routine called by step.F90 
    38    !! dynnept_alloc()     is called only by dyn_nept_init, within this module 
    39    !! div_cur_nept_init   is called only by dyn_nept_init, within this module 
    40    !! dyn_topo_neptunevel is called only by dyn_cor_topo,  within this module 
     38   PUBLIC dyn_nept_cor       ! routine called by step.F90 
     39   !! dynnept_alloc()         is called only by dyn_nept_init, within this module 
     40   !! dyn_nept_div_cur_init   is called only by dyn_nept_init, within this module 
     41   !! dyn_nept_vel            is called only by dyn_nept_cor,  within this module 
    4142 
    4243   !! * Shared module variables 
     
    4647 
    4748 
    48    !! * Namelist nam_dynnept variables 
    49    LOGICAL, PUBLIC  ::  ln_neptsimp    = .FALSE.  ! yes/no simplified neptune 
    50  
    51    LOGICAL          ::  ln_smoothtopo  = .FALSE.  ! yes/no smooth zunep, zvnep 
    52    REAL(wp)         ::  rn_tslse       =  1.2e4   ! value of lengthscale L at the equator 
    53    REAL(wp)         ::  rn_tslsp       =  3.0e3   ! value of lengthscale L at the pole 
     49   !! * Namelist namdyn_nept variables 
     50   LOGICAL, PUBLIC  ::  ln_neptsimp        = .FALSE.  ! yes/no simplified neptune 
     51 
     52   LOGICAL          ::  ln_smooth_neptvel  = .FALSE.  ! yes/no smooth zunep, zvnep 
     53   REAL(wp)         ::  rn_tslse           =  1.2e4   ! value of lengthscale L at the equator 
     54   REAL(wp)         ::  rn_tslsp           =  3.0e3   ! value of lengthscale L at the pole 
    5455!! Specify whether to ramp down the Neptune velocity in shallow 
    5556!! water, and the depth range controlling such ramping down 
    56    LOGICAL          ::  ln_neptramp    = .FALSE.  ! ramp down Neptune velocity in shallow water 
    57    REAL(wp)         ::  rn_htrmin      =  100.0   ! min. depth of transition range 
    58    REAL(wp)         ::  rn_htrmax      =  200.0   ! max. depth of transition range 
     57   LOGICAL          ::  ln_neptramp        = .FALSE.  ! ramp down Neptune velocity in shallow water 
     58   REAL(wp)         ::  rn_htrmin          =  100.0   ! min. depth of transition range 
     59   REAL(wp)         ::  rn_htrmax          =  200.0   ! max. depth of transition range 
    5960 
    6061   !! * Module variables 
     
    117118      REAL(wp) :: hramp                         ! depth over which Neptune vel. is ramped down 
    118119      !! 
    119       NAMELIST/nam_dynnept/ ln_neptsimp,      & 
    120                             ln_smoothtopo,    & 
     120      NAMELIST/namdyn_nept/ ln_neptsimp,      & 
     121                            ln_smooth_neptvel,& 
    121122             rn_tslse,        & 
    122123             rn_tslsp,        & 
     
    131132!!    WRITE(numout,*) ' start dynnept namelist' 
    132133!!    CALL FLUSH(numout) 
    133       REWIND( numnam )                  ! Read Namelist nam_dynnept:  Simplified Neptune 
    134       READ  ( numnam, nam_dynnept ) 
     134      REWIND( numnam )                  ! Read Namelist namdyn_nept:  Simplified Neptune 
     135      READ  ( numnam, namdyn_nept ) 
    135136!!    WRITE(numout,*) ' dynnept namelist done' 
    136137!!    CALL FLUSH(numout) 
     
    140141         WRITE(numout,*) 'dyn_nept_init : Simplified Neptune module enabled' 
    141142         WRITE(numout,*) '~~~~~~~~~~~~~' 
    142          WRITE(numout,*) ' -->   Reading namelist nam_dynnept parameters:' 
     143         WRITE(numout,*) ' -->   Reading namelist namdyn_nept parameters:' 
    143144         WRITE(numout,*) '       ln_neptsimp          = ', ln_neptsimp 
    144145         WRITE(numout,*) 
    145          WRITE(numout,*) '       ln_smoothtopo        = ', ln_smoothtopo 
     146         WRITE(numout,*) '       ln_smooth_neptvel    = ', ln_smooth_neptvel 
    146147         WRITE(numout,*) '       rn_tslse             = ', rn_tslse 
    147148         WRITE(numout,*) '       rn_tslsp             = ', rn_tslsp 
     
    154155      ENDIF 
    155156 
    156       IF( ln_smoothtopo ) THEN 
     157      IF( ln_smooth_neptvel ) THEN 
    157158         IF(lwp) WRITE(numout,*) ' -->   neptune velocities will be smoothed' 
    158159      ELSE 
     
    213214      END DO 
    214215 
    215       IF( ln_smoothtopo ) THEN 
    216          CALL smooth_topo2( htn, ht, .TRUE. ) 
     216      IF( ln_smooth_neptvel ) THEN 
     217         CALL dyn_nept_smooth_vel( htn, ht, .TRUE. ) 
    217218      !! overwrites ht with a smoothed version of htn 
    218219      ELSE 
     
    224225      !! Compute tsp, a stream function for the Neptune velocity, 
    225226      !! with the usual geophysical sign convention 
    226       !! Then zunep = -latitudinal derivative "-d(tsp)/dy" 
    227       !!      zvnep = longitudinal derivative " d(tsp)/dx" 
     227      !! Then zunep = -latitudinal derivative "-(1/H)*d(tsp)/dy" 
     228      !!      zvnep = longitudinal derivative " (1/H)*d(tsp)/dx" 
    228229 
    229230      tsp(:,:)    = 0.0_wp 
     
    235236 
    236237 
    237       IF( ln_smoothtopo ) THEN 
    238          CALL smooth_topo2( hu, hu_n, .TRUE. ) 
     238      IF( ln_smooth_neptvel ) THEN 
     239         CALL dyn_nept_smooth_vel( hu, hu_n, .TRUE. ) 
    239240      !! overwrites hu_n with a smoothed version of hu 
    240241      ELSE 
     
    252253 
    253254 
    254       IF( ln_smoothtopo ) THEN 
    255          CALL smooth_topo2( hv, hv_n, .TRUE. ) 
     255      IF( ln_smooth_neptvel ) THEN 
     256         CALL dyn_nept_smooth_vel( hv, hv_n, .TRUE. ) 
    256257      !! overwrites hv_n with a smoothed version of hv 
    257258      ELSE 
     
    324325      !!  Compute, once and for all, the horizontal divergence (zhdivnep) 
    325326      !!  and the curl (zmrotnep) of the Neptune velocity field (zunep, zvnep) 
    326       CALL div_cur_nept_init 
     327      CALL dyn_nept_div_cur_init 
    327328 
    328329      !! Check the ranges of the computed divergence & vorticity 
     
    353354 
    354355 
    355    SUBROUTINE div_cur_nept_init 
    356       !!---------------------------------------------------------------------- 
    357       !!             ***  ROUTINE div_cur_nept_init  *** 
     356   SUBROUTINE dyn_nept_div_cur_init 
     357      !!---------------------------------------------------------------------- 
     358      !!             ***  ROUTINE dyn_nept_div_cur_init  *** 
    358359      !! 
    359360      !! ** Purpose :   compute the horizontal divergence and the relative 
     
    394395 
    395396      IF(lwp) WRITE(numout,*) 
    396       IF(lwp) WRITE(numout,*) 'div_cur_nept_init :' 
     397      IF(lwp) WRITE(numout,*) 'dyn_nept_div_cur_init :' 
    397398      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~~~' 
    398399      IF(lwp) WRITE(numout,*) 'horizontal velocity divergence and' 
     
    471472      CALL lbc_lnk( zhdivnep, 'T', 1. )   ;   CALL lbc_lnk( zmrotnep , 'F', 1. )     ! lateral boundary cond. (no sign change) 
    472473      ! 
    473    END SUBROUTINE div_cur_nept_init 
    474  
    475  
    476    SUBROUTINE dyn_cor_topo( kt ) 
    477       !!---------------------------------------------------------------------- 
    478       !!                  ***  ROUTINE dyn_cor_topo  *** 
     474   END SUBROUTINE dyn_nept_div_cur_init 
     475 
     476 
     477   SUBROUTINE dyn_nept_cor( kt ) 
     478      !!---------------------------------------------------------------------- 
     479      !!                  ***  ROUTINE dyn_nept_cor  *** 
    479480      !! 
    480481      !! ** Purpose :  Add or subtract the Neptune velocity from the now velocities 
     
    492493         ! 
    493494         IF( lastkt /= kt ) THEN          ! 1st call for this kt: subtract the Neptune velocities zunep, zvnep from un, vn 
    494             CALL dyn_topo_neptunevel( -1 )      ! -1 = subtract 
     495            CALL dyn_nept_vel( -1 )      ! -1 = subtract 
    495496            ! 
    496497         ELSE                              ! 2nd call for this kt: add the Neptune velocities zunep, zvnep to un, vn 
    497             CALL dyn_topo_neptunevel(  1 )      !  1 = add 
     498            CALL dyn_nept_vel(  1 )      !  1 = add 
    498499            ! 
    499500         ENDIF 
     
    503504      ENDIF 
    504505      ! 
    505    END SUBROUTINE dyn_cor_topo 
    506  
    507  
    508    SUBROUTINE dyn_topo_neptunevel( ksign ) 
    509       !!---------------------------------------------------------------------- 
    510       !!                  ***  ROUTINE dyn_topo_neptunevel  *** 
     506   END SUBROUTINE dyn_nept_cor 
     507 
     508 
     509   SUBROUTINE dyn_nept_vel( ksign ) 
     510      !!---------------------------------------------------------------------- 
     511      !!                  ***  ROUTINE dyn_nept_vel  *** 
    511512      !! 
    512513      !! ** Purpose :  Add or subtract the Neptune velocity from the now 
     
    525526      END DO 
    526527      ! 
    527    END SUBROUTINE dyn_topo_neptunevel 
    528  
    529  
    530    SUBROUTINE smooth_topo2( htold, htnew, option ) 
    531  
    532       !!---------------------------------------------------------------------- 
    533       !!                  ***  ROUTINE smooth_topo2  *** 
     528   END SUBROUTINE dyn_nept_vel 
     529 
     530 
     531   SUBROUTINE dyn_nept_smooth_vel( htold, htnew, option ) 
     532 
     533      !!---------------------------------------------------------------------- 
     534      !!                  ***  ROUTINE dyn_nept_smooth_vel  *** 
    534535      !! 
    535536      !! ** Purpose :  Compute smoothed topography field. 
     
    616617      DEALLOCATE( work, nb, iwork ) 
    617618 
    618    END SUBROUTINE smooth_topo2 
     619   END SUBROUTINE dyn_nept_smooth_vel 
    619620 
    620621END MODULE dynnept 
  • branches/2011/dev_r2787_NOCS_NEPTUNE/NEMOGCM/NEMO/OPA_SRC/step.F90

    r2795 r2915  
    221221      IF(  ln_asmiau .AND. & 
    222222         & ln_dyninc       )   CALL dyn_asm_inc( kstp )     ! apply dynamics assimilation increment 
    223       IF( ln_neptsimp )        CALL dyn_cor_topo( kstp )    ! subtract Neptune velocities (simplified) 
     223      IF( ln_neptsimp )        CALL dyn_nept_cor( kstp )    ! subtract Neptune velocities (simplified) 
    224224                               CALL dyn_adv( kstp )         ! advection (vector or flux form) 
    225225                               CALL dyn_vor( kstp )         ! vorticity term including Coriolis 
    226226                               CALL dyn_ldf( kstp )         ! lateral mixing 
    227       IF( ln_neptsimp )        CALL dyn_cor_topo( kstp )    ! add Neptune velocities (simplified) 
     227      IF( ln_neptsimp )        CALL dyn_nept_cor( kstp )    ! add Neptune velocities (simplified) 
    228228#if defined key_agrif 
    229229      IF(.NOT. Agrif_Root())   CALL Agrif_Sponge_dyn        ! momemtum sponge 
Note: See TracChangeset for help on using the changeset viewer.