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 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.
Regardless of how the source locations and conditions are provided (whether by ANA_PSOURCE or a NetCDF file) the point source code is activated by two sets of logical flags in roms.in.
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 for nesting. In nesting all grids run the same executable, 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 takes the divergence of 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+zeta). 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.
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 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.
If passive tracers are being used, by having set #define T_PASSIVE and the associated options in roms.in, then for those tracers to flow in at a point source the LtracerSrc must also be set to True:
- 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 for variable "river_flag(river)" in the input rivers forcing NetCDF file.
In applications where point sources of biogeochemical or sediment tracers are to be imposed there are corresponding logical flags in e.g. bio_Fennel.in and sediment.in.
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) file name.
SSFNAME == ocean_rivers.nc
The format of the rivers NetCDF forcing file 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_time(river_time) ;
river_time:long_name = "river runoff time" ;
river_time:units = "days since 2001-01-01 00:00:00" ;
double river_direction(river) ;
river_direction:long_name = "river runoff direction" ;
river_direction:flag_values = "0, 1" ;
river_direction:flag_meanings = "flow across u-face, flow across v-face" ;
river_direction:LwSrc_True = "flag not used" ;
double river_Xposition(river) ;
river_Xposition:long_name = "river XI-position" ;
river_Xposition:LuvSrc_meaning = "i point index of U or V face source/sink" ;
river_Xposition:LwSrc_meaning = "i point index of RHO center source/sink" ;
double river_Eposition(river) ;
river_Eposition:long_name = "river ETA-position" ;
river_Xposition:LuvSrc_True_meaning = "j point index of U or V face source/sink" ;
river_Xposition:LwSrc_True_meaning = "j point index of RHO center source/sink" ;
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 RHO-cell" ;
river_transport:negative = "LuvSrc=T flow in negative u,v direction, LwSrc=T flow out of RHO-cell" ;
river_transport:time = "river_time" ;
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" ;
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
When using option LuvSrc = T, the variable river_direction defines whether the river enters a cell through the u-face (river_direction = 0) or v-face (river_direction = 1) of the receiving cell.
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 rho 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 rho-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.
LuvSrc river point source position indexing. Case (a) river flows into rho cell i,j from the u-face to the left, here river_Xposition = i and river_Eposition = j. Case (b) river flows into rho cell i,j from the u-face to the right, here river_Xposition = i+1 and river_Eposition = j.
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.
LwSrc option river point source position indexing. Case (c) river flows into rho cell i,j directly. Here river_Xposition = i and river_Eposition = j. Contrast with the figure of LuvSrc to see why for case (b) sources it is necessary to subtract 1 from river_Xposition if changing between LwSrc and LuvSrc options.
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 ot be modified (likewise for v-direction inflows at v_i,j). 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 rho 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 rho 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 over the vertical. If the sum over k is not 1 then the actual volume flux will depart from the value given in river_transport.
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 cells along a coast and use river_Vshape to weight 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)) so as 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.
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 tracers, or tracers for sediment or biology are being used, then the concentrations for those tracers for each river and all vertical levels must also be set in ana_psource.h.
Here are some forum discussions about the vertical profile of the flow.