Nothing
R/Manningrect.R
In iemisc: Irucka Embry's Miscellaneous Functions
#' Rectangular cross-section for the Gauckler-Manning-Strickler equation
#'
#' This function solves for one missing variable in the Gauckler-Manning-
#' Strickler equation for a rectangular cross-section and uniform flow. The
#' \code{\link[stats]{uniroot}} function is used to obtain the missing parameters.
#'
#'
#'
#'
#' Gauckler-Manning-Strickler equation is expressed as
#'
#' \deqn{V = \frac{K_n}{n}R^\frac{2}{3}S^\frac{1}{2}}
#'
#' \describe{
#' \item{\emph{V}}{the velocity (m/s or ft/s)}
#' \item{\emph{n}}{Manning's roughness coefficient (dimensionless)}
#' \item{\emph{R}}{the hydraulic radius (m or ft)}
#' \item{\emph{S}}{the slope of the channel bed (m/m or ft/ft)}
#' \item{\emph{K_n}}{the conversion constant -- 1.0 for SI and
#' 3.2808399 ^ (1 / 3) for English units -- m^(1/3)/s or ft^(1/3)/s}
#' }
#'
#'
#'
#'
#' This equation is also expressed as
#'
#' \deqn{Q = \frac{K_n}{n}\frac{A^\frac{5}{3}}{P^\frac{2}{3}}S^\frac{1}{2}}
#'
#' \describe{
#' \item{\emph{Q}}{the discharge {m^3/s or ft^3/s (cfs)} is VA}
#' \item{\emph{n}}{Manning's roughness coefficient (dimensionless)}
#' \item{\emph{P}}{the wetted perimeters of the channel (m or ft)}
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{S}}{the slope of the channel bed (m/m or ft/ft)}
#' \item{\emph{K_n}}{the conversion constant -- 1.0 for SI and
#' 3.2808399 ^ (1 / 3) for English units -- m^(1/3)/s or ft^(1/3)/s}
#' }
#'
#'
#'
#'
#' Other important equations regarding the rectangular cross-section follow:
#' \deqn{R = \frac{A}{P}}
#'
#' \describe{
#' \item{\emph{R}}{the hydraulic radius (m or ft)}
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{P}}{the wetted perimeters of the channel (m or ft)}
#' }
#'
#'
#'
#'
#' \deqn{A = by}
#'
#' \describe{
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{y}}{the flow depth (normal depth in this function) [m or ft]}
#' \item{\emph{b}}{the bottom width (m or ft)}
#' }
#'
#'
#'
#'
#' \deqn{P = b + 2y}
#'
#' \describe{
#' \item{\emph{P}}{the wetted perimeters of the channel (m or ft)}
#' \item{\emph{y}}{the flow depth (normal depth in this function) [m or ft]}
#' \item{\emph{b}}{the bottom width (m or ft)}
#' }
#'
#'
#'
#'
#' \deqn{B = b}
#'
#' \describe{
#' \item{\emph{B}}{the top width of the channel (m or ft)}
#' \item{\emph{b}}{the bottom width (m or ft)}
#' }
#'
#'
#'
#' \deqn{D = \frac{A}{B}}
#'
#' \describe{
#' \item{\emph{D}}{the hydraulic depth (m or ft)}
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{B}}{the top width of the channel (m or ft)}
#' }
#'
#'
#'
#' \deqn{Z = \frac{\sqrt{2}}{2}my^2.5}
#'
#' \describe{
#' \item{\emph{Z}}{the Section factor (m or ft)}
#' \item{\emph{y}}{the flow depth (normal depth in this function) [m or ft]}
#' \item{\emph{m}}{the horizontal side slope}
#' }
#'
#'
#'
#' \deqn{E = y + \frac{Q^2}{2gA^2}}
#'
#' \describe{
#' \item{\emph{E}}{the Specific Energy (m or ft)}
#' \item{\emph{Q}}{the discharge {m^3/s or ft^3/s (cfs)} is VA}
#' \item{\emph{g}}{gravitational acceleration (m/s^2 or ft/sec^2)}
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{y}}{the flow depth (normal depth in this function) [m or ft]}
#' }
#'
#'
#'
#'
#' \deqn{VH = \frac{V^2}{2g}}
#'
#' \describe{
#' \item{\emph{VH}}{the Velocity Head (m or ft)}
#' \item{\emph{V}}{the velocity (m/s or ft/s)}
#' \item{\emph{g}}{gravitational acceleration (m/s^2 or ft/sec^2)}
#' }
#'
#'
#'
#' A rough turbulent zone check is performed on the water flowing in the
#' channel using the Reynolds number (Re). The Re equation follows:
#'
#' \deqn{Re = \frac{\\rho RV}{\\mu}}
#'
#' \describe{
#' \item{\emph{Re}}{Reynolds number (dimensionless)}
#' \item{\emph{\\rho}}{density (kg/m^3 or slug/ft^3)}
#' \item{\emph{R}}{the hydraulic radius (m or ft)}
#' \item{\emph{V}}{the velocity (m/s or ft/s)}
#' \item{\emph{\\mu}}{dynamic viscosity (* 10^-3 kg/m*s or * 10^-5 lb*s/ft^2)}
#' }
#'
#'
#'
#' A critical flow check is performed on the water flowing in the channel
#' using the Froude number (Fr). The Fr equation follows:
#'
#' \deqn{Fr = \frac{V}{\left(\sqrt{g * D}\right)}}
#'
#' \describe{
#' \item{\emph{Fr}}{the Froude number (dimensionless)}
#' \item{\emph{V}}{the velocity (m/s or ft/s)}
#' \item{\emph{g}}{gravitational acceleration (m/s^2 or ft/sec^2)}
#' \item{\emph{D}}{the hydraulic depth (m or ft)}
#' }
#'
#'
#'
#' @note
#' Assumptions: uniform flow, prismatic channel, and surface water temperature
#' of 20 degrees Celsius (68 degrees Fahrenheit) at atmospheric pressure
#'
#' Note: Units must be consistent
#'
#' Please refer to the iemisc: Manning... Examples using iemiscdata
#' [https://www.ecoccs.com/R_Examples/Manning_iemiscdata_Examples.pdf] and iemisc:
#' Open Channel Flow Examples involving Geometric Shapes with the
#' Gauckler-Manning-Strickler Equation
#' [https://www.ecoccs.com/R_Examples/Open-Channel-Flow_Examples_Geometric_Shapes.pdf]
#' for the cross-section examples using iemiscdata
#'
#'
#'
#'
#' @param Q numeric vector that contains the discharge value {m^3/s or ft^3/s},
#' if known.
#' @param n numeric vector that contains the Manning's roughness coefficient n,
#' if known.
#' @param b numeric vector that contains the bottom width, if known.
#' @param Sf numeric vector that contains the bed slope (m/m or ft/ft),
#' if known.
#' @param y numeric vector that contains the flow depth (m or ft), if known.
#' @param Temp numeric vector that contains the temperature (degrees C or degrees
#' Fahrenheit), if known.
#' @param units character vector that contains the system of units {options are
#' \code{SI} for International System of Units or \code{Eng} for English units
#' (United States Customary System in the United States and Imperial Units in
#' the United Kingdom)}
#'
#' @return the missing parameters (Q, n, b, Sf, or y) & area (A), wetted
#' perimeter (P), velocity (V), top width (B), hydraulic radius (R),
#' Reynolds number (Re), and Froude number (Fr) as a \code{\link[base]{list}}.
#'
#'
#' @references
#' \enumerate{
#' \item Terry W. Sturm, \emph{Open Channel Hydraulics}, 2nd Edition, New York City, New York: The McGraw-Hill Companies, Inc., 2010, page 2, 8, 36, 102, 120, 153-154.
#' \item Dan Moore, P.E., NRCS Water Quality and Quantity Technology Development Team, Portland Oregon, "Using Mannings Equation with Natural Streams", August 2011, \url{https://web.archive.org/web/20210416091858/https://www.wcc.nrcs.usda.gov/ftpref/wntsc/H&H/xsec/manningsNaturally.pdf}. Retrieved thanks to the Internet Archive: Wayback Machine
#' \item Gilberto E. Urroz, Utah State University Civil and Environmental Engineering - OCW, CEE6510 - Numerical Methods in Civil Engineering, Spring 2006 (2006). Course 3. "Solving selected equations and systems of equations in hydraulics using Matlab", August/September 2004, \url{https://digitalcommons.usu.edu/ocw_cee/3/}.
#' \item Tyler G. Hicks, P.E., \emph{Civil Engineering Formulas: Pocket Guide}, 2nd Edition, New York City, New York: The McGraw-Hill Companies, Inc., 2002, page 423, 425.
#' \item Wikimedia Foundation, Inc. Wikipedia, 26 November 2015, "Manning formula", \url{https://en.wikipedia.org/wiki/Manning_formula}.
#' \item John C. Crittenden, R. Rhodes Trussell, David W. Hand, Kerry J. Howe, George Tchobanoglous, \emph{MWH's Water Treatment: Principles and Design}, Third Edition, Hoboken, New Jersey: John Wiley & Sons, Inc., 2012, page 1861-1862.
#' \item Andrew Chadwick, John Morfett and Martin Borthwick, \emph{Hydraulics in Civil and Environmental Engineering}, Fourth Edition, New York City, New York: Spon Press, Inc., 2004, page 133.
#' \item Robert L. Mott and Joseph A. Untener, \emph{Applied Fluid Mechanics}, Seventh Edition, New York City, New York: Pearson, 2015, page 376, 379-380.
#' \item Ven Te Chow, Ph.D., \emph{Open-Channel Hydraulics}, McGraw-Hill Classic Textbook Reissue, New York City, New York: McGraw-Hill Book Company, 1988, pages 21, 40-41.
#' \item Gary P. Merkley, "BIE6300 - Irrigation & Conveyance Control Systems, Spring 2004", 2004, Biological and Irrigation Engineering - OCW. Course 2, \url{https://digitalcommons.usu.edu/ocw_bie/2/}.
#' \item The NIST Reference on Constants, Units, and Uncertainty, Fundamental Constants Data Center of the NIST Physical Measurement Laboratory, "standard acceleration of gravity g_n", \url{https://physics.nist.gov/cgi-bin/cuu/Value?gn}.
#' \item Wikimedia Foundation, Inc. Wikipedia, 15 May 2019, "Conversion of units", \url{https://en.wikipedia.org/wiki/Conversion_of_units}.
#' }
#'
#'
#'
#' @author Irucka Embry
#'
#'
#'
#' @encoding UTF-8
#'
#'
#'
#'
#'
#'
#'
#'
#'
#' @seealso \code{\link{Manningtrap}} for a trapezoidal cross-section, \code{\link{Manningtri}}
#' for a triangular cross-section, \code{\link{Manningpara}} for a parabolic
#' cross-section, and \code{\link{Manningcirc}} for a circular cross-section.
#'
#'
#'
#'
#'
#'
#'
#'
#' @importFrom fpCompare %==%
#' @importFrom assertthat assert_that
#' @importFrom checkmate qtest
#' @importFrom units set_units make_units drop_units
#' @importFrom round round_r3
#' @importFrom stats uniroot
#'
#' @export
Manningrect <- function (Q = NULL, n = NULL, b = NULL, Sf = NULL, y = NULL, Temp = NULL, units = c("SI", "Eng")) {
K <- NULL
# due to NSE notes in R CMD check
checks <- c(Q, n, b, Sf, y)
units <- units
# Check
assert_that(!any(qtest(checks, "N+(0,)") == FALSE), msg = "Either Q, n, b, Sf, or y is 0, NA, NaN, Inf, -Inf, empty, or a string. Please try again.")
# only process with finite values and provide an error message if the check fails
assert_that(qtest(units, "S==1"), msg = "There is not an unit type or more than 1 unit type. Please specify either 'SI' or 'Eng'.")
# only process with enough known variables and provide an error message if the check fails
assert_that(isTRUE(any(c("SI", "Eng") %in% units)), msg = "The unit system has not been identified correctly as either 'SI' or 'Eng'. Please try again.")
# only process with a specified unit and provide a stop warning if not
# units
if (units == "SI") {
# use the temperature to determine the density & absolute and kinematic viscosities
Temp <- ifelse(is.null(Temp), 20, Temp) # degrees C
assert_that(qtest(Temp, "N+(0,)"), msg = "Either Temp is equal to or less than 0 C, NA, NaN, Inf, -Inf, empty, or a string. The equation is valid only for temperatures greater than 0 C / 32 F. Please try again.")
# only process with specified, finite values and provide an error message if the check fails
# create a numeric vector with the units of degrees Celsius
T_C <- set_units(Temp, "degree_C")
# create a numeric vector to convert from degrees Celsius to Kelvin
T_K <- T_C
# create a numeric vector with the units of Kelvin
units(T_K) <- make_units(K)
# create viscosities based on temperatures
# saturated liquid density at given temperature in degrees Celsius (SI units)
rho_SI <- density_water(Temp, "SI")
# absolute or dynamic viscosity at given temperature in degrees Celsius and density of rho (SI units)
mu_SI <- dyn_visc_water(Temp, "SI")
# kinematic viscosity at given temperature in degrees Celsius and density of rho (SI units)
nu_SI <- kin_visc_water(rho_SI, mu_SI, rho_units = "kg/m^3", mu_units = "Pa*s or kg/m/s")
k <- 1
g <- 9.80665 # m / s^2
rho <- rho_SI
nu <- nu_SI
mu <- mu_SI
# unit weight of water at given temperature in degrees Celsius and density of rho (SI units)
gamma <- unit_wt(rho = rho, units = "SI")
density_water_units <- "kg/m^3"
dyn_visc_water_units <- "Pa * s or kg/m*s"
kin_visc_water_units <- "m^2/s"
result_units <- c("m", "m^2", "m", "m", "m", "m", "m", "m/s", "m^3/s", "dimensionless", "m/m", "degrees Celsius", "Kelvin", density_water_units, dyn_visc_water_units, kin_visc_water_units, "dimensionless", "dimensionless", "m/m", "m/m", "m/m", "m", "m", "m", "m", "m^3/s", "m", "m", "pascal (N/m^2)", "pascal (N/m^2)")
} else if (units == "Eng") {
# use the temperature to determine the density & absolute and kinematic viscosities
Temp <- ifelse(is.null(Temp), 68, Temp) # degrees F
assert_that(qtest(Temp, "N+(32,)"), msg = "Either Temp is equal to or less than 32 F, NA, NaN, Inf, -Inf, empty, or a string. The equation is valid only for temperatures greater than 32 F / 0 C. Please try again.")
# only process with specified, finite values and provide an error message if the check fails
T_F <- Temp
# create a numeric vector with the units of degrees Fahrenheit
T_F <- set_units(T_F, "degree_F")
# create a numeric vector to convert from degrees Fahrenheit to Kelvin
T_K <- T_F
# create a numeric vector with the units of Kelvin
units(T_K) <- make_units(K)
# create viscosities based on temperatures
# saturated liquid density at given temperature in degrees Fahrenheit (US Customary units)
rho_Eng <- density_water(Temp, "Eng", Eng_units = "slug/ft^3")
# absolute or dynamic viscosity at given temperature in degrees Fahrenheit and density of rho (US Customary units)
mu_Eng <- dyn_visc_water(Temp, "Eng", Eng_units = "slug/ft/s")
# kinematic viscosity at given temperature in degrees Fahrenheit and density of rho (US Customary units)
nu_Eng <- kin_visc_water(rho_Eng, mu_Eng, rho_units = "slug/ft^3", mu_units = "slug/ft/s")
k <- 3.2808399 ^ (1 / 3)
g <- 9.80665 * (3937 / 1200) # ft / sec^2
gc <- 9.80665 * (3937 / 1200) # lbm-ft/lbf-sec^2
rho <- rho_Eng
nu <- nu_Eng
mu <- mu_Eng
# unit weight of water at given temperature in degrees Fahrenheit and density of rho (US Customary units)
gamma <- unit_wt(rho = rho, units = "Eng", Eng_units = "slug/ft^3")
density_water_units <- "slug/ft^3"
dyn_visc_water_units <- "slug/ft*s"
kin_visc_water_units <- "ft^2/s"
result_units <- c("ft", "ft^2", "ft", "ft", "ft", "ft", "ft", "ft/sec (fps)", "ft^3/sec (cfs)", "dimensionless", "ft/ft", "degrees Fahrenheit", "Kelvin", density_water_units, dyn_visc_water_units, kin_visc_water_units, "dimensionless", "dimensionless", "ft/ft", "ft/ft", "ft/ft", "ft", "ft", "ft", "ft", "ft^3/sec (cfs)", "ft", "ft", "lb/ft^2", "lb/ft^2")
}
if (missing(Q)) {
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
Qfun <- function(Q) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
Quse <- uniroot(Qfun, interval = c(0.0000001, 200), extendInt = "yes")
Q <- Quse$root
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(Q = Q, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
} else if (missing(n)) {
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
nfun <- function(n) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
nuse <- uniroot(nfun, interval = c(0.0000001, 200), extendInt = "yes")
n <- nuse$root
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(n = n, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
} else if (missing(b)) {
bfun <- function(b) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
buse <- uniroot(bfun, interval = c(0.0000001, 200), extendInt = "yes")
b <- buse$root
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(b = b, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
} else if (missing(y)) {
yfun <- function(y) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
yuse <- uniroot(yfun, interval = c(0.0000001, 200), extendInt = "yes")
y <- yuse$root
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(y = y, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
} else if (missing(Sf)) {
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
Sffun <- function(Sf) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
Sfuse <- uniroot(Sffun, interval = c(0.0000001, 200), extendInt = "yes")
Sf <- Sfuse$root
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(Sf = Sf, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
}
}
Try the iemisc package in your browser
Any scripts or data that you put into this service are public.
iemisc documentation built on Sept. 25, 2023, 5:09 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
#' Rectangular cross-section for the Gauckler-Manning-Strickler equation
#'
#' This function solves for one missing variable in the Gauckler-Manning-
#' Strickler equation for a rectangular cross-section and uniform flow. The
#' \code{\link[stats]{uniroot}} function is used to obtain the missing parameters.
#'
#'
#'
#'
#' Gauckler-Manning-Strickler equation is expressed as
#'
#' \deqn{V = \frac{K_n}{n}R^\frac{2}{3}S^\frac{1}{2}}
#'
#' \describe{
#' \item{\emph{V}}{the velocity (m/s or ft/s)}
#' \item{\emph{n}}{Manning's roughness coefficient (dimensionless)}
#' \item{\emph{R}}{the hydraulic radius (m or ft)}
#' \item{\emph{S}}{the slope of the channel bed (m/m or ft/ft)}
#' \item{\emph{K_n}}{the conversion constant -- 1.0 for SI and
#' 3.2808399 ^ (1 / 3) for English units -- m^(1/3)/s or ft^(1/3)/s}
#' }
#'
#'
#'
#'
#' This equation is also expressed as
#'
#' \deqn{Q = \frac{K_n}{n}\frac{A^\frac{5}{3}}{P^\frac{2}{3}}S^\frac{1}{2}}
#'
#' \describe{
#' \item{\emph{Q}}{the discharge {m^3/s or ft^3/s (cfs)} is VA}
#' \item{\emph{n}}{Manning's roughness coefficient (dimensionless)}
#' \item{\emph{P}}{the wetted perimeters of the channel (m or ft)}
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{S}}{the slope of the channel bed (m/m or ft/ft)}
#' \item{\emph{K_n}}{the conversion constant -- 1.0 for SI and
#' 3.2808399 ^ (1 / 3) for English units -- m^(1/3)/s or ft^(1/3)/s}
#' }
#'
#'
#'
#'
#' Other important equations regarding the rectangular cross-section follow:
#' \deqn{R = \frac{A}{P}}
#'
#' \describe{
#' \item{\emph{R}}{the hydraulic radius (m or ft)}
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{P}}{the wetted perimeters of the channel (m or ft)}
#' }
#'
#'
#'
#'
#' \deqn{A = by}
#'
#' \describe{
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{y}}{the flow depth (normal depth in this function) [m or ft]}
#' \item{\emph{b}}{the bottom width (m or ft)}
#' }
#'
#'
#'
#'
#' \deqn{P = b + 2y}
#'
#' \describe{
#' \item{\emph{P}}{the wetted perimeters of the channel (m or ft)}
#' \item{\emph{y}}{the flow depth (normal depth in this function) [m or ft]}
#' \item{\emph{b}}{the bottom width (m or ft)}
#' }
#'
#'
#'
#'
#' \deqn{B = b}
#'
#' \describe{
#' \item{\emph{B}}{the top width of the channel (m or ft)}
#' \item{\emph{b}}{the bottom width (m or ft)}
#' }
#'
#'
#'
#' \deqn{D = \frac{A}{B}}
#'
#' \describe{
#' \item{\emph{D}}{the hydraulic depth (m or ft)}
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{B}}{the top width of the channel (m or ft)}
#' }
#'
#'
#'
#' \deqn{Z = \frac{\sqrt{2}}{2}my^2.5}
#'
#' \describe{
#' \item{\emph{Z}}{the Section factor (m or ft)}
#' \item{\emph{y}}{the flow depth (normal depth in this function) [m or ft]}
#' \item{\emph{m}}{the horizontal side slope}
#' }
#'
#'
#'
#' \deqn{E = y + \frac{Q^2}{2gA^2}}
#'
#' \describe{
#' \item{\emph{E}}{the Specific Energy (m or ft)}
#' \item{\emph{Q}}{the discharge {m^3/s or ft^3/s (cfs)} is VA}
#' \item{\emph{g}}{gravitational acceleration (m/s^2 or ft/sec^2)}
#' \item{\emph{A}}{the cross-sectional area (m^2 or ft^2)}
#' \item{\emph{y}}{the flow depth (normal depth in this function) [m or ft]}
#' }
#'
#'
#'
#'
#' \deqn{VH = \frac{V^2}{2g}}
#'
#' \describe{
#' \item{\emph{VH}}{the Velocity Head (m or ft)}
#' \item{\emph{V}}{the velocity (m/s or ft/s)}
#' \item{\emph{g}}{gravitational acceleration (m/s^2 or ft/sec^2)}
#' }
#'
#'
#'
#' A rough turbulent zone check is performed on the water flowing in the
#' channel using the Reynolds number (Re). The Re equation follows:
#'
#' \deqn{Re = \frac{\\rho RV}{\\mu}}
#'
#' \describe{
#' \item{\emph{Re}}{Reynolds number (dimensionless)}
#' \item{\emph{\\rho}}{density (kg/m^3 or slug/ft^3)}
#' \item{\emph{R}}{the hydraulic radius (m or ft)}
#' \item{\emph{V}}{the velocity (m/s or ft/s)}
#' \item{\emph{\\mu}}{dynamic viscosity (* 10^-3 kg/m*s or * 10^-5 lb*s/ft^2)}
#' }
#'
#'
#'
#' A critical flow check is performed on the water flowing in the channel
#' using the Froude number (Fr). The Fr equation follows:
#'
#' \deqn{Fr = \frac{V}{\left(\sqrt{g * D}\right)}}
#'
#' \describe{
#' \item{\emph{Fr}}{the Froude number (dimensionless)}
#' \item{\emph{V}}{the velocity (m/s or ft/s)}
#' \item{\emph{g}}{gravitational acceleration (m/s^2 or ft/sec^2)}
#' \item{\emph{D}}{the hydraulic depth (m or ft)}
#' }
#'
#'
#'
#' @note
#' Assumptions: uniform flow, prismatic channel, and surface water temperature
#' of 20 degrees Celsius (68 degrees Fahrenheit) at atmospheric pressure
#'
#' Note: Units must be consistent
#'
#' Please refer to the iemisc: Manning... Examples using iemiscdata
#' [https://www.ecoccs.com/R_Examples/Manning_iemiscdata_Examples.pdf] and iemisc:
#' Open Channel Flow Examples involving Geometric Shapes with the
#' Gauckler-Manning-Strickler Equation
#' [https://www.ecoccs.com/R_Examples/Open-Channel-Flow_Examples_Geometric_Shapes.pdf]
#' for the cross-section examples using iemiscdata
#'
#'
#'
#'
#' @param Q numeric vector that contains the discharge value {m^3/s or ft^3/s},
#' if known.
#' @param n numeric vector that contains the Manning's roughness coefficient n,
#' if known.
#' @param b numeric vector that contains the bottom width, if known.
#' @param Sf numeric vector that contains the bed slope (m/m or ft/ft),
#' if known.
#' @param y numeric vector that contains the flow depth (m or ft), if known.
#' @param Temp numeric vector that contains the temperature (degrees C or degrees
#' Fahrenheit), if known.
#' @param units character vector that contains the system of units {options are
#' \code{SI} for International System of Units or \code{Eng} for English units
#' (United States Customary System in the United States and Imperial Units in
#' the United Kingdom)}
#'
#' @return the missing parameters (Q, n, b, Sf, or y) & area (A), wetted
#' perimeter (P), velocity (V), top width (B), hydraulic radius (R),
#' Reynolds number (Re), and Froude number (Fr) as a \code{\link[base]{list}}.
#'
#'
#' @references
#' \enumerate{
#' \item Terry W. Sturm, \emph{Open Channel Hydraulics}, 2nd Edition, New York City, New York: The McGraw-Hill Companies, Inc., 2010, page 2, 8, 36, 102, 120, 153-154.
#' \item Dan Moore, P.E., NRCS Water Quality and Quantity Technology Development Team, Portland Oregon, "Using Mannings Equation with Natural Streams", August 2011, \url{https://web.archive.org/web/20210416091858/https://www.wcc.nrcs.usda.gov/ftpref/wntsc/H&H/xsec/manningsNaturally.pdf}. Retrieved thanks to the Internet Archive: Wayback Machine
#' \item Gilberto E. Urroz, Utah State University Civil and Environmental Engineering - OCW, CEE6510 - Numerical Methods in Civil Engineering, Spring 2006 (2006). Course 3. "Solving selected equations and systems of equations in hydraulics using Matlab", August/September 2004, \url{https://digitalcommons.usu.edu/ocw_cee/3/}.
#' \item Tyler G. Hicks, P.E., \emph{Civil Engineering Formulas: Pocket Guide}, 2nd Edition, New York City, New York: The McGraw-Hill Companies, Inc., 2002, page 423, 425.
#' \item Wikimedia Foundation, Inc. Wikipedia, 26 November 2015, "Manning formula", \url{https://en.wikipedia.org/wiki/Manning_formula}.
#' \item John C. Crittenden, R. Rhodes Trussell, David W. Hand, Kerry J. Howe, George Tchobanoglous, \emph{MWH's Water Treatment: Principles and Design}, Third Edition, Hoboken, New Jersey: John Wiley & Sons, Inc., 2012, page 1861-1862.
#' \item Andrew Chadwick, John Morfett and Martin Borthwick, \emph{Hydraulics in Civil and Environmental Engineering}, Fourth Edition, New York City, New York: Spon Press, Inc., 2004, page 133.
#' \item Robert L. Mott and Joseph A. Untener, \emph{Applied Fluid Mechanics}, Seventh Edition, New York City, New York: Pearson, 2015, page 376, 379-380.
#' \item Ven Te Chow, Ph.D., \emph{Open-Channel Hydraulics}, McGraw-Hill Classic Textbook Reissue, New York City, New York: McGraw-Hill Book Company, 1988, pages 21, 40-41.
#' \item Gary P. Merkley, "BIE6300 - Irrigation & Conveyance Control Systems, Spring 2004", 2004, Biological and Irrigation Engineering - OCW. Course 2, \url{https://digitalcommons.usu.edu/ocw_bie/2/}.
#' \item The NIST Reference on Constants, Units, and Uncertainty, Fundamental Constants Data Center of the NIST Physical Measurement Laboratory, "standard acceleration of gravity g_n", \url{https://physics.nist.gov/cgi-bin/cuu/Value?gn}.
#' \item Wikimedia Foundation, Inc. Wikipedia, 15 May 2019, "Conversion of units", \url{https://en.wikipedia.org/wiki/Conversion_of_units}.
#' }
#'
#'
#'
#' @author Irucka Embry
#'
#'
#'
#' @encoding UTF-8
#'
#'
#'
#'
#'
#'
#'
#'
#'
#' @seealso \code{\link{Manningtrap}} for a trapezoidal cross-section, \code{\link{Manningtri}}
#' for a triangular cross-section, \code{\link{Manningpara}} for a parabolic
#' cross-section, and \code{\link{Manningcirc}} for a circular cross-section.
#'
#'
#'
#'
#'
#'
#'
#'
#' @importFrom fpCompare %==%
#' @importFrom assertthat assert_that
#' @importFrom checkmate qtest
#' @importFrom units set_units make_units drop_units
#' @importFrom round round_r3
#' @importFrom stats uniroot
#'
#' @export
Manningrect <- function (Q = NULL, n = NULL, b = NULL, Sf = NULL, y = NULL, Temp = NULL, units = c("SI", "Eng")) {
K <- NULL
# due to NSE notes in R CMD check
checks <- c(Q, n, b, Sf, y)
units <- units
# Check
assert_that(!any(qtest(checks, "N+(0,)") == FALSE), msg = "Either Q, n, b, Sf, or y is 0, NA, NaN, Inf, -Inf, empty, or a string. Please try again.")
# only process with finite values and provide an error message if the check fails
assert_that(qtest(units, "S==1"), msg = "There is not an unit type or more than 1 unit type. Please specify either 'SI' or 'Eng'.")
# only process with enough known variables and provide an error message if the check fails
assert_that(isTRUE(any(c("SI", "Eng") %in% units)), msg = "The unit system has not been identified correctly as either 'SI' or 'Eng'. Please try again.")
# only process with a specified unit and provide a stop warning if not
# units
if (units == "SI") {
# use the temperature to determine the density & absolute and kinematic viscosities
Temp <- ifelse(is.null(Temp), 20, Temp) # degrees C
assert_that(qtest(Temp, "N+(0,)"), msg = "Either Temp is equal to or less than 0 C, NA, NaN, Inf, -Inf, empty, or a string. The equation is valid only for temperatures greater than 0 C / 32 F. Please try again.")
# only process with specified, finite values and provide an error message if the check fails
# create a numeric vector with the units of degrees Celsius
T_C <- set_units(Temp, "degree_C")
# create a numeric vector to convert from degrees Celsius to Kelvin
T_K <- T_C
# create a numeric vector with the units of Kelvin
units(T_K) <- make_units(K)
# create viscosities based on temperatures
# saturated liquid density at given temperature in degrees Celsius (SI units)
rho_SI <- density_water(Temp, "SI")
# absolute or dynamic viscosity at given temperature in degrees Celsius and density of rho (SI units)
mu_SI <- dyn_visc_water(Temp, "SI")
# kinematic viscosity at given temperature in degrees Celsius and density of rho (SI units)
nu_SI <- kin_visc_water(rho_SI, mu_SI, rho_units = "kg/m^3", mu_units = "Pa*s or kg/m/s")
k <- 1
g <- 9.80665 # m / s^2
rho <- rho_SI
nu <- nu_SI
mu <- mu_SI
# unit weight of water at given temperature in degrees Celsius and density of rho (SI units)
gamma <- unit_wt(rho = rho, units = "SI")
density_water_units <- "kg/m^3"
dyn_visc_water_units <- "Pa * s or kg/m*s"
kin_visc_water_units <- "m^2/s"
result_units <- c("m", "m^2", "m", "m", "m", "m", "m", "m/s", "m^3/s", "dimensionless", "m/m", "degrees Celsius", "Kelvin", density_water_units, dyn_visc_water_units, kin_visc_water_units, "dimensionless", "dimensionless", "m/m", "m/m", "m/m", "m", "m", "m", "m", "m^3/s", "m", "m", "pascal (N/m^2)", "pascal (N/m^2)")
} else if (units == "Eng") {
# use the temperature to determine the density & absolute and kinematic viscosities
Temp <- ifelse(is.null(Temp), 68, Temp) # degrees F
assert_that(qtest(Temp, "N+(32,)"), msg = "Either Temp is equal to or less than 32 F, NA, NaN, Inf, -Inf, empty, or a string. The equation is valid only for temperatures greater than 32 F / 0 C. Please try again.")
# only process with specified, finite values and provide an error message if the check fails
T_F <- Temp
# create a numeric vector with the units of degrees Fahrenheit
T_F <- set_units(T_F, "degree_F")
# create a numeric vector to convert from degrees Fahrenheit to Kelvin
T_K <- T_F
# create a numeric vector with the units of Kelvin
units(T_K) <- make_units(K)
# create viscosities based on temperatures
# saturated liquid density at given temperature in degrees Fahrenheit (US Customary units)
rho_Eng <- density_water(Temp, "Eng", Eng_units = "slug/ft^3")
# absolute or dynamic viscosity at given temperature in degrees Fahrenheit and density of rho (US Customary units)
mu_Eng <- dyn_visc_water(Temp, "Eng", Eng_units = "slug/ft/s")
# kinematic viscosity at given temperature in degrees Fahrenheit and density of rho (US Customary units)
nu_Eng <- kin_visc_water(rho_Eng, mu_Eng, rho_units = "slug/ft^3", mu_units = "slug/ft/s")
k <- 3.2808399 ^ (1 / 3)
g <- 9.80665 * (3937 / 1200) # ft / sec^2
gc <- 9.80665 * (3937 / 1200) # lbm-ft/lbf-sec^2
rho <- rho_Eng
nu <- nu_Eng
mu <- mu_Eng
# unit weight of water at given temperature in degrees Fahrenheit and density of rho (US Customary units)
gamma <- unit_wt(rho = rho, units = "Eng", Eng_units = "slug/ft^3")
density_water_units <- "slug/ft^3"
dyn_visc_water_units <- "slug/ft*s"
kin_visc_water_units <- "ft^2/s"
result_units <- c("ft", "ft^2", "ft", "ft", "ft", "ft", "ft", "ft/sec (fps)", "ft^3/sec (cfs)", "dimensionless", "ft/ft", "degrees Fahrenheit", "Kelvin", density_water_units, dyn_visc_water_units, kin_visc_water_units, "dimensionless", "dimensionless", "ft/ft", "ft/ft", "ft/ft", "ft", "ft", "ft", "ft", "ft^3/sec (cfs)", "ft", "ft", "lb/ft^2", "lb/ft^2")
}
if (missing(Q)) {
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
Qfun <- function(Q) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
Quse <- uniroot(Qfun, interval = c(0.0000001, 200), extendInt = "yes")
Q <- Quse$root
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(Q = Q, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
} else if (missing(n)) {
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
nfun <- function(n) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
nuse <- uniroot(nfun, interval = c(0.0000001, 200), extendInt = "yes")
n <- nuse$root
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(n = n, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
} else if (missing(b)) {
bfun <- function(b) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
buse <- uniroot(bfun, interval = c(0.0000001, 200), extendInt = "yes")
b <- buse$root
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(b = b, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
} else if (missing(y)) {
yfun <- function(y) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
yuse <- uniroot(yfun, interval = c(0.0000001, 200), extendInt = "yes")
y <- yuse$root
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(y = y, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
} else if (missing(Sf)) {
A <- b * y
P <- b + 2 * y
B <- b
R <- A / P
D <- A / B
Sffun <- function(Sf) {Q - ((((b * y)) ^ (5 / 3) * sqrt(Sf)) * (k / n) / ((b + 2 * y)) ^ (2 / 3))}
Sfuse <- uniroot(Sffun, interval = c(0.0000001, 200), extendInt = "yes")
Sf <- Sfuse$root
V <- Q / A
Re <- (rho * R * V) / mu
if (Re > 2000) {
cat("\nFlow IS in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is acceptable to use.\n\n")
} else {
cat("\nFlow is NOT in the rough turbulent zone so the Gauckler-Manning-Strickler equation\n is not acceptable to use.\n\n")
}
Fr <- V / (sqrt(g * D))
if (Fr %==% 1) {
cat("\nThis is critical flow.\n\n")
} else if (Fr < 1) {
cat("\nThis is subcritical flow.\n\n")
} else if (Fr > 1) {
cat("\nThis is supercritical flow.\n\n")
}
return(list(Sf = Sf, V = V, A = A, P = P, R = R, B = B, D = D, Re = Re, Fr = Fr))
}
}
Try the iemisc package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.