jiang
5966d039de
Reorganize algorithm packages by business responsibility, move orchestration into services, and keep database access behind pooled repositories. Harden analysis API validation, remove unsafe legacy simulation endpoints, and add regression and architecture boundary coverage. BREAKING CHANGE: legacy algorithm module paths and obsolete simulation endpoints are removed.
392 lines
12 KiB
Python
392 lines
12 KiB
Python
from pathlib import Path
|
|
|
|
import numpy as np
|
|
import pytest
|
|
import wntr
|
|
from scipy.sparse import csr_matrix, isspmatrix_csr
|
|
from scipy.sparse.csgraph import dijkstra
|
|
|
|
from app.algorithms.pressure_sensor_placement import sensitivity_placement as sensitivity
|
|
|
|
|
|
def _build_test_network() -> wntr.network.WaterNetworkModel:
|
|
wn = wntr.network.WaterNetworkModel()
|
|
wn.options.time.duration = 0
|
|
wn.add_reservoir("R1", base_head=100.0, coordinates=(-1.0, 0.0))
|
|
wn.add_junction("J0", elevation=5.0, coordinates=(0.0, 0.0))
|
|
|
|
for index in range(1, 13):
|
|
wn.add_junction(
|
|
f"J{index}",
|
|
base_demand=0.001 + index * 0.00001,
|
|
elevation=5.0 + index * 0.05,
|
|
coordinates=(float(index % 4), float(index // 4)),
|
|
)
|
|
|
|
wn.add_pipe("P0", "R1", "J0", length=100.0, diameter=0.4, roughness=110)
|
|
wn.add_pipe("P1", "J0", "J1", length=100.0, diameter=0.4, roughness=110)
|
|
for index in range(1, 11):
|
|
wn.add_pipe(
|
|
f"P{index + 1}",
|
|
f"J{index}",
|
|
f"J{index + 1}",
|
|
length=80.0 + index,
|
|
diameter=0.3,
|
|
roughness=105,
|
|
)
|
|
# J12 is connected only through a small pipe, while J11 also touches P11.
|
|
wn.add_pipe("P12", "J11", "J12", length=90.0, diameter=0.1, roughness=105)
|
|
wn.add_pipe("PX1", "J2", "J6", length=120.0, diameter=0.3, roughness=105)
|
|
wn.add_pipe("PX2", "J5", "J9", length=120.0, diameter=0.3, roughness=105)
|
|
return wn
|
|
|
|
|
|
def _prepared_selection_network(
|
|
coordinates: np.ndarray,
|
|
edges: list[tuple[int, int, float]],
|
|
) -> sensitivity._PreparedNetwork:
|
|
node_count = len(coordinates)
|
|
rows: list[int] = []
|
|
columns: list[int] = []
|
|
weights: list[float] = []
|
|
for start, end, weight in edges:
|
|
rows.extend((start, end))
|
|
columns.extend((end, start))
|
|
weights.extend((weight, weight))
|
|
coverage_graph = csr_matrix(
|
|
(weights, (rows, columns)),
|
|
shape=(node_count, node_count),
|
|
)
|
|
return sensitivity._PreparedNetwork(
|
|
node_names=tuple(f"N{index:04d}" for index in range(node_count)),
|
|
full_node_indices=np.arange(node_count, dtype=np.int64),
|
|
candidate_indices=np.arange(node_count, dtype=np.int64),
|
|
coordinates=np.asarray(coordinates, dtype=np.float64),
|
|
incidence=csr_matrix((node_count, 1), dtype=np.float64),
|
|
conductance=np.ones(1, dtype=np.float64),
|
|
roughness_response=np.ones(1, dtype=np.float64),
|
|
distance_graph=csr_matrix((node_count, node_count), dtype=np.float64),
|
|
coverage_graph=coverage_graph,
|
|
)
|
|
|
|
|
|
def test_algorithm_is_deterministic_and_runs_epanet_once(monkeypatch, tmp_path):
|
|
wn = _build_test_network()
|
|
original_run_sim = wntr.sim.EpanetSimulator.run_sim
|
|
prefixes: list[str] = []
|
|
|
|
def counted_run_sim(simulator, *args, **kwargs):
|
|
prefixes.append(str(kwargs["file_prefix"]))
|
|
return original_run_sim(simulator, *args, **kwargs)
|
|
|
|
monkeypatch.setattr(wntr.sim.EpanetSimulator, "run_sim", counted_run_sim)
|
|
monkeypatch.chdir(tmp_path)
|
|
|
|
first = sensitivity.optimize_sensor_placement(wn, sensor_num=4, min_diameter=0)
|
|
second = sensitivity.optimize_sensor_placement(wn, sensor_num=4, min_diameter=0)
|
|
|
|
assert first == second
|
|
assert len(first) == len(set(first)) == 4
|
|
assert len(prefixes) == 2
|
|
assert all(not Path(prefix).parent.exists() for prefix in prefixes)
|
|
assert not list(tmp_path.glob("temp.*"))
|
|
|
|
|
|
def test_hydraulic_simulation_keeps_only_initial_state_and_restores_duration():
|
|
wn = _build_test_network()
|
|
wn.options.time.duration = 24 * 60 * 60
|
|
|
|
results = sensitivity._run_hydraulic_simulation(wn)
|
|
|
|
assert len(results.node["head"].index) == 1
|
|
assert wn.options.time.duration == 24 * 60 * 60
|
|
|
|
|
|
def test_preparation_keeps_network_matrices_sparse():
|
|
wn = _build_test_network()
|
|
results = sensitivity._run_hydraulic_simulation(wn)
|
|
|
|
prepared = sensitivity._prepare_network(wn, results, min_diameter=0)
|
|
|
|
assert isspmatrix_csr(prepared.incidence)
|
|
assert isspmatrix_csr(prepared.distance_graph)
|
|
assert isspmatrix_csr(prepared.coverage_graph)
|
|
assert prepared.incidence.nnz <= 2 * prepared.incidence.shape[1]
|
|
assert prepared.distance_graph.nnz <= wn.num_pipes
|
|
assert prepared.coverage_graph.nnz <= 2 * wn.num_links
|
|
assert (prepared.coverage_graph != prepared.coverage_graph.T).nnz == 0
|
|
dense_incidence_bytes = int(np.prod(prepared.incidence.shape)) * 8
|
|
sparse_payload_bytes = (
|
|
prepared.incidence.data.nbytes
|
|
+ prepared.incidence.indices.nbytes
|
|
+ prepared.incidence.indptr.nbytes
|
|
)
|
|
assert sparse_payload_bytes < dense_incidence_bytes
|
|
|
|
|
|
def test_minimum_diameter_filters_installation_candidates_in_millimetres():
|
|
wn = _build_test_network()
|
|
results = sensitivity._run_hydraulic_simulation(wn)
|
|
prepared = sensitivity._prepare_network(wn, results, min_diameter=300)
|
|
candidate_names = {
|
|
prepared.node_names[index] for index in prepared.candidate_indices
|
|
}
|
|
|
|
assert "J12" not in candidate_names
|
|
assert "J11" in candidate_names
|
|
|
|
selected = sensitivity.optimize_sensor_placement(
|
|
wn,
|
|
sensor_num=4,
|
|
min_diameter=300,
|
|
)
|
|
assert set(selected) <= candidate_names
|
|
|
|
with pytest.raises(ValueError, match="候选节点少于"):
|
|
sensitivity.optimize_sensor_placement(
|
|
wn,
|
|
sensor_num=len(candidate_names) + 1,
|
|
min_diameter=300,
|
|
)
|
|
|
|
|
|
def test_sparse_estimate_preserves_dense_reference_placement_quality():
|
|
wn = _build_test_network()
|
|
results = sensitivity._run_hydraulic_simulation(wn)
|
|
prepared = sensitivity._prepare_network(wn, results, min_diameter=0)
|
|
|
|
approximate_log_sensitivity = sensitivity._estimate_log_pressure_sensitivity(
|
|
prepared
|
|
)
|
|
approximate_distance = sensitivity._estimate_hydraulic_distance_sums(prepared)
|
|
approximate_selected = sensitivity._select_sensor_nodes(
|
|
prepared,
|
|
approximate_log_sensitivity,
|
|
approximate_distance,
|
|
sensor_num=4,
|
|
)
|
|
|
|
incidence = prepared.incidence.toarray()
|
|
laplacian = (
|
|
prepared.incidence.multiply(prepared.conductance)
|
|
@ prepared.incidence.T
|
|
).toarray()
|
|
diagonal_scale = float(np.max(np.abs(np.diag(laplacian))))
|
|
laplacian += np.eye(laplacian.shape[0]) * (
|
|
diagonal_scale * np.sqrt(np.finfo(np.float64).eps)
|
|
)
|
|
response = np.linalg.solve(
|
|
laplacian,
|
|
incidence * prepared.roughness_response,
|
|
)
|
|
exact_sensitivity = np.abs(response).sum(axis=1)
|
|
|
|
exact_distances = dijkstra(
|
|
prepared.distance_graph.transpose().tocsr(),
|
|
directed=True,
|
|
indices=prepared.full_node_indices,
|
|
return_predecessors=False,
|
|
)[:, prepared.full_node_indices]
|
|
exact_distances[~np.isfinite(exact_distances)] = 0.0
|
|
exact_distance = exact_distances.sum(axis=0)
|
|
exact_selected = sensitivity._select_sensor_nodes(
|
|
prepared,
|
|
np.log(np.maximum(exact_sensitivity, np.finfo(np.float64).tiny)),
|
|
exact_distance,
|
|
sensor_num=4,
|
|
)
|
|
|
|
exact_score = exact_sensitivity * exact_distance
|
|
node_index = {
|
|
node_name: index for index, node_name in enumerate(prepared.node_names)
|
|
}
|
|
approximate_objective = sum(
|
|
exact_score[node_index[node_name]] for node_name in approximate_selected
|
|
)
|
|
exact_objective = sum(
|
|
exact_score[node_index[node_name]] for node_name in exact_selected
|
|
)
|
|
|
|
assert approximate_objective / exact_objective >= 0.95
|
|
|
|
|
|
def test_mixed_coverage_avoids_candidate_density_bias():
|
|
dense_west = np.linspace(0.0, 2.0, 200)
|
|
sparse_east = np.linspace(3.0, 10.0, 20)
|
|
x_coordinates = np.concatenate((dense_west, sparse_east))
|
|
coordinates = np.column_stack(
|
|
(x_coordinates, np.zeros(len(x_coordinates), dtype=np.float64))
|
|
)
|
|
ordered = np.argsort(x_coordinates)
|
|
edges = [
|
|
(
|
|
int(start),
|
|
int(end),
|
|
float(x_coordinates[end] - x_coordinates[start]),
|
|
)
|
|
for start, end in zip(ordered[:-1], ordered[1:])
|
|
]
|
|
prepared = _prepared_selection_network(coordinates, edges)
|
|
log_scores = np.linspace(4.0, 0.0, len(coordinates))
|
|
distance_sums = np.ones(len(coordinates), dtype=np.float64)
|
|
|
|
selected = sensitivity._select_sensor_nodes(
|
|
prepared,
|
|
log_scores,
|
|
distance_sums,
|
|
sensor_num=6,
|
|
)
|
|
name_to_position = {
|
|
name: position for position, name in enumerate(prepared.node_names)
|
|
}
|
|
selected_positions = [name_to_position[name] for name in selected]
|
|
new_metrics = sensitivity._geographic_coverage_metrics(
|
|
coordinates,
|
|
selected_positions,
|
|
)
|
|
|
|
legacy_labels, _centers = sensitivity._cluster_labels(
|
|
coordinates,
|
|
6,
|
|
random_seed=sensitivity._RANDOM_SEED + 2,
|
|
)
|
|
legacy_positions: list[int] = []
|
|
represented: set[int] = set()
|
|
for position in np.argsort(-log_scores):
|
|
label = int(legacy_labels[position])
|
|
if label in represented:
|
|
continue
|
|
represented.add(label)
|
|
legacy_positions.append(int(position))
|
|
legacy_metrics = sensitivity._geographic_coverage_metrics(
|
|
coordinates,
|
|
legacy_positions,
|
|
)
|
|
|
|
assert new_metrics[0] <= legacy_metrics[0] * 0.6
|
|
assert new_metrics[2] >= legacy_metrics[2] * 1.5
|
|
assert max(x_coordinates[selected_positions]) >= 9.0
|
|
|
|
|
|
def test_disconnected_components_each_receive_a_sensor_when_slots_allow():
|
|
coordinates = np.asarray(
|
|
[
|
|
(0.0, 0.0),
|
|
(1.0, 0.0),
|
|
(0.0, 0.01),
|
|
(1.0, 0.01),
|
|
]
|
|
)
|
|
prepared = _prepared_selection_network(
|
|
coordinates,
|
|
[(0, 1, 1.0), (2, 3, 1.0)],
|
|
)
|
|
|
|
selected = sensitivity._select_sensor_nodes(
|
|
prepared,
|
|
np.asarray([10.0, 9.0, 8.0, 7.0]),
|
|
np.ones(4),
|
|
sensor_num=2,
|
|
)
|
|
|
|
assert len(set(selected) & {"N0000", "N0001"}) == 1
|
|
assert len(set(selected) & {"N0002", "N0003"}) == 1
|
|
|
|
|
|
def test_overlapping_components_respect_global_geographic_spacing():
|
|
coordinates = np.asarray(
|
|
[
|
|
(0.0, 0.0),
|
|
(10.0, 0.0),
|
|
(0.1, 0.0),
|
|
(10.1, 0.0),
|
|
]
|
|
)
|
|
prepared = _prepared_selection_network(
|
|
coordinates,
|
|
[(0, 1, 10.0), (2, 3, 10.0)],
|
|
)
|
|
|
|
selected = sensitivity._select_sensor_nodes(
|
|
prepared,
|
|
np.asarray([10.0, 1.0, 9.0, 0.0]),
|
|
np.ones(4),
|
|
sensor_num=2,
|
|
)
|
|
selected_positions = [prepared.node_names.index(name) for name in selected]
|
|
minimum_gap = sensitivity._geographic_coverage_metrics(
|
|
coordinates,
|
|
selected_positions,
|
|
)[2]
|
|
|
|
assert minimum_gap >= 0.9
|
|
|
|
|
|
def test_component_quota_prefers_longer_networks_when_slots_are_limited():
|
|
coordinates = np.asarray(
|
|
[
|
|
(0.0, 0.0),
|
|
(10.0, 0.0),
|
|
(20.0, 0.0),
|
|
(25.0, 0.0),
|
|
(30.0, 0.0),
|
|
(31.0, 0.0),
|
|
]
|
|
)
|
|
prepared = _prepared_selection_network(
|
|
coordinates,
|
|
[(0, 1, 10.0), (2, 3, 5.0), (4, 5, 1.0)],
|
|
)
|
|
|
|
selected = sensitivity._select_sensor_nodes(
|
|
prepared,
|
|
np.asarray([1.0, 1.0, 2.0, 2.0, 100.0, 100.0]),
|
|
np.ones(6),
|
|
sensor_num=2,
|
|
)
|
|
|
|
assert set(selected) <= {"N0000", "N0001", "N0002", "N0003"}
|
|
assert len(set(selected) & {"N0000", "N0001"}) == 1
|
|
assert len(set(selected) & {"N0002", "N0003"}) == 1
|
|
|
|
|
|
def test_duplicate_coordinates_use_topology_and_return_exact_count():
|
|
coordinates = np.zeros((6, 2), dtype=np.float64)
|
|
prepared = _prepared_selection_network(
|
|
coordinates,
|
|
[(index, index + 1, 1.0) for index in range(5)],
|
|
)
|
|
log_scores = np.linspace(6.0, 1.0, 6)
|
|
|
|
first = sensitivity._select_sensor_nodes(
|
|
prepared,
|
|
log_scores,
|
|
np.ones(6),
|
|
sensor_num=4,
|
|
)
|
|
second = sensitivity._select_sensor_nodes(
|
|
prepared,
|
|
log_scores,
|
|
np.ones(6),
|
|
sensor_num=4,
|
|
)
|
|
|
|
assert first == second
|
|
assert len(first) == len(set(first)) == 4
|
|
|
|
|
|
@pytest.mark.parametrize(
|
|
("sensor_num", "min_diameter", "message"),
|
|
[
|
|
(0, 0, "监测点数量必须大于 0"),
|
|
(1, -1, "最小管径不能小于 0"),
|
|
],
|
|
)
|
|
def test_algorithm_rejects_invalid_parameters(sensor_num, min_diameter, message):
|
|
with pytest.raises(ValueError, match=message):
|
|
sensitivity.optimize_sensor_placement(
|
|
_build_test_network(),
|
|
sensor_num=sensor_num,
|
|
min_diameter=min_diameter,
|
|
)
|