💻 Fortran Source Code Library
We currently offer 172 open-source, production-grade Fortran codes for offline testing. Run calculations locally on your own machine, view code structure, read technical explanations, and download compilation packages including sample input files.
Manning Equation — Open Channel
Core Numerical Engine in Fortran 90 • 60 total downloads
manning_equation.f90
! =========================================================================
! Source File: manning_equation.f90
! =========================================================================
! ============================================================
! ThermoFluidCalc — Manning's Equation Open Channel Flow
! Solve for Q (discharge) or y_n (normal depth)
! Supports: rectangular, trapezoidal, triangular, circular
! ============================================================
program manning_equation
implicit none
real(8), parameter :: g = 9.81d0
real(8), parameter :: PI = 3.141592653589793d0
integer :: shape, solve_mode
real(8) :: b, z, D, mann_n, S0, known_val
real(8) :: A, P_wet, R_h, T, y_h, V, Q, Fr, E, y_c, y_n, Re
real(8) :: nu_water
integer :: i, n_points
real(8) :: y_plot, Q_plot, y_max
character(len=20) :: shape_name, solve_name, regime_name
! Read inputs
read(*,*) shape ! 1=rect, 2=trap, 3=tri, 4=circular
read(*,*) b ! bottom width [m] (rect/trap)
read(*,*) z ! side slope H:V (trap/tri)
read(*,*) D ! diameter [m] (circular)
read(*,*) mann_n ! Manning's n
read(*,*) S0 ! bed slope
read(*,*) solve_mode ! 1=find Q from y, 2=find y from Q
read(*,*) known_val ! depth y [m] or discharge Q [m³/s]
nu_water = 1.0d-6 ! kinematic viscosity of water [m²/s]
! Map shape names
select case (shape)
case (1)
shape_name = 'Rectangular'
case (2)
shape_name = 'Trapezoidal'
case (3)
shape_name = 'Triangular'
case (4)
shape_name = 'Circular'
case default
shape_name = 'Unknown'
end select
! Map solve mode names
if (solve_mode == 1) then
solve_name = 'Q from Depth y_n'
else
solve_name = 'Depth y_n from Q'
end if
if (solve_mode == 1) then
! === MODE 1: Given depth y, find Q ===
y_n = known_val
call geometry(shape, b, z, D, y_n, A, P_wet, R_h, T)
if (A <= 0.0d0 .or. P_wet <= 0.0d0) then
write(*,'(A)') "ERROR: Invalid geometry or zero depth"
stop
end if
Q = (1.0d0/mann_n) * A * R_h**(2.0d0/3.0d0) * sqrt(S0)
V = Q / A
y_h = A / T
Fr = V / sqrt(g * y_h)
E = y_n + V**2 / (2.0d0*g)
Re = V * R_h / nu_water
! Critical depth by bisection
call critical_depth(shape, b, z, D, Q, y_c)
else
! === MODE 2: Given Q, find normal depth y_n by bisection ===
Q = known_val
call normal_depth(shape, b, z, D, mann_n, S0, Q, y_n)
call geometry(shape, b, z, D, y_n, A, P_wet, R_h, T)
V = Q / A
y_h = A / T
Fr = V / sqrt(g * y_h)
E = y_n + V**2 / (2.0d0*g)
Re = V * R_h / nu_water
call critical_depth(shape, b, z, D, Q, y_c)
end if
! Determine flow regime
if (Fr < 0.95d0) then
regime_name = 'Subcritical'
else if (Fr > 1.05d0) then
regime_name = 'Supercritical'
else
regime_name = 'Critical'
end if
! === OUTPUT ASCII REPORT ===
write(*,'(A)') '============================================================'
write(*,'(A)') ' MANNINGS EQUATION OPEN CHANNEL FLOW CALCULATOR'
write(*,'(A)') '============================================================'
write(*,*)
write(*,'(A)') '--- INPUT CONFIGURATION ---'
write(*,'(A,A)') ' Channel Shape = ', trim(shape_name)
select case (shape)
case (1)
write(*,'(A,F12.4,A)') ' Bottom Width (b) = ', b, ' m'
case (2)
write(*,'(A,F12.4,A)') ' Bottom Width (b) = ', b, ' m'
write(*,'(A,F12.4)') ' Side Slope (z) = ', z
case (3)
write(*,'(A,F12.4)') ' Side Slope (z) = ', z
case (4)
write(*,'(A,F12.4,A)') ' Pipe Diameter (D) = ', D, ' m'
end select
write(*,'(A,F12.4)') ' Mannings Roughness (n) = ', mann_n
write(*,'(A,F12.6,A)') ' Bed Slope (S0) = ', S0, ' m/m'
write(*,'(A,A)') ' Calculation Mode = ', trim(solve_name)
write(*,*)
write(*,'(A)') '--- HYDRAULIC RESULTS ---'
write(*,'(A,F12.4,A)') ' Normal Depth (y_n) = ', y_n, ' m'
write(*,'(A,F12.4,A)') ' Flow Area (A) = ', A, ' m2'
write(*,'(A,F12.4,A)') ' Wetted Perimeter (P) = ', P_wet, ' m'
write(*,'(A,F12.4,A)') ' Hydraulic Radius (R_h) = ', R_h, ' m'
write(*,'(A,F12.4,A)') ' Top Width (T) = ', T, ' m'
write(*,'(A,F12.4,A)') ' Hydraulic Depth (y_h) = ', y_h, ' m'
write(*,'(A,F12.4,A)') ' Mean Velocity (V) = ', V, ' m/s'
write(*,'(A,F12.4,A)') ' Discharge (Q) = ', Q, ' m3/s'
write(*,'(A,F12.4)') ' Froude Number (Fr) = ', Fr
write(*,'(A,F12.4,A)') ' Specific Energy (E) = ', E, ' m'
write(*,'(A,ES12.4)') ' Reynolds Number (Re) = ', Re
write(*,'(A,F12.4,A)') ' Critical Depth (y_c) = ', y_c, ' m'
write(*,'(A,A)') ' Flow Regime = ', trim(regime_name)
write(*,*)
write(*,'(A)') '--- EQUATIONS USED ---'
write(*,'(A)') ' Manning Equation: Q = (1/n) * A * R_h^(2/3) * sqrt(S0)'
write(*,'(A)') ' Hydraulic Radius: R_h = A / P'
write(*,'(A)') ' Froude Number: Fr = V / sqrt(g * A / T)'
write(*,'(A)') ' Specific Energy: E = y + V^2 / (2g)'
write(*,'(A)') ' Critical Depth: Solved by A^3 / T = Q^2 / g'
write(*,'(A)') '============================================================'
write(*,*)
! === CHART DATA (Rating curve Q vs y) ===
write(*,'(A)') "--- CHART_DATA ---"
if (shape == 4) then
y_max = D * 0.95d0
else
y_max = max(y_n * 2.0d0, y_c * 2.0d0)
if (y_max < 0.1d0) y_max = 0.5d0
end if
n_points = 30
do i = 1, n_points
y_plot = y_max * real(i, 8) / real(n_points, 8)
call geometry(shape, b, z, D, y_plot, A, P_wet, R_h, T)
if (A > 0.0d0 .and. P_wet > 0.0d0) then
Q_plot = (1.0d0/mann_n) * A * R_h**(2.0d0/3.0d0) * sqrt(S0)
else
Q_plot = 0.0d0
end if
write(*,'(F10.4,A,F10.4)') y_plot, ",", Q_plot
end do
write(*,'(A)') "--- END_CHART_DATA ---"
contains
! ---- Cross-section geometry ----
subroutine geometry(shape, b, z, D, y, A, P, R, T)
integer, intent(in) :: shape
real(8), intent(in) :: b, z, D, y
real(8), intent(out) :: A, P, R, T
real(8) :: theta
select case (shape)
case (1) ! Rectangular
A = b * y
P = b + 2.0d0*y
T = b
case (2) ! Trapezoidal
A = (b + z*y) * y
P = b + 2.0d0*y*sqrt(1.0d0 + z**2)
T = b + 2.0d0*z*y
case (3) ! Triangular
A = z * y**2
P = 2.0d0*y*sqrt(1.0d0 + z**2)
T = 2.0d0*z*y
case (4) ! Circular
if (y >= D) then
theta = 2.0d0*PI
else
theta = 2.0d0 * acos(1.0d0 - 2.0d0*y/D)
end if
A = (D**2/8.0d0) * (theta - sin(theta))
P = D * theta / 2.0d0
T = D * sin(theta/2.0d0)
case default
A = b * y; P = b + 2.0d0*y; T = b
end select
if (P > 0.0d0) then
R = A / P
else
R = 0.0d0
end if
end subroutine geometry
! ---- Normal depth by bisection ----
subroutine normal_depth(shape, b, z, D, n_mann, S0, Q_target, yn)
integer, intent(in) :: shape
real(8), intent(in) :: b, z, D, n_mann, S0, Q_target
real(8), intent(out) :: yn
real(8) :: y_lo, y_hi, y_mid, Q_mid, A_, P_, R_, T_
integer :: iter
y_lo = 1.0d-4
if (shape == 4) then
y_hi = D * 0.99d0
else
y_hi = 50.0d0
end if
do iter = 1, 100
y_mid = (y_lo + y_hi) / 2.0d0
call geometry(shape, b, z, D, y_mid, A_, P_, R_, T_)
if (A_ > 0.0d0 .and. P_ > 0.0d0) then
Q_mid = (1.0d0/n_mann) * A_ * R_**(2.0d0/3.0d0) * sqrt(S0)
else
Q_mid = 0.0d0
end if
if (Q_mid < Q_target) then
y_lo = y_mid
else
y_hi = y_mid
end if
if ((y_hi - y_lo) < 1.0d-6) exit
end do
yn = (y_lo + y_hi) / 2.0d0
end subroutine normal_depth
! ---- Critical depth by bisection: A³/T = Q²/g ----
subroutine critical_depth(shape, b, z, D, Q, yc)
integer, intent(in) :: shape
real(8), intent(in) :: b, z, D, Q
real(8), intent(out) :: yc
real(8) :: y_lo, y_hi, y_mid, A_, P_, R_, T_, f_mid
integer :: iter
y_lo = 1.0d-4
if (shape == 4) then
y_hi = D * 0.99d0
else
y_hi = 50.0d0
end if
do iter = 1, 100
y_mid = (y_lo + y_hi) / 2.0d0
call geometry(shape, b, z, D, y_mid, A_, P_, R_, T_)
if (T_ > 0.0d0) then
f_mid = A_**3 / T_ - Q**2 / g
else
f_mid = -Q**2 / g
end if
if (f_mid < 0.0d0) then
y_lo = y_mid
else
y_hi = y_mid
end if
if ((y_hi - y_lo) < 1.0d-6) exit
end do
yc = (y_lo + y_hi) / 2.0d0
end subroutine critical_depth
end program manning_equation
Solver Description
Open channel flow under gravity is modeled by Manning's Equation. It relates flow velocity ($V$) and volumetric flow ($Q$) to wetted perimeter roughness and hydraulic gradient:
Key Numerical Methods & Architecture
- Input Redirection: Reads parameters sequentially from standard input (`stdin`) using Fortran sequential read (`read(*,*)`), ensuring modular integration.
- Modular Design: Formulated using pure mathematical routines, separation of equations from output formatting, and precise numerical solvers (e.g. bisection, Newton-Raphson).
- Standard Compliant: Written in clean, standards-compliant Fortran 90 to ensure cross-compiler compatibility.
🛠️ Local Compilation
To test this code on your machine, compile the source code file(s) using a standard Fortran compiler (e.g., `gfortran`).
Compilation Command:
gfortran -O3 manning_equation.f90 -o manning_equation_calc
Execution Command:
Execute the program by feeding the sample input file into the program using stdin redirection:
manning_equation_calc < input.txt
📥 Downloads & Local Files
Preview of the required input file (input.txt):
! Channel shape (1=Trapezoidal, 2=Circular, 3=Rectangular, 4=Triangular)
1
! Bottom width b [m]
2.0
! Side slope z (z:1 horizontal:vertical)
1.5
! Diameter D [m]
0.0
! Manning roughness coefficient n
0.025
! Channel bed slope S0
0.001
! Solve mode (1=Solve Normal Depth, 2=Solve Discharge Q)
1
! Known parameter (Discharge Q [m3/s] or Depth y [m])
5.0
1
! Bottom width b [m]
2.0
! Side slope z (z:1 horizontal:vertical)
1.5
! Diameter D [m]
0.0
! Manning roughness coefficient n
0.025
! Channel bed slope S0
0.001
! Solve mode (1=Solve Normal Depth, 2=Solve Discharge Q)
1
! Known parameter (Discharge Q [m3/s] or Depth y [m])
5.0