River (Point Sources) Input
To configure point sources (typically rivers) or sinks and their associated concentrations of tracers (temperature, salt, inert tracers, ecosystem,, and sediment constituents, etc.), it is necessary to specify the information on the locations of the sources, their transport as a volume flux, how the volume flux is to be distributed over vertical s-coordinate levels, flow direction, and tracer concentrations. The information can be provided analytically or in a forcing NetCDF file.
To provide the information analytically, you need to define the CPP option ANA_PSOURCE and configure the sources in the ana_psource.h Functional. This procedure is described later below. If ANA_PSOURCE is not defined, the default is to look for the information in a source/sink forcing NetCDF file.
Note: In versions of ROMS prior to November 2013 point sources were activated with CPP options UV_PSOURCE (or Q_PSOURCE) and TS_PSOURCE. Those options have been eliminated (see ticket https://www.myroms.org/projects/src/ticket/616) to allow nesting. All grids run the same executable in nesting, so code that might execute in one grid but not all must be activated selectively using logical flags.
The relevant section of roms.in is:
- ! Logical switches (TRUE/FALSE) to activate horizontal momentum transport
! point Sources/Sinks (like river runoff transport) and mass point
! Sources/Sinks (like volume vertical influx), [1:Ngrids].
LuvSrc == T ! horizontal momentum transport
LwSrc == F ! volume vertical influx
The LuvSrc option imposes sources/sinks via the horizontal velocity that crosses a u-face or v-face of a cell. The advection operator then diverges the associated volume and tracer fluxes.
The LwSrc option imposes sources by directly adding to the volume and tracer divergence terms.
Only one of LuvSrc or LwSrc can be true and this selection applies to all sources in the model set-up.
Prior to July 2020 the LwSrc option was coded incorrectly and the majority of ROMS users opted for the LuvSrc approach when specifying river inflows. The LwSrc option is now corrected (see ticket https://www.myroms.org/projects/src/ticket/860) and both options are valid choices.
The principal difference between these options is that LuvSrc adds horizontal momentum to the system at the inflow cell, whereas LwSrc does not.
With LuvSrc it can happen that at very high river discharge the inflow velocity violates CFL stability because the inflow cell area is small (narrow width and/or small height h+ζ). In this situation, LwSrc should still be stable because it does not affect the horizontal momentum equations.
With LuvSrc, the rivers must be placed at cells on the land/sea boundary. With LwSrc, they can in principle be located at any ocean cell.
Note: You may not be able to simply switch between LuvSrc and LwSrc options in roms.in without modifying the river_Xposition and river_Eposition information. This is explained below.
The introduction of tracer fluxes at point sources/sinks is activated with a second set of logical flags.
- ! Logical switches (TRUE/FALSE) to activate tracers point Sources/Sinks
! (like river runoff) and to specify which tracer variables to consider:
! [1:NAT+NPT,Ngrids]. See glossary below for details.
LtracerSrc == F T ! temperature, salinity, inert ...
This example says to impose salinity but not temperature conditions in the river inflow. This requires that variable river_salt is present in the river forcing NetCDF, or is set in ana_psource.h.
- LtracerSrc == F T F T ! temperature, salinity, inert ...
This example says impose the concentration of passive tracer 2, but not passive tracer 1, in all the rivers. (Remember, the first two flags are for temp and salt). In this case, variable river_dye_02 must be in the river forcing NetCDF (or set in ana_psource.h). If you want dye_02 to enter through some rivers but not all, then you set the concentration of river_dye_02 to zero in the appropriate rivers in the forcing file.
The LtracerSrc logical switches replace and eliminate the variable river_flag(river) in the input rivers forcing NetCDF file.
For example, in bio_Fennel.in:
- ! Logical switches (TRUE/FALSE) to activate biological tracers point
! Sources/Sinks (like river runoff) and to specify which tracer variables
! to consider: [NBT,Ngrids] values are expected. See glossary below for
LtracerSrc == T T 10*F
This example says impose concentrations of the first two Fennel model biogeochemical tracers (dissolved nitrate and ammonium) in the rivers but do not impose inflows for any of the other biological model tracers. In this case, variables river_NO3 and river_NH4 must be in the river forcing NetCDF file (or set in ana_psource.h).
The name of the river forcing file is set by its own keyword SSFNAME in roms.in.
- ! Input Sources/Sinks forcing (like river runoff) filename.
SSFNAME == ocean_rivers.nc
Note:This is a change from versions before 11/25/2013 in which the river file was included in the FRCNAME list with other forcing files (meteorology, tides, etc.).
The rivers NetCDF forcing file format is presented in the template frc_rivers.cdl file (in the Data/ROMS/CDL section of the svn distribution).
River NetCDF template
s_rho = 30 ;
river = 4 ;
river_time = UNLIMITED ; // (0 currently)
double river(river) ;
river:long_name = "river runoff identification number" ;
double river_direction(river) ;
river_direction:long_name = "river runoff grid-cell face flag" ;
river_direction:flag_values = "0, 1, 2" ;
river_direction:flag_meanings = "flow across u-face, flow across v-face, flow across w-face" ;
double river_Xposition(river) ;
river_Xposition:long_name = "river XI-position" ;
river_Xposition:LuvSrc_meaning = "i-index grid-cell of u- or v-face source/sink" ;
river_Xposition:LwSrc_meaning = "i-index grid-cell of w-face source/sink" ;
double river_Eposition(river) ;
river_Eposition:long_name = "river ETA-position" ;
river_Xposition:LuvSrc_True_meaning = "j-index grid-cell of u- or v-face source/sink" ;
river_Xposition:LwSrc_True_meaning = "j-index grid-cell of w-face source/sink" ;
double river_Vshape(s_rho, river) ;
river_Vshape:long_name = "river runoff mass transport vertical profile" ;
river_Vshape:requires = "must sum to 1 over s_rho coordinate" ;
double river_time(river_time) ;
river_time:long_name = "river runoff time" ;
river_time:units = "days since 2001-01-01 00:00:00" ;
double river_transport(river_time, river) ;
river_transport:long_name = "river runoff vertically integrated mass transport" ;
river_transport:units = "meter3 second-1" ;
river_transport:positive = "LuvSrc=T flow in positive u,v direction, LwSrc=T flow into grid-cell" ;
river_transport:negative = "LuvSrc=T flow in negative u,v direction, LwSrc=T flow out of grid-cell" ;
river_transport:time = "river_time" ;
double river_temp(river_time, s_rho, river) ;
river_temp:long_name = "river runoff potential temperature" ;
river_temp:units = "Celsius" ;
river_temp:time = "river_time" ;
double river_salt(river_time, s_rho, river) ;
river_salt:long_name = "river runoff salinity" ;
river_salt:time = "river_time" ;
// global attributes:
:rivers = "(1) Connecticut River at Hartford, CT, (2) Hudson River at Green Island NY, (3) Penobscot River at Eddington, ME, (4) Delaware River at Trenton NJ " ;
Creating a river NetCDF file
One approach to generating a rivers NetCDF file is to edit the CDL template above to reflect the number of rivers you have, the number of vertical levels (s_rho), and any other supporting metadata for your convenience. Then the unix command
will create an empty netcdf file you can subsequently fill with data (e.g. inside Matlab). Another tool for creating rivers files is in Python as part of the pyroms package under the examples/rivers directory. Other users can suggest other tools or approaches.
Function and interpretation of river forcing file variables
The NetCDF variable river_direction is a flag indicating the grid-cell face where the point source or sink is specified for a particular application. Its values are read into ROMS variable Dsrc. A value of river_direction=0 corresponds for flow across the u-face, river_direction=1 corresponds for flow across the v-face, and river_direction=2 corresponds for flow across the w-face, as shown in the diagram below.
There are two different methods to impose the point Source/Sink forcing in ROMS:
- As a flux tern to the horizontal advection operator using input switch LuvSrc = T. The point source enters the cell through the u-face or v-face at the specified (i,j,k) location.
- As a volumetric flux term to horizontal mass divergence and horizontal tracer divergence operators using input switch LwSrc = T. The point source enters the cell through the w-face at the specified (i,j,k) location.
Notice that advection and divergence operators center the point Source/Sink effect at the ρ-point (cell center) for all values of Dsrc, as shown in the above figure.
Starting update src:ticket:905 (Jan 2022), both LuvSrc and LwSrc methods are allowed in an application. However, applying both methods to the same grid location is incorrect.
river_Xposition and river_Eposition
When using option LuvSrc = T, river_Xposition and river_Eposition refer to the (i,j) index of the u-face or v-face the flow crosses - NOT the (i,j) index of the ρ-cell it flows into. The (i,j) values must follow ROMS Fortran numbering convention for the appropriate u-point or v-point on the ROMS staggered grid.
This numbering convention is shown in the figure below for flow crossing a u-face into a cell from either the left or the right. This makes it more obvious why the index of the u-face must be specified because to give the (i,j) indices of the receiving ρ-cell would be ambiguous.
The u-face or v-face should be a land/sea mask boundary (i.e., a coastline). If the cell face is placed wholly in the land you get nothing because there is no wet cell for the flow to enter. If the face is in the middle of open water you have a situation where the flow at that cell face computed by the advection algorithm is REPLACED, not augmented, by the source.
When using option LwSrc = T, river_Xposition and river_Eposition refer to the (i,j) index of the rho cell it flows into. The (i,j) values must follow ROMS Fortran numbering convention appropriate to rho-points on the ROMS staggered grid. This is illustrated in the figure below.
Comparing the two figures shows that you may not be able to simply switch between LuvSrc and LwSrc options in roms.in without modifying the river_Xposition and river_Eposition information. If all the inflows are Case (a), and enter the cell from the left then it so happens that the (i,j) index of the receiving rho cell is the same, and river_Xposition does not need to be modified (likewise for v-direction inflows at ). But if the inflow is Case (b) and enters from the right, then the i-index of the rho-cell will be one less and the rivers file needs modification. A change to river_transport is also required.
It is very easy to misconfigure source/sink positions so caution and careful checking is required.
This is very important and very confusing:
The direction of the flow across a u-face or v-face source is controlled by the sign of the variable river_transport, not by the variable river_direction. In the figure above for LuvSrc, if the flow is through the left-hand u-face so as to add water to the ρ-cell immediately to the right, then that flow is in the positive u-direction and river_transport should be positive.
If the flow is through the right-hand u-face so as to add water to the ρ-cell immediately to the left, then that flow is in the negative u-direction and river_transport should be negative.
Likewise for a v-face source.
This is a confusing sign convention, but it exists because in principle one could use the contrary sign for river_transport to create a sink that removes volume from a cell. There is not a lot of experience using this option, but it's there.
For LwSrc, an inflowing source always has positive river_transport.
So if you want to modify a river-forcing file to switch from using LuvSrc to LwSrc, then if you had sources flowing in the negative u- or v-direction, in addition to the changes to river_Xposition or river_Eposition you will need to change river_transport of those sources to be positive.
Variable river_Vshape sets the fractional distribution of the river_transport among the vertical cells and must sum to 1.0 over the vertical. If the sum over k is not 1.0, then the actual volume flux will depart from the value given in river_transport.
Notice that river_Vshape variable can be also used to specify the vertical location of the point Source/Sink. The point Source/Sink can be distributed over the entire water column or be confined to a particular vertical level like surface, bottom, or other levels in between. The only constraint is that the non-dimensional shape function must sum to 1.0 over the vertical range k=1:N for each Source/Sink.
If the river source is in shallow water several grid cells up a river tributary then a simple uniform vertical distribution (all river_Vshape values = 1/s_rho) usually works well. ROMS generates an estuarine circulation to mix the river water with the open ocean. You can manufacture a short tributary in the land/sea mask to do this.
Some users set inflows directly into a deep cell along a coast and use river_Vshape to weigh the inflow toward the sea surface. Use caution if you do this because it can generate 2-grid-point noise at the source and lead to tracer concentration overshoot and undershoot, which can introduce negative salinities.
Configuring point sources/sinks with ana_psource.h
To configure point sources and sinks analytically instead of via a NetCDF file, set the following variables in ana_psource.h
Dsrc(is) ! river_direction for each source is = 1:Nsrc
Isrc(is) ! river_Xposition i-coordinate (u or v point) for each source
Jsrc(is) ! river_Eposition j-coordinate (u or v point) for each source
Qbar(is) ! river_transport for each source following sign conventions noted above
Qshape(is,k) ! river_Vshape distribution of the flow for each source,
! and for every s-coordinate level k = 1:N
The example code in ana_psource.h for the test cases RIVERPLUME2 and SED_TEST1 illustrate how to compute Qshape using the vertical layer thicknesses (by differencing z_w(i,j,k)) to achieve an inflow velocity profile that is uniform with depth.
Tsrc(is,k,isalt) ! salinity values (usually zero for rivers), for each source and all levels.
To include this code for analytical point sources/sinks, you must define ANA_PSOURCE. However, we discourage such practice in realistic applications because of shared-memory or distributed-memory parallel bugs. We recommend always creating the point source/sink NetCDF file. It is easy to check and modify and not subject to parallel bugs.
The logical flags that activate the analytical sources/sink code at runtime (LuvSrc or LwSrc, and LtracerSrc) must still be set in roms.in, sediment.in and/or the appropriate biological model parameter file.
If passive (biology, inert, sediment) tracers are being used, then their concentrations for each river and all vertical levels must also be set in ana_psource.h or in the input NetCDF file.
Here are some forum discussions about the vertical profile of the flow.