Manning's n Bulk Sensitivity Analysis¶
This notebook demonstrates bulk Manning's n sensitivity analysis for HEC-RAS 2D models with spatially variable roughness. The workflow:
- Extracts a 2D project with Manning's n land cover regions
- Defines minimum/maximum parameter ranges for each roughness class
- Creates modified geometries with min and max Manning's n values
- Executes HEC-RAS for current, minimum, and maximum scenarios
- Extracts water surface elevation at a point of interest
- Quantifies sensitivity range and visualizes results
Manning's n Physical Context¶
| Land Cover | Typical Range |
|---|---|
| Paved/open space | 0.020 - 0.060 |
| Mowed grass/parks | 0.030 - 0.080 |
| Residential | 0.050 - 0.120 |
| Urban (mixed) | 0.060 - 0.150 |
| Forest/trees | 0.080 - 0.200 |
Reference¶
- Chow (1959) Open-Channel Hydraulics: Table 5-6
- HEC-RAS 2D User Manual: Section 3.7
- USGS TWI 3-A1: Roughness Characteristics of Natural Channels
Setup¶
# Development mode toggle
USE_LOCAL_SOURCE = False
if USE_LOCAL_SOURCE:
import sys
from pathlib import Path
local_path = str(Path.cwd().parent)
if local_path not in sys.path:
sys.path.insert(0, local_path)
print(f"LOCAL SOURCE MODE: Loading from {local_path}/ras_commander")
else:
from pathlib import Path
print("PIP PACKAGE MODE: Loading installed ras-commander")
from ras_commander import (
init_ras_project, RasExamples, RasPlan, RasCmdr,
RasGeo, HdfMesh, HdfResultsMesh
)
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from shapely.geometry import Point
import ras_commander
print(f"Loaded: {ras_commander.__file__}")
Parameters¶
Configure these values to customize for your project.
# Project Configuration
PROJECT_NAME = "Muncie"
RAS_VERSION = "7.0"
# Template plan with 2D Manning's n regions
TEMPLATE_PLAN = "04" # "Unsteady Run with 2D 50ft User n Value Regions"
# Point of interest for result extraction (State Plane Indiana East, ft)
POINT_OF_INTEREST = (408350.0, 1802550.0)
# Which Manning's categories to adjust
INCLUDE_REGIONAL_OVERRIDES = True
INCLUDE_BASE_OVERRIDES = True
# Execution settings
NUM_CORES = 2
Define Manning's n Value Ranges¶
These ranges are based on literature values for the Muncie model's urban land cover
classification. The building class uses a high obstruction value and is excluded
from sensitivity variation.
def create_manning_minmax_df():
"""
Create min/max Manning's n ranges for Muncie's urban land cover types.
The 'building' class is an obstruction (n=100) and is held constant.
All other classes vary based on published literature ranges.
"""
manning_data = [
{"Land Cover Name": "building", "min_n": 100.0, "max_n": 100.0},
{"Land Cover Name": "medium density residential", "min_n": 0.050, "max_n": 0.120},
{"Land Cover Name": "open space", "min_n": 0.020, "max_n": 0.060},
{"Land Cover Name": "park", "min_n": 0.030, "max_n": 0.080},
{"Land Cover Name": "trees", "min_n": 0.080, "max_n": 0.200},
{"Land Cover Name": "urban", "min_n": 0.060, "max_n": 0.150},
]
df = pd.DataFrame(manning_data)
df['mid_n'] = (df['min_n'] + df['max_n']) / 2
print(f"Manning's n ranges for {len(df)} land cover types:")
print(df.to_string(index=False))
return df
manning_minmax_df = create_manning_minmax_df()
Extract Project and Initialize¶
project_folder = RasExamples.extract_project(PROJECT_NAME, suffix="710")
ras = init_ras_project(project_folder, RAS_VERSION)
print(f"\nProject: {ras.project_name}")
print(f"Plans:")
print(ras.plan_df[['plan_number', 'Plan Title', 'geometry_number']].to_string())
Inspect Current Manning's n Values¶
template_geom = ras.plan_df.loc[
ras.plan_df['plan_number'] == TEMPLATE_PLAN, 'geometry_number'
].values[0]
geom_path = ras.geom_df.loc[
ras.geom_df['geom_number'] == template_geom, 'full_path'
].values[0]
original_base = RasGeo.get_mannings_baseoverrides(geom_path)
original_region = RasGeo.get_mannings_regionoverrides(geom_path)
print("=== Base Manning's n Overrides ===")
print(original_base.to_string())
print(f"\n=== Regional Manning's n Overrides ===")
if not original_region.empty:
print(original_region.to_string())
else:
print("None found")
Execute Template Plan (Baseline)¶
print(f"Executing template plan {TEMPLATE_PLAN} (current Manning's values)...")
result = RasCmdr.compute_plan(TEMPLATE_PLAN, num_cores=NUM_CORES, clear_geompre=True)
print(f"Result: {result}")
Create Min/Max Scenarios¶
Clone the template plan and geometry, then modify Manning's n values to the minimum and maximum of the defined ranges.
def create_modified_scenario(name, shortid, manning_minmax_df, use_min=False, use_max=False):
"""
Clone template plan/geometry and apply modified Manning's n values.
Args:
name: Scenario name
shortid: Plan short identifier
manning_minmax_df: DataFrame with min_n and max_n columns
use_min: Apply minimum values
use_max: Apply maximum values
Returns:
dict with scenario metadata
"""
new_plan = RasPlan.clone_plan(TEMPLATE_PLAN, new_plan_shortid=shortid, ras_object=ras)
new_geom = RasPlan.clone_geom(template_geom, ras_object=ras)
RasPlan.set_geom(new_plan, new_geom, ras_object=ras)
new_geom_path = ras.geom_df.loc[
ras.geom_df['geom_number'] == new_geom, 'full_path'
].values[0]
# Modify base overrides
if INCLUDE_BASE_OVERRIDES:
modified_base = original_base.copy()
for idx, row in modified_base.iterrows():
match = manning_minmax_df[
manning_minmax_df['Land Cover Name'] == row['Land Cover Name']
]
if not match.empty:
if use_min:
modified_base.loc[idx, 'Base Mannings n Value'] = match['min_n'].values[0]
elif use_max:
modified_base.loc[idx, 'Base Mannings n Value'] = match['max_n'].values[0]
RasGeo.set_mannings_baseoverrides(new_geom_path, modified_base)
# Modify regional overrides
if INCLUDE_REGIONAL_OVERRIDES and not original_region.empty:
modified_region = original_region.copy()
for idx, row in modified_region.iterrows():
match = manning_minmax_df[
manning_minmax_df['Land Cover Name'] == row['Land Cover Name']
]
if not match.empty:
if use_min:
modified_region.loc[idx, 'MainChannel'] = match['min_n'].values[0]
elif use_max:
modified_region.loc[idx, 'MainChannel'] = match['max_n'].values[0]
RasGeo.set_mannings_regionoverrides(new_geom_path, modified_region)
print(f" Created: {name} (plan={new_plan}, geom={new_geom})")
return {'name': name, 'plan_number': new_plan, 'geom_number': new_geom, 'shortid': shortid}
print("Creating sensitivity scenarios...")
min_scenario = create_modified_scenario("Minimum", "Min_n", manning_minmax_df, use_min=True)
max_scenario = create_modified_scenario("Maximum", "Max_n", manning_minmax_df, use_max=True)
scenarios = [
{'name': 'Current', 'plan_number': TEMPLATE_PLAN, 'geom_number': template_geom, 'shortid': 'Current'},
min_scenario,
max_scenario,
]
print(f"\n{len(scenarios)} scenarios ready.")
Execute Min/Max Plans¶
plans_to_run = [min_scenario['plan_number'], max_scenario['plan_number']]
print(f"Executing plans: {plans_to_run}")
results = RasCmdr.compute_parallel(
plan_number=plans_to_run,
max_workers=2,
num_cores=NUM_CORES,
clear_geompre=True
)
for plan, res in results.items():
print(f" Plan {plan}: {res}")
Extract Results at Point of Interest¶
poi = Point(POINT_OF_INTEREST[0], POINT_OF_INTEREST[1])
# Find nearest mesh cell
mesh_cells = HdfMesh.get_mesh_cell_points(TEMPLATE_PLAN, ras_object=ras)
distances = mesh_cells.geometry.apply(lambda g: g.distance(poi))
nearest_idx = distances.idxmin()
nearest_cell = mesh_cells.loc[nearest_idx]
cell_id = nearest_cell['cell_id']
mesh_name = nearest_cell['mesh_name']
print(f"Point of Interest: {POINT_OF_INTEREST}")
print(f"Nearest cell: id={cell_id}, distance={distances[nearest_idx]:.1f} ft")
print(f"Mesh: {mesh_name}")
# Extract water surface time series for each scenario
all_results = {}
for scenario in scenarios:
plan_num = scenario['plan_number']
name = scenario['name']
try:
results_xr = HdfResultsMesh.get_mesh_cells_timeseries(plan_num, ras_object=ras)
ws_data = results_xr[mesh_name]['Water Surface'].sel(cell_id=int(cell_id))
ws_df = pd.DataFrame({
'time': ws_data.time.values,
'water_surface': ws_data.values
})
max_ws = ws_df['water_surface'].max()
all_results[name] = {'df': ws_df, 'max_ws': max_ws, 'plan': plan_num}
print(f" {name}: Max WSE = {max_ws:.2f} ft")
except Exception as e:
print(f" {name}: ERROR - {e}")
print(f"\nSuccessfully extracted {len(all_results)}/{len(scenarios)} scenarios.")
Sensitivity Analysis Results¶
if len(all_results) == 3:
current_ws = all_results['Current']['max_ws']
min_ws = all_results['Minimum']['max_ws']
max_ws = all_results['Maximum']['max_ws']
sensitivity_range = max_ws - min_ws
print("=" * 50)
print("MANNING'S n BULK SENSITIVITY SUMMARY")
print("=" * 50)
print(f" Minimum n scenario WSE: {min_ws:.2f} ft")
print(f" Current n scenario WSE: {current_ws:.2f} ft")
print(f" Maximum n scenario WSE: {max_ws:.2f} ft")
print(f"")
print(f" Sensitivity Range: {sensitivity_range:.2f} ft")
print(f" Current vs Min: +{current_ws - min_ws:.2f} ft")
print(f" Max vs Current: +{max_ws - current_ws:.2f} ft")
if sensitivity_range > 0:
position = (current_ws - min_ws) / sensitivity_range * 100
print(f" Current position in range: {position:.0f}%")
print("=" * 50)
else:
print(f"Only {len(all_results)} scenarios extracted. Check execution results above.")
Visualize Results¶
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
# Time series plot
ax1 = axes[0]
colors = {'Current': 'black', 'Minimum': 'blue', 'Maximum': 'red'}
styles = {'Current': '-', 'Minimum': '--', 'Maximum': '--'}
for name, result in all_results.items():
df = result['df']
ax1.plot(df['time'], df['water_surface'],
label=f"{name} (max={result['max_ws']:.2f} ft)",
color=colors[name], linestyle=styles[name], linewidth=1.5)
ax1.set_xlabel('Time')
ax1.set_ylabel('Water Surface Elevation (ft)')
ax1.set_title(f"WSE Sensitivity to Manning's n at Cell {cell_id}")
ax1.legend()
ax1.grid(True, alpha=0.3)
ax1.tick_params(axis='x', rotation=45)
# Bar chart
ax2 = axes[1]
names = list(all_results.keys())
values = [all_results[n]['max_ws'] for n in names]
bar_colors = [colors[n] for n in names]
bars = ax2.bar(names, values, color=bar_colors, alpha=0.7, edgecolor='black')
for bar, val in zip(bars, values):
ax2.text(bar.get_x() + bar.get_width()/2, val + 0.02,
f'{val:.2f}', ha='center', va='bottom', fontsize=10)
ax2.set_ylabel('Maximum WSE (ft)')
ax2.set_title("Peak WSE by Manning's n Scenario")
ax2.grid(axis='y', alpha=0.3)
# Set y-axis to show differences clearly
all_vals = values
y_margin = max(0.5, (max(all_vals) - min(all_vals)) * 0.3)
ax2.set_ylim(min(all_vals) - y_margin, max(all_vals) + y_margin)
plt.tight_layout()
plt.savefig(Path(project_folder) / "mannings_bulk_sensitivity.png", dpi=150, bbox_inches='tight')
plt.show()
print(f"Plot saved to: {Path(project_folder) / 'mannings_bulk_sensitivity.png'}")
Manning's n Comparison Table¶
# Show how Manning's values changed across scenarios
comparison = manning_minmax_df[manning_minmax_df['Land Cover Name'] != 'building'].copy()
comparison = comparison.rename(columns={'min_n': 'Min Scenario', 'max_n': 'Max Scenario', 'mid_n': 'Mid'})
# Add current values
current_vals = []
for _, row in comparison.iterrows():
match = original_base[original_base['Land Cover Name'] == row['Land Cover Name']]
if not match.empty:
current_vals.append(match['Base Mannings n Value'].values[0])
else:
current_vals.append(np.nan)
comparison['Current'] = current_vals
comparison['Range'] = comparison['Max Scenario'] - comparison['Min Scenario']
print("Manning's n values by scenario:")
print(comparison[['Land Cover Name', 'Min Scenario', 'Current', 'Max Scenario', 'Range']].to_string(index=False))
Interpretation¶
The bulk sensitivity analysis shows the total envelope of uncertainty in water surface elevation due to Manning's n selection. Key takeaways:
- Range magnitude: The difference between min and max WSE indicates overall model sensitivity to roughness parameters
- Current position: Where the current calibration sits within the uncertainty envelope
- Hydraulic significance: If the range exceeds 0.5 ft, Manning's n selection materially affects flood mapping results
For individual land cover sensitivity (one-at-a-time analysis), see
notebook 711_mannings_sensitivity_multi_interval.ipynb.