Scientific Investigations Report 2007–5185
U.S. GEOLOGICAL SURVEY
Scientific Investigations Report 2007–5185
USGS personnel obtained the Willamette River temperature TMDL models and model source code from ODEQ staff. ODEQ used CE‑QUAL‑W2 version 3.12 from August 19, 2003, for each of the submodels, with one modification by the PSU model development team. That “Tmax” modification enabled the model to internally calculate and output the daily maximum water temperature for every segment in the model in three different ways: (a) surface-water temperature only, (b) volume-weighted water temperature for the entire cross section, and (c) flow-weighted water temperature for the entire cross section. ODEQ used this “Tmax” version for almost all of its Willamette TMDL modeling. The Tmax version was not used on the North Santiam and Santiam Rivers; the version used there did not have the Tmax code additions, but the model computations were otherwise identical.
When USGS staff needed to run the Willamette River models, several code changes were necessary to make the models easier to use or to eliminate problems. Great restraint was exercised in keeping code changes to a minimum, so as to ensure that the model calibration would remain unchanged. All code changes were checked extensively and are described and justified in this appendix. In the end, USGS used its “generic Tmax3” version for all Willamette River temperature models except the lower Willamette River, where the “generic Tmax3a” version was applied. The details behind these versions, and the reasons for their use, are described below.
The models as provided by ODEQ to USGS contained a “visual interface” such that when they were run under Microsoft© Windows, a new window opened to display selected model inputs and results during the model run. When the model run was completed, however, the display window did not close itself and did not return control to the operating system. Therefore, the model version that contains this visual interface was not amenable to unattended operation in which batch files are used to run the submodels in sequence, passing the output of one model as input to the next model downstream without active human intervention.
To enable unattended linkages between submodels, USGS staff removed the visual interface from the source code obtained from ODEQ. The code modifications were designed to create a “generic” version of CE‑QUAL‑W2. When CE‑QUAL‑W2 is distributed to the public by the PSU model development team, two versions typically are provided: the “cvf” version with the visual interface, and a generic version without the visual interface. This generic version of the model code created by USGS was different from the version with the visual interface in exactly the same way that the PSU cvf and generic versions differ.
This code change has no effect on the model’s computations and therefore has absolutely no effect on the model results. This was checked by running the cvf and generic versions of the Willamette models and comparing their output. No significant differences were observed. The minor differences that did exist were easily attributable to the use of slightly different compilation options.
The version of CE‑QUAL‑W2 used in the Willamette River temperature models includes a clever bit of coding by the PSU model development team that helps sloping river models to continue running under a wide range of flow conditions. In each of the model’s “water bodies,” the model keeps track of a “surface-layer index.” The water surface, though, can actually reside in a model layer that is above the surface-layer index. To enable the model to continue running when a model segment still has flow, but its bottom is above the surface layer index, the PSU team added some code that tells the model that the river bottom is lower in that segment than originally specified. This “trick” works fairly well, but was only partially implemented by the PSU model development team at the time the Willamette temperature models were being developed.
In running one of the Willamette River models under slightly different conditions, USGS staff observed that the model could crash because this “trick” was not sufficiently robust. The trick allows the surface layer index to drop below the bottom-most active layer in a model segment, assigning a small nonzero width to the next lower layer to keep the water flowing in a sloping reach. This originally was not allowed to occur if the segment had the largest Z value in the water body – if the water surface in that segment was the lowest in the water body, relative to the layer boundaries and irrespective of slope.
The code was changed so that all segments in a water body could have their bottoms “lowered,” as long as they are not lowered past the bottom-most active layer in that water body. Without this fix, the original trick is insufficient to keep the river from “drying up” under certain conditions, and the model could crash unnecessarily. This modification was tested originally by the USGS on another river system (the Tualatin River), and it was determined to work as intended. The code changes were shared with the PSU model development team at that time and have since been incorporated into the latest public release version of CE‑QUAL‑W2. So, in the interest of keeping the model from crashing unnecessarily, this change was made to the code. In the model code, these modifications are denoted with the comments “!SR 08/01/05” and “!SR 08/08/05,” which are the dates of modification for when these fixes were made to USGS versions of the code. This change does not affect the computations of the model. It simply keeps it from crashing.
As stated previously, great care was taken to minimize the number of code changes that were made to the Willamette models, in order to preserve the model calibration and ensure that the models were essentially the same as those used to develop the TMDL. In a few instances, though, it made sense to correct a couple of small coding errors.
First, in the LATERAL_WITHDRAWAL subroutine, several fixes were applied. The statement:
IF (KBOT > KB(ID)) KBOT = KB(ID)
incorrectly referenced the most-downstream segment in the branch (ID), when it should have referenced the segment where the withdrawal is located (I). The line was changed to
IF (KBOT > KB(I)) KBOT = KB(I)
The original statement could result in the improper allocation of a withdrawal flow rate as a function of model layer. This error was fixed by the PSU model development team in later versions, but was never changed in the Willamette models. Next, the following code was added:
IF (QWD(JWD) == 0.0) THEN
KTW(JWD) = KTOP
KBW(JWD) = KBOT
RETURN
END IF
This helps to avoid problems with the calculation of KTOP and KBOT, and streamlines the code slightly. This might have been the source of problems in which the original code would hang if the withdrawal flow rate was zero. In any case, the change explicitly enables the model to handle withdrawal flow rates of zero.
Elsewhere in the LATERAL_WITHDRAWAL subroutine, the distribution of the withdrawal flow rate among the relevant model layers utilized values of BHR() and BHRKT2() in the original code. These were carry-overs from the DOWNSTREAM_WITHDRAWAL routine and were not appropriate for a lateral withdrawal. These variables were changed to BH() and BHKT2(). These fixes were cleared with the PSU development team in 2004. In addition, two further slight code modifications were added to help avoid potential divide-by-zero errors in this routine. Those fixes are labeled with the comment “!SR 09/16/04,” which was the date when these errors were originally fixed in other versions of the code.
In the DOWNSTREAM_WITHDRAWAL subroutine, two snippets of code were added in an effort to avoid divide-by-zero errors, mirroring code that was added to the LATERAL_WITHDRAWAL routine. Look for comments with the label “!SR 09/16/04.”
Finally, a problem was fixed in the setting of the BHKT1(IU-1) variable when segments are being subtracted after a layer is subtracted. There was an extra set of parentheses in that computation that would cause the water depth in that cell to be computed incorrectly. Look for the comment “!SR 07/29/05.” This change is inconsequential, because segments are never subtracted in the Willamette temperature model applications. Still, this was an obvious coding error that has since been corrected by the PSU model development team.
All changes up to this point were incorporated into the USGS model version named “generic Tmax3,” or “Tmax3” for short. The effect of these changes was tested, and no significant changes to the model output resulted. Most of the coding changes were minor, would rarely be invoked, or were “defensive” in nature, guarding against divide-by-zero errors, for example, and therefore would not normally be expected to cause a change in the model results. Mainly, these bug fixes serve to keep the model running under extraordinary circumstances, and to keep its computations more correct. Many more bug fixes and upgrades could be applied to the Willamette temperature models (the public release version of CE‑QUAL‑W2 is at version 3.5 as of April of 2007, and the Willamette models were built on version 3.12), but a major upgrade to the model code would require a re-evaluation of the model calibration, which is beyond the scope and need of this investigation. The Tmax3 version was used by USGS to run all Willamette submodels except the lower Willamette River, as explained below.
In version Tmax3 and in the original CE‑QUAL‑W2 version obtained from ODEQ, the calculations that determine the daily maximum surface temperature, the daily maximum volume-weighted temperature, and the daily maximum flow-weighted temperature for every segment only calculated those quantities when the “time series” output was printed. So, if the user specified that the model should write time-series output 20 times per day (time series output frequency (TSRF) = 0.05), then the modeled temperatures would be queried only 20 times per day to find the daily maxima. At times, this may not be frequent enough. It particularly can be a problem when the model time steps are large, such that different model runs might not calculate daily maxima using information from the exact same times. If the user-specified maximum time step is, say, 360 seconds, as it was in the lower Willamette River model, then the typical time step might be several minutes, and different model runs (no point sources versus all point sources, for example) might end up using information from several minutes apart in determining the daily maxima, simply due to variations in the model’s variable time steps. As a result, two different model runs might show daily maximum temperature differences of several hundredths of a degree or more, when such differences are only an artifact of the number of times per day that the model checked to determine whether a daily maximum had occurred. This “timing artifact” is not real and can be minimized or eliminated in one of several ways.
Simply increasing the number of times per day that the model checks for a daily maximum helps to minimize this timing artifact. As a test, the lower Willamette River model was run with a time-series output frequency of 200 points per day (TSRF = 0.005) and the results were compared to those from a model run using a time-series output frequency of only 20 points per day. The magnitude of the artifact was decreased, but not eliminated. Similarly, tests were run using a user-specified maximum time step (DLTMAX) of 60 seconds rather than 360 seconds for the lower Willamette River model, and this also decreased the magnitude of the timing artifact. In fact, decreasing the maximum time step to 60 seconds was more effective than checking more frequently.
In an attempt to minimize this artifact to the greatest extent possible, two things were changed for the lower Willamette River model. First, the maximum time step was decreased from 360 seconds to 60 seconds. That did not entail a code change. Second, a change was made to the model code so that the model would determine the daily maximum water temperatures (surface, volume-weighted, flow-weighted) using information from each and every time step, rather than only a certain number of times per day. This change means that the daily maximum temperature calculation is no longer tied to the TSRF. This new version of the model was named “generic Tmax3a” and was used only for the lower Willamette River submodel, because that model had the longest time steps and was the most susceptible to this timing artifact.
Tests showed that the new Tmax3a code had smaller timing artifacts in the computed daily maximum water temperatures. Separate tests showed that the Tmax3a version did not produce results that were any different from the Tmax3 version for the upper Willamette model, where the typical time step is much smaller.
U.S. Department of the Interior | U.S. Geological Survey
URL: https://pubs.usgs.gov/sir/2007/5185
Page Contact Information: Publications Team
Page Last Modified: Thursday, 01-Dec-2016 19:50:25 EST
