I’m trying to set up a new regional configuration with open boundaries, using release 4.0.5. (Probably switching to 4.2_RC as soon as available.)
After quite some testing there are remaining questions, for which I hope to get some input from the community. In the following, I try to summarize, what I did so far, what works and what doesn’t.
First of all, installing the code and reproducing the test configuration AMM12 works!
The new configuration has three straight open boundaries (N,W,S - coast in the east) that I define using the nambdy_index option with automatic definition of the boundary locations by setting nbdyind = -1 in all three nambdy_index. From ocean.output, I see that these boundaries are correctly identified.
Running with closed boundaries works.
My aim is to prescribe temporally varying boundary conditions for 3D velocities, temperature, salinity and SSH (tides might follow at a later stage - I leave them out for now).
The core question I am struggling with is: What is the correct structure of the boundary conditions to be provided?
I tried several options from boundary data provided on the native model grid to using all the (seemingly) available options of interpolation from the original grid of the boundary data (weights, land/sea mask and vertical interpolation), and everything in between, all failing at some (different) point.
Additionally, I came across two particular problems:
The velocities are full velocities, so I set ln_full_vel = .true.. Doing so leads to segmentation faults. This seems to be connected to undefined dta_alias%u2d and dta_alias%v2d in subroutine bdy_dta: IF( bf_alias(jp_bdyu3d)%ltotvel). I can overcome this by setting cn_dyn2d to something else than none (although none would be the logical choice) and nn_dyn2d_dta = 0 . Is this a bug or am I misunderstanding something here? [I saw similar things being done in older configurations in NEMO3.6]
I’d also like to prescribe SSH but this is not possible with full velocities - correct?
I’ve been using NEMO for quite a while now but am new to BDY and appreciate any input that might help to get me through this process.
I’m glad to start over with this new forum. Hope exchanges will be faster.
boundary data online interpolation: Only vertical interpolation is supported (and provided you have added depths and thicknesses at T/U/V points in the boundary data file). So, it means you have to provide boundary data on your model grid (you can’t use fldread weights functionality).
In the case of straight open boundary segments (basically, your case with ln_coords_file=.false.), arrays in boundary data files should be ordered with along boundary dimension first, and cross boundary dimension second. The latter dimension should be in any case greater or equal to your nn_rimwidth namelist parameter.
cn_dyn2d and ln_full_vel=.true. issue: I don’t understand why you claim none would be the logical choice. cn_dyn2d is the open boundary scheme you choose for barotropic variables (e.g. barotropic transports and sea level anomaly) whatever ln_full_vel is. In practice, we use whether 'frs' (flow relaxation), or 'flather', the latter making use of both internal (model) and external sea level data to prescribe barotropic velocities at the boundaries. The fact that the combination of cn_dyn2d='none' and ln_full_vel=.true. does not work, may, however, be a bug. This situation in real cases is nevertheless unlikely.
Prescribing SSH: this is not possible in NEMO. The only use of external SSH data is when you use the Flather scheme (see item above). The Flather scheme has the good property of maintaining the nested model close to your forcing dataset, hence avoiding possible volume drift. That’s also an excellent gravity wave absorber, all of which making it the ideal scheme for the barotropic mode (the best way to force tides also).
Jerome, I am not sure to see your point about ln_full_vel=.true.. For you, when ln_full_vel=.true. this just means that we must remove the barotropic component and use the one coming from cn_dyn2d?
In my mind ln_full_vel=.true. and as it is coded in bdydta.F90, is working as Franziska is expected: the barotropic component is computed from the 3D array and overwrite u2d which is not read.
For me (say what we had in mind a long time ago when this option ln_full_vel was introduced), ln_full_vel has nothing to do with the BDY scheme you use. It’s just a way to indicate that the barotropic velocity, needed or not, is deduced from the 3d velocity data.
The discussion around cn_dyn2d in combination with the example namelist seems to explain some of the behaviour I observed but didn’t understand so far: Setting ln_full_vel=T, nn_dyn2d_dta = 1 in combination with cn_dyn2d = 'flather' will not try to read [uv]2d but only ssh. That’s exactly what I’d like it to be.
Concerning interpolation: vertical interpolation (including the respective grid information in the bdy-files) seems to work but leads to e.g. negative salinity if ocean values are not extrapolated into the sea floor. Is there a way to mask/flag “dry” grid cells to be excluded from the interpolation?
Are there any attempts to implement the weights functionality of fldread for BDY? This would reduce the currently mandatory pre-processing a lot.
I’ll generate a new set of boundary conditions according to the above described structure and give it another try.
Did you define a _FillValue attribute in your input BDY data and fill land values with it ? This is how the vertical interpolation finds the maximum depth of your input field and eventually extrapolate.
That would be a nice add-on indeed.
You may start by posting an “enhancement” ticket on NEMO website, describing your solution, the NEMO version it applies to, and eventually providing the modified routines as a tar file.
Then, if you belong to one of the consortium members, you can “push” your development into the trunk, an associated development branch being scheduled in the yearly workplan.
Thanks ! I can’t see the “tickets” section, I will read the documentation and create one in the next few days. There is no “enhancement” ticket type, I will use “request” but I can’t see how to attach a tar file…
I realized that this definition was misleading. It does only apply to South and North open boundaries. For West and East boundaries that’s the opposite: cross boundary dimension first, along boundary dimension second. In any cases, the data is ordered in the cross boundary dimension starting from the edge.
I’m still struggling with the structure of the boundary data files, getting warnings pointing to “cluster of rim 1 points” for the northern and southern boundaries. I tried to exclude parts over land, so far without success.
So maybe we can go through it using an example:
Assuming my model domain is rectangular with 200 grid cells in i-direction, 100 in j-direction and 20 depth levels with open ocean in the north, south and west. The coast line in the east is tilted, starting at let’s say i=150 at the northern boundary, ending at i=199 at the southern boundary.
The automatically defined boundary segments will be located at j=2 (south), j=99 (north) and i=2 (west) (corners at 2,2 and 2,99)
Further assuming a rimwidth of 5 and a boundary dataset covering 1 year at daily resolution. How should the dimensions of the individual boundary fields look like?
So far I understood, the order should be (time (unlimited), depth, x, y) (where x is along boundary for south and north and cross boundary for west (and east)). This would result in
south: (365,20,199,5) [or (365,20,197,5) if excluding the land part of the northern boundary]
north: (365,20,199,5) [or (365,20,148,5)]
Is there any mistake? Which order/dimensions would be expected?
A couple of “ncdumps” may worth a thousand words. Here the (working) case of a North East Atlantic domain with 4 open boundaries. nn_rimwidth=30 except for the eastern bdy (15), 1 time dump per file, 50 vertical levels: