I am not exactly sure whether this is a real bug, but I would consider the following problem at least as an avoidable inconsistency.

I use ROMS version 854 with the Fennel biology on the northern Gulf of Mexico shelf. The Mississippi and Atchafalaya rivers are implemented as horizontal point sources at various locations on the land-sea interface. At most of the point source locations the freshwater (and tracer) input is distributed over the entire water column. However, at 3 locations the river input is applied to the upper 8 layers only (i.e., k = 13,N(ng)) because of the great bottom depth at these locations. Here, the river input is defined as point sources in Xi direction.

I tried calculating nitrate mass balances for individual grid cells based on the

**DIAGNOSTICS_TS**and

**DIAGNOSTICS_BIO**output (after implementing nitrification as additional

**DIAGNOSTICS_BIO**output). The mass balances work well (in the order of what can be expected from the

**DIAGNOSTICS**output) for all grid cells of the domain . However, in those deeper grid cells (k = 1,12) of the water columns with river input only in the upper layers the horizontal advection across the land-sea interface is non-zero, which from my perspective should simply not be the case.

I was able to trace back this issue to the step3d_uv.F routine, more precisely to the code lines 1188-1207:

Code: Select all

```
!
! Compute correct mass flux, Hz*u/n.
!
DO k=N(ng),1,-1
DO i=IstrP,IendT
Huon(i,j,k)=0.5_r8*(Huon(i,j,k)+u(i,j,k,nnew)*DC(i,k))
# ifdef NEARSHORE_MELLOR
Huon(i,j,k)=Huon(i,j,k)+0.5_r8*u_stokes(i,j,k)*DC(i,k)
# endif
FC(i,0)=FC(i,0)+Huon(i,j,k)
# ifdef DIAGNOSTICS_UV
DiaU3wrk(i,j,k,M3rate)=u(i,j,k,nnew)-DiaU3wrk(i,j,k,M3rate)
# endif
END DO
END DO
DO i=IstrP,IendT
FC(i,0)=DC(i,0)*(FC(i,0)-DU_avg2(i,j)) ! recursive
END DO
DO k=1,N(ng)
DO i=IstrP,IendT
Huon(i,j,k)=Huon(i,j,k)-DC(i,k)*FC(i,0)
END DO
END DO
```

**DC(i,k)*FC(i,0)**in the last loop. Before that,

**Huon**is zero in the layers below the river input cells.

In the following I will consider my case of the river mouth location with

**i = i_river**. To my understanding the following happens:

After the first loop over k and i,

**FC(i_river,0)**equals the overall freshwater discharge by the river.

After the second loop over i (that one labeled as

**! recursive**)

**FC(i_river,0)**equals the difference between the river freshwater input and water column integrated U transport, divided by the water column depth, i.e., it is the water column averaged U transport difference in m2/s.

In the last loop this transport difference is distributed over the entire water column using the individual layer thicknesses as weights, and

**Huon(i,j,k)**is corrected by the corresponding values.

In my case

**FC(i_river,0)-DU_avg2(i_river,j)**is negative which causes the corrected horizontal transport to become positive, i.e., from sea (West) to land (East).

I suppose that the "non-zeroness" of the difference results from the time-stepping of the 2D equations and the incorporation of the river discharge, without knowing where and exactly this mismatch is generated.

I let my simulation run over a period of 4 months in order to check whether this problem only occurs at the very beginning of the simulation. However, this is not the case and the imbalance persists throughout the entire period showing no convergence towards zero. The values range between +/-5.0e-8 mmol N m-3 s-1 which is about +/-4.e4 mmol N d-1 in the considered volume. Indeed, this is only a minor difference and will not change any qualitative results. Nevertheless, advection across the land-sea interface should be generally avoided.

In my opinion, this problem could be solved by applying the volume transport correction at river input locations only to the actual river input grid cells (i.e., by using the rivers Vshape). Though, I am not sure whether this is the best/correct solution or not.

Best,

Fabian

EDIT: My proposed solution may not work as - to my understanding - the proportion of water column depth attributed to the different sigma layers is fixed. That means that a change in water column depth, e.g., due to point source volume discharge, always has to be distributed over the entire water column according to these fixed proportions. Otherwise the model would diverge from these fixed proportions.