NEMO with s-coordinate

Tuesday March 01 2011 06:34:18 pm

Hello,

I am currently attempting to set up a new NEMO configuration using s-levels in the vertical (ln_sco = T) using the latest code v3.3 (out of the SVN repository). I am encountering several problems when using s-coordinates. I was wondering who else is using NEMO with s-coordinates and what code modifications they have made to adapt the code to s-coordinates.


To give you an idea what I'm trying to achieve with my setup, I copy some parts of my namelist below.
!-----------------------------------------------------------------------
&namzgr ! vertical coordinate
!-----------------------------------------------------------------------
ln_zco = .false. ! z-coordinate - full steps (T/F) ("key_zco" may also be defined)
ln_zps = .false. ! z-coordinate - partial steps (T/F)
ln_sco = .true. ! s- or hybrid z-s-coordinate (T/F)
/
!-----------------------------------------------------------------------
&namzgr_sco ! s-coordinate or hybrid z-s-coordinate
!-----------------------------------------------------------------------
rn_sbot_min = 0.001 ! minimum depth of s-bottom surface (>0) (m)
rn_sbot_max = 0.32 ! maximum depth of s-bottom surface (= ocean depth) (>0) (m)
rn_theta = 7.4 ! surface control parameter (0<=rn_theta<=20)
rn_thetb = 1.0 ! bottom control parameter (0<=rn_thetb<= 1)
rn_rmax = 1.0 ! maximum cut-off r-value allowed (0<rn_max<1)
ln_s_sigma = .true. ! hybrid s-sigma coordinates
rn_bb = 1.0 ! stretching with s-sigma
rn_hc = 0.0001 ! critical depth with s-sigma
/

So here are some of my observations/comments/questions:

I want to provide all initial conditions on s-levels, so I can keep the logic of mapping z- to s-levels within the realm of pre-processing. So I deactivated the code which attempts to interpolate the initial conditions (tem_dta/sal_dta) from z-levels to s-levels. The code I mean is in dtatem.f90 where it says "DO jj = 1, jpj ! interpolation of temperatures" and the same again in dtasal.f90 where it says "DO jj = 1, jpj ! interpolation of salinites".

I was wondering why the s-coordinate code refers to the gdept and gdept_0 arrays and why "zgr_z" is still called in domzgr.f90 (line 110)? Is there a deeper reason for this internal reference z-coordinate, other than trying to simplify the pre/post-processing for users who want to deal with input and output files on z-levels? I would like all input to be specified on s-levels and all output to be written on s-levels (my post-processing scripts can easily deal with the inverse mapping from s-levels to z-levels).

This is not a criticism of an internal reference z-coordinate, as I am sure the pre/post-processing effort would be reduced for some users. For me, however, this is not the case and a switch to make this handling optional would help me out in the future. Currently it seems I'll modify the code myself.


Talking of the code... I encountered portions of code which are executed in case of "ln_s_sigma = .true." in domzgr.f90 around line 1394, which has this nested loop:
DO ji = 1, jpi
DO jj = 1, jpj
DO jk = 1, jpk
and within the loop it refers to hbatt(ji+1,jj) and hbatt(ji,jj+1). This is very puzzling as hbatt is allocated to size (jpi,jpj) in dom_oce.f90 (line 162). So, the code will always fail in its current state (the index will always overrun i or j by 1).

Should I change the dimensions of hbatt? Should it have a halo of 1 gridpoint in i and j direction, such that ji+1/jj+1 will not overrun its dimensions?
Or should the above mentioned loop in domzgr.f90 be limited to jpi-i/jpj-1?

Any input from the experts is very much appreciated.
Thank you very much,

Fred Wobus