Update with WMH's code
This commit is contained in:
DingZQ
13
globals.py
13
globals.py
@@ -2,7 +2,6 @@
|
||||
# reservoir basic height
|
||||
RESERVOIR_BASIC_HEIGHT = float(250.35)
|
||||
PATTERN_TIME_STEP = None # 浮点数
|
||||
|
||||
# 实时数据类：element_id和api_query_id对应
|
||||
reservoirs_id = {}
|
||||
tanks_id = {}
|
||||
@@ -11,22 +10,18 @@ variable_pumps_id = {}
|
||||
pressure_id = {}
|
||||
demand_id = {}
|
||||
quality_id = {}
|
||||
|
||||
# 实时数据类：pattern_id和api_query_id对应
|
||||
source_outflow_pattern_id = {}
|
||||
realtime_pipe_flow_pattern_id = {}
|
||||
pipe_flow_region_patterns = {} # 根据realtime的pipe_flow，对non_realtime的demand进行分区
|
||||
|
||||
# 分区查询
|
||||
source_outflow_region = {} # 以绑定的管段作为value
|
||||
source_outflow_region_id = {} # 以api_query_id作为value
|
||||
source_outflow_region_patterns = {} # 以associated_pattern作为value
|
||||
# 非实时数据的pattern
|
||||
non_realtime_region_patterns = {} # 基于source_outflow_region进行区分
|
||||
|
||||
realtime_region_pipe_flow_and_demand_id = {} # 基于source_outflow_region搜索该分区中的实时pipe_flow和demand的api_query_id，后续用region的流量 - 实时流量计的流量
|
||||
realtime_region_pipe_flow_and_demand_patterns = {} # 基于source_outflow_region搜索该分区中的实时pipe_flow和demand的associated_pattern，后续用region的流量 - 实时流量计的流量
|
||||
|
||||
# ---------------------------------------------------------
|
||||
# influxdb_api.py中的全局变量
|
||||
# 全局变量，用于存储不同类型的realtime api_query_id
|
||||
@@ -39,11 +34,9 @@ pipe_flow_realtime_ids = []
|
||||
pressure_realtime_ids = []
|
||||
demand_realtime_ids = []
|
||||
quality_realtime_ids = []
|
||||
|
||||
# transmission_frequency的最大值
|
||||
transmission_frequency = None
|
||||
hydraulic_timestep = None # 时间字符串
|
||||
|
||||
reservoir_liquid_level_non_realtime_ids = []
|
||||
tank_liquid_level_non_realtime_ids = []
|
||||
fixed_pump_non_realtime_ids = []
|
||||
@@ -54,3 +47,9 @@ pressure_non_realtime_ids = []
|
||||
demand_non_realtime_ids = []
|
||||
quality_non_realtime_ids = []
|
||||
|
||||
# api_query_id和associated_element_id对应，不包含液位和泵
|
||||
scheme_source_outflow_ids = {}
|
||||
scheme_pipe_flow_ids = {}
|
||||
scheme_pressure_ids = {}
|
||||
scheme_demand_ids = {}
|
||||
scheme_quality_ids = {}
|
||||
|
||||
737
influxdb_api.py
737
influxdb_api.py
File diff suppressed because it is too large
Load Diff
@@ -56,11 +56,11 @@ ZBBDFQJL000035,demand,ZBBDFQJL000035,ZhongYangXinDu,,ZBBGXSZW000377,P17021,,,,95
|
||||
J06186,demand,J06186,XinHaiJiaYuan,,ZBBGXSZW000377,P17021,,,,9525,non_realtime,6:00:00,106.4067,29.8278
|
||||
ZBBDFQJL000052,demand,ZBBDFQJL000052,DongFengJie,,ZBBGXSZW000377,P17021,,,,9526,non_realtime,6:00:00,106.3678,29.766
|
||||
ZBBDFQJL000049,demand,ZBBDFQJL000049,DingYaXinYu,,ZBBGXSZW000377,P17021,,,,9500,non_realtime,6:00:00,106.3674,29.7654
|
||||
GSD2302160925534D50CD17540D_E,demand,GSD2302160925534D50CD17540D_E,ZiYunTai,,ZBBGXSZW000377,P17021,,,,9529,non_realtime,6:00:00,106.387,29.7884
|
||||
GSD2302160925534D50CD17540D_E,demand,GSD2302160925534D50CD17540D_End,ZiYunTai,,ZBBGXSZW000377,P17021,,,,9529,non_realtime,6:00:00,106.387,29.7884
|
||||
ZBBDFQJL000050,demand,ZBBDFQJL000050,XieMaGuangChang,,ZBBGXSZW000377,P17021,,,,9477,non_realtime,6:00:00,106.3652,29.7657
|
||||
GSD2307192058578BCE265C6EA8,demand,GSD2307192058578BCE265C6EA8,YongJinFu,,ZBBGXSZW000377,P17021,,,,9534,non_realtime,6:00:00,106.3831,29.781
|
||||
ZBBDFQJL000028,demand,ZBBDFQJL000028,PanXiMingDu,,P16504,,,,,9514,non_realtime,6:00:00,106.4223,29.821
|
||||
GSD230113095617514F4A2066D7_E,demand,GSD230113095617514F4A2066D7_E,WanKeJinYuHuaFuGaoCeng,,P16504,,,,,9528,non_realtime,6:00:00,106.4228,29.8113
|
||||
GSD230113095617514F4A2066D7_E,demand,GSD230113095617514F4A2066D7_End,WanKeJinYuHuaFuGaoCeng,,P16504,,,,,9528,non_realtime,6:00:00,106.4228,29.8113
|
||||
ZBBDFQJL000014,demand,ZBBDFQJL000014,KeJiXiao,,P16504,,,,,9503,non_realtime,6:00:00,106.4332,29.8258
|
||||
ZBBDFQJL000016,demand,ZBBDFQJL000016,LuGouQiao,,P16504,,,,,9527,non_realtime,6:00:00,106.437,29.8285
|
||||
ZBBDFQJL000057,demand,ZBBDFQJL000057,LongJiangHuaYuan,,P16504,,,,,9495,non_realtime,6:00:00,106.4316,29.8309
|
||||
@@ -74,7 +74,7 @@ ZBBDFQJL000005,demand,ZBBDFQJL000005,TaiJiBinJiangYiQi,,P16504,,,,,9484,non_real
|
||||
ZBBDFQJL000006,demand,ZBBDFQJL000006,TianQiHuaYuan,,P16504,,,,,9523,non_realtime,6:00:00,106.4361,29.8189
|
||||
ZBBDFQJL000004,demand,ZBBDFQJL000004,TaiJiBinJiangErQi,,P16504,,,,,9515,non_realtime,6:00:00,106.4443,29.8297
|
||||
ZBBDFQJL000039,demand,ZBBDFQJL000039,122Zhong,,P16504,,,,,9516,non_realtime,6:00:00,106.4284,29.8109
|
||||
GSD230112144241F42EF6065148_E,demand,GSD230112144241F42EF6065148_E,WanKeJinYuHuaFuYangFang,,P16504,,,,,9530,non_realtime,6:00:00,106.4263,29.8154
|
||||
GSD230112144241F42EF6065148_E,demand,GSD230112144241F42EF6065148_End,WanKeJinYuHuaFuYangFang,,P16504,,,,,9530,non_realtime,6:00:00,106.4263,29.8154
|
||||
J06194,pressure,J06194,,,,,,,,2514,realtime,0:01:00,106.4195183,29.83782213
|
||||
J06190,pressure,J06190,,,,,,,,2510,realtime,0:01:00,106.4194644,29.83781489
|
||||
ZBBDFQJL000025_p,pressure,ZBBDFQJL000025,,,,,,,,9536,non_realtime,6:00:00,106.4120554,29.82753087
|
||||
|
||||
|
557
sensor_placement.py
Normal file
557
sensor_placement.py
Normal file
@@ -0,0 +1,557 @@
|
||||
# 改进灵敏度法
|
||||
import networkx
|
||||
import numpy as np
|
||||
import pandas
|
||||
import wntr
|
||||
import pandas as pd
|
||||
import copy
|
||||
import matplotlib.pyplot as plt
|
||||
import networkx as nx
|
||||
from sklearn.cluster import KMeans
|
||||
from wntr.epanet.toolkit import EpanetException
|
||||
from numpy.linalg import slogdet
|
||||
import random
|
||||
from tjnetwork import *
|
||||
|
||||
|
||||
# 2025/03/12
|
||||
# Step1: 获取节点坐标
|
||||
def getCoor(wn: wntr.network.WaterNetworkModel) -> pandas.DataFrame:
|
||||
"""
|
||||
获取管网模型的节点坐标
|
||||
:param wn: 由wntr生成的模型
|
||||
:return: 节点坐标
|
||||
"""
|
||||
# site: pandas.Series
|
||||
# index：节点名称（wn.node_name_list）
|
||||
# values：每个节点的坐标，格式为 tuple（如 (x, y) 或 (x, y, z)）
|
||||
site = wn.query_node_attribute('coordinates')
|
||||
# Coor: pandas.Series
|
||||
# index：与site相同（节点名称）。
|
||||
# values：坐标转换为numpy.ndarray（如array([10.5, 20.3])）
|
||||
Coor = site.apply(lambda x: np.array(x)) # 将节点坐标转换为numpy数组
|
||||
# x, y: list[float]
|
||||
x = [] # 存储所有节点的 x 坐标
|
||||
y = [] # 存储所有节点的 y 坐标
|
||||
for i in range(0, len(Coor)):
|
||||
x.append(Coor.values[i][0]) # 将 x 坐标存入 x 列表。
|
||||
y.append(Coor.values[i][1]) # 将 y 坐标存入 y 列表
|
||||
# xy: dict[str, list], x、y 坐标的字典
|
||||
xy = {'x': x, 'y': y}
|
||||
# Coor_node: pandas.DataFrame, 存储节点 x, y 坐标的 DataFrame
|
||||
Coor_node = pd.DataFrame(xy, index=wn.node_name_list, columns=['x', 'y'])
|
||||
return Coor_node
|
||||
|
||||
|
||||
# 2025/03/12
|
||||
# Step2: KMeans 聚类
|
||||
# 将节点用kmeans根据坐标分为k组，存入字典g
|
||||
def kgroup(coor: pandas.DataFrame, knum: int) -> dict[int, list[str]]:
|
||||
"""
|
||||
使用KMeans聚类，将节点坐标分组
|
||||
:param coor: 存储所有节点的坐标数据
|
||||
:param knum: 需要分成的聚类数
|
||||
:return: 聚类结果字典
|
||||
"""
|
||||
g = {}
|
||||
# estimator: sklearn.cluster.KMeans,KMeans 聚类模型
|
||||
estimator = KMeans(n_clusters=knum)
|
||||
estimator.fit(coor)
|
||||
# label_pred: numpy.ndarray（int）,每个点的类别标签
|
||||
label_pred = estimator.labels_
|
||||
for i in range(0, knum):
|
||||
g[i] = coor[label_pred == i].index.tolist()
|
||||
return g
|
||||
|
||||
|
||||
# 2025/03/12
|
||||
# Step3: wn_func类，水力计算
|
||||
# wn_func 主要用于计算：
|
||||
# 水力距离（hydraulic length）：即节点之间的水力阻力。
|
||||
# 灵敏度分析（sensitivity analysis）：用于优化测压点的布置。
|
||||
# 一些与水力相关的函数，包括 CtoS：求水力距离,stafun：求状态函数F
|
||||
# # diff：求F对P的导数，返回灵敏度矩阵A
|
||||
# # sensitivity：返回灵敏度和总灵敏度
|
||||
class wn_func(object):
|
||||
|
||||
# Step3.1: 初始化
|
||||
def __init__(self, wn: wntr.network.WaterNetworkModel):
|
||||
"""
|
||||
获取管网模型信息
|
||||
:param wn: 由wntr生成的模型
|
||||
"""
|
||||
# self.results: wntr.sim.results.SimulationResults,仿真结果，包含压力、流量、水头等数据
|
||||
self.results = wntr.sim.EpanetSimulator(wn).run_sim() # 存储运行结果
|
||||
self.wn = wn
|
||||
# self.q：pandas.DataFrame,管道流量，索引为时间步长，列为管道名称
|
||||
self.q = self.results.link['flowrate']
|
||||
# ReservoirIndex / Tankindex: list[str],水库 / 水箱节点名称列表
|
||||
ReservoirIndex = wn.reservoir_name_list
|
||||
Tankindex = wn.tank_name_list
|
||||
# 删除水库节点，删除与直接水库相连的虚拟管道
|
||||
# self.pipes: list[str],所有管道的名称
|
||||
self.pipes = wn.pipe_name_list
|
||||
# self.nodes: list[str],所有节点的名称
|
||||
self.nodes = wn.node_name_list
|
||||
# self.coordinates：pandas.Series,节点坐标，索引为节点名，值为 (x, y) 坐标的 tuple
|
||||
self.coordinates = wn.query_node_attribute('coordinates')
|
||||
# allpumps / allvalves: list[str],所有泵/阀门名称列表
|
||||
allpumps = wn.pump_name_list
|
||||
allvalves = wn.valve_name_list
|
||||
# pumpstnode / pumpednode / valvestnode / valveednode: list[str],存储泵和阀门 起终点节点的名称
|
||||
pumpstnode = []
|
||||
pumpednode = []
|
||||
valvestnode = []
|
||||
valveednode = []
|
||||
# Reservoirpipe / Reservoirednode: list[str],记录与水库相关的管道和节点
|
||||
Reservoirpipe = []
|
||||
Reservoirednode = []
|
||||
for pump in allpumps:
|
||||
pumpstnode.append(wn.links[pump].start_node.name)
|
||||
pumpednode.append(wn.links[pump].end_node.name)
|
||||
for valve in allvalves:
|
||||
valvestnode.append(wn.links[valve].start_node.name)
|
||||
valveednode.append(wn.links[valve].end_node.name)
|
||||
for pipe in self.pipes:
|
||||
if wn.links[pipe].start_node.name in ReservoirIndex:
|
||||
Reservoirpipe.append(pipe)
|
||||
Reservoirednode.append(wn.links[pipe].end_node.name)
|
||||
if wn.links[pipe].start_node.name in Tankindex:
|
||||
Reservoirpipe.append(pipe)
|
||||
Reservoirednode.append(wn.links[pipe].end_node.name)
|
||||
if wn.links[pipe].end_node.name in Tankindex:
|
||||
Reservoirpipe.append(pipe)
|
||||
Reservoirednode.append(wn.links[pipe].start_node.name)
|
||||
# 泵的起终点、tank、reservoir
|
||||
# self.delnodes: list[str],需要删除的节点（包括水库、泵、阀门连接的节点）
|
||||
self.delnodes = list(
|
||||
set(ReservoirIndex).union(Tankindex, pumpstnode, pumpednode, valvestnode, valveednode, Reservoirednode))
|
||||
# 泵、起终点为tank、reservoir的管道
|
||||
# self.delpipes: list[str],需要删除的管道（包括水库、泵、阀门连接的管道）
|
||||
self.delpipes = list(set(wn.pump_name_list).union(wn.valve_name_list).union(Reservoirpipe))
|
||||
self.pipes = [pipe for pipe in wn.pipe_name_list if pipe not in self.delpipes]
|
||||
# self.L: list[float],所有管道的长度（以米为单位）
|
||||
self.L = wn.query_link_attribute('length')[self.pipes].tolist()
|
||||
self.n = len(self.nodes)
|
||||
self.m = len(self.pipes)
|
||||
# self.unit_headloss: list[float],单位水头损失（headloss 数据的第一行，单位：米/km）
|
||||
self.unit_headloss = self.results.link['headloss'].iloc[0, :].tolist()
|
||||
##
|
||||
self.delnodes1 = list(set(ReservoirIndex).union(Tankindex))
|
||||
|
||||
# Step3.2: 计算水力距离
|
||||
def CtoS(self):
|
||||
"""
|
||||
计算水力距离矩阵
|
||||
:return:
|
||||
"""
|
||||
# 水力距离：当行索引对应的节点为控制点时，列索引对应的节点距离控制点的（路径*水头损失）的最小值
|
||||
# nodes：list[str]（节点名称）
|
||||
nodes = copy.deepcopy(self.nodes)
|
||||
# pipes：list[str]（管道名称）
|
||||
pipes = self.pipes
|
||||
wn = self.wn
|
||||
# n / m：int（节点数 / 管道数）
|
||||
n = self.n
|
||||
m = self.m
|
||||
s1 = [0] * m
|
||||
q = self.q
|
||||
L = self.L
|
||||
# H1：pandas.DataFrame,水头数据，索引为时间步长，列为节点名
|
||||
H1 = self.results.node['head'].T
|
||||
# hh：list[float],计算管道两端水头之差
|
||||
hh = []
|
||||
# 水头损失
|
||||
for p in pipes:
|
||||
h1 = self.wn.links[p].start_node.name
|
||||
h1 = H1.loc[str(h1)]
|
||||
h2 = self.wn.links[p].end_node.name
|
||||
h2 = H1.loc[str(h2)]
|
||||
hh.append(abs(h1 - h2))
|
||||
hh = np.array(hh)
|
||||
# headloss：pandas.DataFrame,管道水头损失矩阵
|
||||
headloss = pd.DataFrame(hh, index=pipes).T
|
||||
# s1:管道阻力系数，s2：将管道阻力系数与管道的起始节点和终止节点对应
|
||||
hf = pd.DataFrame(np.array([0] * (n ** 2)).reshape(n, n), index=nodes, columns=nodes, dtype=float)
|
||||
weightL = pd.DataFrame(np.array([0] * (n ** 2)).reshape(n, n), index=nodes, columns=nodes, dtype=float)
|
||||
# s2为对应管道起始节点与终止节点的粗糙度系数矩阵，index代表起始节点，columns代表终止节点
|
||||
G = nx.DiGraph()
|
||||
for i in range(0, m):
|
||||
pipe = pipes[i]
|
||||
a = wn.links[pipe].start_node.name
|
||||
b = wn.links[pipe].end_node.name
|
||||
if q.loc[0, pipe] > 0:
|
||||
hf.loc[a, b] = headloss.loc[0, pipe]
|
||||
weightL.loc[a, b] = headloss.loc[0, pipe] * L[i]
|
||||
G.add_weighted_edges_from([(a, b, weightL.loc[a, b])])
|
||||
|
||||
else:
|
||||
hf.loc[b, a] = headloss.loc[0, pipe]
|
||||
weightL.loc[b, a] = headloss.loc[0, pipe] * L[i]
|
||||
G.add_weighted_edges_from([(b, a, weightL.loc[b, a])])
|
||||
|
||||
hydraulicL = pd.DataFrame(np.array([0] * (n ** 2)).reshape(n, n), index=nodes, columns=nodes, dtype=float)
|
||||
|
||||
for a in nodes:
|
||||
if a in G.nodes:
|
||||
d = nx.shortest_path_length(G, source=a, weight='weight')
|
||||
for b in list(d.keys()):
|
||||
hydraulicL.loc[a, b] = d[b]
|
||||
|
||||
hydraulicL = hydraulicL.drop(self.delnodes)
|
||||
hydraulicL = hydraulicL.drop(self.delnodes, axis=1)
|
||||
|
||||
# 求加权水力距离
|
||||
return hydraulicL, G
|
||||
|
||||
# Step3.3: 计算灵敏度矩阵
|
||||
# 获取关系矩阵
|
||||
def get_Conn(self):
|
||||
"""
|
||||
计算管网连接关系矩阵
|
||||
:return:
|
||||
"""
|
||||
m = self.wn.num_links
|
||||
n = self.wn.num_nodes
|
||||
p = self.wn.num_pumps
|
||||
v = self.wn.num_valves
|
||||
|
||||
self.nonjunc_index = []
|
||||
self.non_link_index = []
|
||||
for r in self.wn.reservoirs():
|
||||
self.nonjunc_index.append(r[0])
|
||||
for t in self.wn.tanks():
|
||||
self.nonjunc_index.append(t[0])
|
||||
# Conn：numpy.matrix，节点-管道连接矩阵，起点 -1，终点 1
|
||||
Conn = np.mat(np.zeros([n, m - p - v])) # 节点和管道的关系矩阵，行为节点，列为管道，起点为-1，终点为1
|
||||
# NConn：numpy.matrix，节点-节点连接矩阵，有管道相连的地方设为 1
|
||||
NConn = np.mat(np.zeros([n, n])) # 节点之间的关系，之间有管道为1，反之为0
|
||||
# pipes：list[str]，去除泵和阀门的管道列表
|
||||
pipes = [pipe for pipe in self.wn.pipes() if pipe not in self.wn.pumps() and pipe not in self.wn.valves()]
|
||||
for pipe_name, pipe in pipes:
|
||||
start = self.wn.node_name_list.index(pipe.start_node_name)
|
||||
end = self.wn.node_name_list.index(pipe.end_node_name)
|
||||
p_index = self.wn.link_name_list.index(pipe_name)
|
||||
Conn[start, p_index] = -1
|
||||
Conn[end, p_index] = 1
|
||||
NConn[start, end] = 1
|
||||
NConn[end, start] = 1
|
||||
self.A = Conn
|
||||
link_name_list = [link for link in self.wn.link_name_list if
|
||||
link not in self.wn.pump_name_list and link not in self.wn.valve_name_list]
|
||||
self.A2 = pd.DataFrame(self.A, index=self.wn.node_name_list, columns=link_name_list)
|
||||
self.A2 = self.A2.drop(self.delnodes)
|
||||
for pipe in self.delpipes:
|
||||
if pipe not in self.wn.pump_name_list and pipe not in self.wn.valve_name_list:
|
||||
self.A2 = self.A2.drop(columns=pipe)
|
||||
self.junc_list = self.A2.index
|
||||
self.A2 = np.mat(self.A2) # 节点管道关系
|
||||
self.A3 = NConn
|
||||
|
||||
def Jaco(self, hL: pandas.DataFrame):
|
||||
"""
|
||||
计算灵敏度矩阵（节点压力对粗糙度变化的响应）
|
||||
:param hL: 水力距离矩阵
|
||||
:return:
|
||||
"""
|
||||
# global result
|
||||
# A：numpy.matrix, 节点-管道关系矩阵
|
||||
A = self.A2
|
||||
wn = self.wn
|
||||
|
||||
try:
|
||||
result = wntr.sim.EpanetSimulator(wn).run_sim()
|
||||
except EpanetException:
|
||||
pass
|
||||
finally:
|
||||
h = result.link['headloss'][self.pipes].values[0]
|
||||
q = result.link['flowrate'][self.pipes].values[0]
|
||||
l = self.wn.query_link_attribute('length')[self.pipes]
|
||||
C = self.wn.query_link_attribute('roughness')[self.pipes]
|
||||
# headloss：numpy.ndarray,水头损失数组
|
||||
headloss = np.array(h)
|
||||
# 调整流量方向
|
||||
for i in range(0, len(q)):
|
||||
if q[i] < 0:
|
||||
A[:, i] = -A[:, i]
|
||||
# q：numpy.ndarray,流量数组
|
||||
q = np.abs(q)
|
||||
# 两个灵敏度矩阵
|
||||
# B / S：numpy.matrix,灵敏度计算的中间矩阵
|
||||
B = np.mat(np.diag(q / ((1.852 * headloss) + 1e-10)))
|
||||
S = np.mat(np.diag(q / C))
|
||||
# X：numpy.matrix, 灵敏度矩阵
|
||||
X = A * B * A.T
|
||||
try:
|
||||
det = np.linalg.det(X)
|
||||
except RuntimeError as e:
|
||||
sign, logdet = slogdet(X) # 防止溢出
|
||||
det = sign * np.exp(logdet)
|
||||
if det != 0:
|
||||
J_H_Cw = X.I * A * S
|
||||
# J_H_Q = -X.I
|
||||
J_q_Cw = S - B * A.T * X.I * A * S # 去掉了delnodes和delpipes
|
||||
# J_q_Q = B * A.T * X.I
|
||||
else: # 当X不可逆
|
||||
J_H_Cw = np.linalg.pinv(X) @ A @ S
|
||||
# J_H_Q = -np.linalg.pinv(X)
|
||||
J_q_Cw = S - B * A.T * np.linalg.pinv(X) * A * S
|
||||
# J_q_Q = B * A.T * np.linalg.pinv(X)
|
||||
|
||||
Sen_pressure = []
|
||||
S_pressure = np.abs(J_H_Cw).sum(axis=1).tolist() # 修改为绝对值
|
||||
for ss in S_pressure:
|
||||
Sen_pressure.append(ss[0])
|
||||
# 求总灵敏度
|
||||
SS_pressure = copy.deepcopy(hL)
|
||||
for i in range(0, len(Sen_pressure)):
|
||||
SS_pressure.iloc[i, :] = SS_pressure.iloc[i, :] * Sen_pressure[i]
|
||||
SS = copy.deepcopy(hL)
|
||||
for i in range(0, len(Sen_pressure)):
|
||||
SS.iloc[i, :] = SS.iloc[i, :] * Sen_pressure[i]
|
||||
# SS[i,j]:节点nodes[i]的灵敏度*该节点到nodes[j]的水力距离
|
||||
return SS
|
||||
|
||||
|
||||
# 2025/03/12
|
||||
# Step4: 传感器布置优化
|
||||
# Sensorplacement
|
||||
# weight：分配权重
|
||||
# sensor：传感器布置的位置
|
||||
class Sensorplacement(wn_func):
|
||||
"""
|
||||
Sensorplacement 类继承了 wn_func 类，并且用于计算和优化传感器布置的位置。
|
||||
"""
|
||||
def __init__(self, wn: wntr.network.WaterNetworkModel, sensornum: int):
|
||||
"""
|
||||
|
||||
:param wn: 由wntr生成的模型
|
||||
:param sensornum: 传感器的数量
|
||||
"""
|
||||
wn_func.__init__(self, wn)
|
||||
self.sensornum = sensornum
|
||||
|
||||
# 1.某个节点到所有节点的加权距离之和
|
||||
# 2.某个节点到该组内所有节点的加权距离之和
|
||||
def sensor(self, SS: pandas.DataFrame, G: networkx.Graph, group: dict[int, list[str]]):
|
||||
"""
|
||||
sensor 方法是用来根据灵敏度矩阵 SS 和加权图 G 来确定传感器布置位置的
|
||||
:param SS: 灵敏度矩阵，每个节点的行和列代表不同节点，矩阵元素表示节点间的灵敏度。SS.iloc[i, :] 表示第 i 行对应节点 i 到所有其他节点的灵敏度
|
||||
:param G: 加权图，表示管网的拓扑结构，每个节点通过管道连接。图的边的权重通常是根据水力距离或者流量等计算的
|
||||
:param group: 节点分组，字典的键是分组编号，值是该组的节点名称列表
|
||||
:return:
|
||||
"""
|
||||
# 传感器布置个数以及位置
|
||||
# W = self.weight()
|
||||
n = self.n - len(self.delnodes)
|
||||
nodes = copy.deepcopy(self.nodes)
|
||||
for node in self.delnodes:
|
||||
nodes.remove(node)
|
||||
# sumSS：list[float],每个节点到其他节点的灵敏度之和。SS.iloc[i, :] 返回第 i 个节点与所有其他节点的灵敏度值，sum(SS.iloc[i, :]) 计算这些灵敏度值的总和。
|
||||
sumSS = []
|
||||
for i in range(0, n):
|
||||
sumSS.append(sum(SS.iloc[i, :]))
|
||||
# 一个整数范围，表示每个节点的索引，用作sumSS_ DataFrame的索引
|
||||
indices = range(0, n)
|
||||
# sumSS_：pandas.DataFrame,将 sumSS 转换成 DataFrame 格式，并且将节点的总灵敏度保存到 CSV 文件 sumSS_data.csv 中
|
||||
sumSS_ = pd.DataFrame(np.array(sumSS), index=indices)
|
||||
sumSS_.to_csv('sumSS_data.csv') # 存储节点总灵敏度
|
||||
# sumSS：pandas.DataFrame,sumSS 被转换为 DataFrame 类型，并且按总灵敏度（即灵敏度之和）降序排列。此时，sumSS 是按节点的灵敏度之和排序的 DataFrame
|
||||
sumSS = pd.DataFrame(np.array(sumSS), index=nodes)
|
||||
sumSS = sumSS.sort_values(by=[0], ascending=[False])
|
||||
# sensorindex：list[str],用于存储根据灵敏度排序选出的传感器位置的节点名称,存储根据总灵敏度排序的节点列表，用于传感器布置
|
||||
sensorindex = []
|
||||
# sensorindex_2：list[str],用于存储每组内根据灵敏度排序选出的传感器位置的节点名称,存储每个组内根据灵敏度排序选择的传感器节点
|
||||
sensorindex_2 = []
|
||||
# group_S：dict[int, pandas.DataFrame],存储每个组内的灵敏度矩阵
|
||||
group_S = {}
|
||||
# group_sumSS：dict[int, list[float]],存储每个组内节点的总灵敏度，值为每个组内节点灵敏度之和的列表
|
||||
group_sumSS = {}
|
||||
for i in range(0, len(group)):
|
||||
for node in self.delnodes:
|
||||
# 这里的group[i]是每个组的节点列表，代码首先去除已经被标记为删除的节点self.delnodes
|
||||
if node in group[i]:
|
||||
group[i].remove(node)
|
||||
group_S[i] = SS.loc[group[i], group[i]]
|
||||
# 对每个组内的节点，计算组内节点的总灵敏度（group_sumSS[i]）。它将每个组内节点的灵敏度值相加，并且按灵敏度降序排序
|
||||
group_sumSS[i] = []
|
||||
for j in range(0, len(group[i])):
|
||||
group_sumSS[i].append(sum(group_S[i].iloc[j, :]))
|
||||
group_sumSS[i] = pd.DataFrame(np.array(group_sumSS[i]), index=group[i])
|
||||
group_sumSS[i] = group_sumSS[i].sort_values(by=[0], ascending=[False])
|
||||
pass
|
||||
|
||||
# 1.选sumSS最大的节点，然后把这个节点所在的那个组删掉，就可以不再从这个组选点。再重新排序选sumSS最大的;
|
||||
# 2.在每组内选group_sumSS最大的节点
|
||||
# 在这个循环中，首先选择灵敏度最高的节点Smaxnode并添加到sensorindex。然后根据灵敏度排序，删除已选的节点并继续选择下一个灵敏度最大的节点。这个过程用于选择传感器的位置
|
||||
sensornum = self.sensornum
|
||||
for i in range(0, sensornum):
|
||||
# Smaxnode：str,最大灵敏度节点，sumSS.index[0] 表示灵敏度最高的节点
|
||||
Smaxnode = sumSS.index[0]
|
||||
sensorindex.append(Smaxnode)
|
||||
sensorindex_2.append(group_sumSS[i].index[0])
|
||||
|
||||
for key, value in group.items():
|
||||
if Smaxnode in value:
|
||||
sumSS = sumSS.drop(index=group[key])
|
||||
continue
|
||||
|
||||
sumSS = sumSS.sort_values(by=[0], ascending=[False])
|
||||
|
||||
return sensorindex, sensorindex_2
|
||||
|
||||
|
||||
# 2025/03/13
|
||||
def get_sensor_coord(name: str, sensor_num: int) -> dict[str, float]:
|
||||
"""
|
||||
获取布置测压点的坐标,初始测压点布置根据灵敏度来布置，计算初始情况下的校准过程的error
|
||||
:param name: 数据库名称
|
||||
:param sensor_num: 测压点数目
|
||||
:return: 测压点坐标字典
|
||||
"""
|
||||
# inp_file_real：str,输入文件名，表示原始水力模型文件的路径,该文件格式为 EPANET 输入文件（.inp），包含管网的结构信息、节点、管道、泵等数据
|
||||
inp_file_real = f'./db_inp/{name}.db.inp'
|
||||
# sensornum：int,需要布置的传感器数量
|
||||
# sensornum = sensor_num
|
||||
# wn_real：wntr.network.WaterNetworkModel,加载 EPANET 水力模型
|
||||
wn_real = wntr.network.WaterNetworkModel(inp_file_real) # 真实粗糙度的原始管网
|
||||
# sim_real：wntr.sim.EpanetSimulator,创建一个水力仿真器对象
|
||||
sim_real = wntr.sim.EpanetSimulator(wn_real)
|
||||
# results_real：wntr.sim.results.SimulationResults,运行仿真并返回结果
|
||||
results_real = sim_real.run_sim()
|
||||
|
||||
# real_C：list[float]，包含所有管道粗糙度的列表
|
||||
real_C = wn_real.query_link_attribute('roughness').tolist()
|
||||
# wn_fun1：wn_func（继承自 object），创建 wn_func 类的实例，传入 wn_real 水力模型对象。wn_func 用于计算管网相关的水力属性，比如水力距离、灵敏度等
|
||||
wn_fun1 = wn_func(wn_real)
|
||||
# nodes：list[str]，管网的节点名称列表
|
||||
nodes = wn_fun1.nodes
|
||||
# delnodes：list[str]，被删除的节点（如水库、泵、阀门连接的节点等）
|
||||
delnodes = wn_fun1.delnodes
|
||||
# Coor_node：pandas.DataFrame
|
||||
Coor_node = getCoor(wn_real)
|
||||
Coor_node = Coor_node.drop(wn_fun1.delnodes)
|
||||
nodes = [node for node in wn_fun1.nodes if node not in delnodes]
|
||||
# coordinates：pandas.Series，存储所有节点的坐标，类型为 Series，索引为节点名称，值为 (x, y) 坐标对
|
||||
coordinates = wn_fun1.coordinates
|
||||
|
||||
# 随机产生监测点
|
||||
# junctionnum：int，nodes 的长度，表示节点的数量
|
||||
junctionnum = len(nodes)
|
||||
# random_numbers：list[int]，使用 random.sample 随机选择 sensornum（20）个节点的编号。它返回一个不重复的随机编号列表
|
||||
# random_numbers = random.sample(range(junctionnum), sensor_num)
|
||||
# for i in range(sensor_num):
|
||||
# # print(random_numbers[i])
|
||||
|
||||
wn_fun1.get_Conn()
|
||||
# hL：pandas.DataFrame，水力距离矩阵，表示每个节点到其他节点的水力阻力
|
||||
# G：networkx.DiGraph，加权有向图，表示管网的拓扑结构，节点之间的边带有权重
|
||||
hL, G = wn_fun1.CtoS()
|
||||
# SS：pandas.DataFrame，灵敏度矩阵，表示每个节点对管网变化（如粗糙度、流量等）的响应
|
||||
SS = wn_fun1.Jaco(hL)
|
||||
# group：dict[int, list[str]]，使用 kgroup 函数将节点按坐标分成若干组，每组包含的节点数不一定相同。group 是一个字典，键为分组编号，值为节点名列表
|
||||
group = kgroup(Coor_node, sensor_num)
|
||||
# wn_fun：Sensorplacement（继承自wn_func）
|
||||
# 创建Sensorplacement类的实例，传入水力网络模型wn_real和传感器数量sensornum。Sensorplacement用于计算和布置传感器
|
||||
wn_fun = Sensorplacement(wn_real, sensor_num)
|
||||
wn_fun.__dict__.update(wn_fun1.__dict__)
|
||||
# sensorindex：list[str]，初始传感器布置位置的节点名称
|
||||
# sensorindex_2：list[str]，根据分组选择的传感器位置
|
||||
sensorindex, sensorindex_2 = wn_fun.sensor(SS, G, group) # 初始的sensorindex
|
||||
# print(str(sensor_num), "个测压点，测压点位置：", sensorindex)
|
||||
sensor_coord = {}
|
||||
# 重新打开数据库
|
||||
if is_project_open(name=name):
|
||||
close_project(name=name)
|
||||
open_project(name=name)
|
||||
for node_id in sensorindex:
|
||||
sensor_coord[node_id] = get_node_coord(name=name, node_id=node_id)
|
||||
close_project(name=name)
|
||||
# print(sensor_coord)
|
||||
return sensor_coord
|
||||
|
||||
|
||||
if __name__ == '__main__':
|
||||
sensor_coord = get_sensor_coord(name='bb', sensor_num=20)
|
||||
print(sensor_coord)
|
||||
# '''
|
||||
# 初始测压点布置根据灵敏度来布置，计算初始情况下的校准过程的error
|
||||
# '''
|
||||
#
|
||||
# # inp_file_real：str,输入文件名，表示原始水力模型文件的路径,该文件格式为 EPANET 输入文件（.inp），包含管网的结构信息、节点、管道、泵等数据
|
||||
# inp_file_real = './db_inp/bb.db.inp'
|
||||
# # sensornum：int,需要布置的传感器数量
|
||||
# sensornum = 20
|
||||
# # wn_real：wntr.network.WaterNetworkModel,加载 EPANET 水力模型
|
||||
# wn_real = wntr.network.WaterNetworkModel(inp_file_real) # 真实粗糙度的原始管网
|
||||
# # sim_real：wntr.sim.EpanetSimulator,创建一个水力仿真器对象
|
||||
# sim_real = wntr.sim.EpanetSimulator(wn_real)
|
||||
# # results_real：wntr.sim.results.SimulationResults,运行仿真并返回结果
|
||||
# results_real = sim_real.run_sim()
|
||||
#
|
||||
# # real_C：list[float]，包含所有管道粗糙度的列表
|
||||
# real_C = wn_real.query_link_attribute('roughness').tolist()
|
||||
# # wn_fun1：wn_func（继承自 object），创建 wn_func 类的实例，传入 wn_real 水力模型对象。wn_func 用于计算管网相关的水力属性，比如水力距离、灵敏度等
|
||||
# wn_fun1 = wn_func(wn_real)
|
||||
# # nodes：list[str]，管网的节点名称列表
|
||||
# nodes = wn_fun1.nodes
|
||||
# # delnodes：list[str]，被删除的节点（如水库、泵、阀门连接的节点等）
|
||||
# delnodes = wn_fun1.delnodes
|
||||
# # Coor_node：pandas.DataFrame
|
||||
# Coor_node = getCoor(wn_real)
|
||||
# Coor_node = Coor_node.drop(wn_fun1.delnodes)
|
||||
# nodes = [node for node in wn_fun1.nodes if node not in delnodes]
|
||||
# # coordinates：pandas.Series，存储所有节点的坐标，类型为 Series，索引为节点名称，值为 (x, y) 坐标对
|
||||
# coordinates = wn_fun1.coordinates
|
||||
#
|
||||
# # 随机产生监测点
|
||||
# # junctionnum：int，nodes 的长度，表示节点的数量
|
||||
# junctionnum = len(nodes)
|
||||
# # random_numbers：list[int]，使用 random.sample 随机选择 sensornum（20）个节点的编号。它返回一个不重复的随机编号列表
|
||||
# random_numbers = random.sample(range(junctionnum), sensornum)
|
||||
# for i in range(sensornum):
|
||||
# print(random_numbers[i])
|
||||
#
|
||||
# wn_fun1.get_Conn()
|
||||
# # hL：pandas.DataFrame，水力距离矩阵，表示每个节点到其他节点的水力阻力
|
||||
# # G：networkx.DiGraph，加权有向图，表示管网的拓扑结构，节点之间的边带有权重
|
||||
# hL, G = wn_fun1.CtoS()
|
||||
# # SS：pandas.DataFrame，灵敏度矩阵，表示每个节点对管网变化（如粗糙度、流量等）的响应
|
||||
# SS = wn_fun1.Jaco(hL)
|
||||
# # group：dict[int, list[str]]，使用 kgroup 函数将节点按坐标分成若干组，每组包含的节点数不一定相同。group 是一个字典，键为分组编号，值为节点名列表
|
||||
# group = kgroup(Coor_node, sensornum)
|
||||
# # wn_fun：Sensorplacement（继承自wn_func）
|
||||
# # 创建Sensorplacement类的实例，传入水力网络模型wn_real和传感器数量sensornum。Sensorplacement用于计算和布置传感器
|
||||
# wn_fun = Sensorplacement(wn_real, sensornum)
|
||||
# wn_fun.__dict__.update(wn_fun1.__dict__)
|
||||
# # sensorindex：list[str]，初始传感器布置位置的节点名称
|
||||
# # sensorindex_2：list[str]，根据分组选择的传感器位置
|
||||
# sensorindex, sensorindex_2 = wn_fun.sensor(SS, G, group) # 初始的sensorindex
|
||||
# print(str(sensornum), "个测压点，测压点位置：", sensorindex)
|
||||
|
||||
# # 分区画图
|
||||
# colorlist = ['lightpink', 'coral', 'rosybrown', 'olive', 'powderblue', 'lightskyblue', 'steelblue', 'peachpuff','brown','silver','indigo','lime','gold','violet','maroon','navy','teal','magenta','cyan',
|
||||
# 'burlywood', 'tan', 'slategrey', 'thistle', 'lightseagreen', 'lightgreen', 'red','blue','yellow','orange','purple','grey','green','pink','lightblue','beige','chartreuse','turquoise','lavender','fuchsia','coral']
|
||||
# G = wn_real.to_graph()
|
||||
# G = G.to_undirected() # 变为无向图
|
||||
# pos = nx.get_node_attributes(G, 'pos')
|
||||
# pass
|
||||
# for i in range(0, sensornum):
|
||||
# ax = plt.gca()
|
||||
# ax.set_title(inp_file_real + str(sensornum))
|
||||
# nodes = nx.draw_networkx_nodes(G, pos, nodelist=group[i], node_color=colorlist[i], node_size=20)
|
||||
# nodes = nx.draw_networkx_nodes(G, pos,
|
||||
# nodelist=sensorindex_2, node_color='black', node_size=70, node_shape='*'
|
||||
# )
|
||||
# edges = nx.draw_networkx_edges(G, pos)
|
||||
# ax.spines['top'].set_visible(False)
|
||||
# ax.spines['right'].set_visible(False)
|
||||
# ax.spines['bottom'].set_visible(False)
|
||||
# ax.spines['left'].set_visible(False)
|
||||
# plt.savefig(inp_file_real + str(sensornum) + ".png")
|
||||
# plt.show()
|
||||
#
|
||||
# wntr.graphics.plot_network(wn_real, node_attribute=sensorindex_2, node_size=50, node_labels=False,
|
||||
# title=inp_file_real + '_Projetion' + str(sensornum))
|
||||
# plt.savefig(inp_file_real + '_S' + str(sensornum) + ".png")
|
||||
# plt.show()
|
||||
@@ -656,11 +656,15 @@ def run_simulation(name: str, simulation_type: str, modify_pattern_start_time: s
|
||||
# 修改工频泵的pattern
|
||||
fixed_pump_SCADA_data_dict = influxdb_api.query_SCADA_data_by_device_ID_and_time(
|
||||
query_ids_list=list(globals.fixed_pumps_id.values()), query_time=modify_pattern_start_time)
|
||||
# print(fixed_pump_SCADA_data_dict)
|
||||
fixed_pump_dict = {key: fixed_pump_SCADA_data_dict[value] for key, value in globals.fixed_pumps_id.items()}
|
||||
# print(fixed_pump_dict)
|
||||
for fixed_pump_name, value in fixed_pump_dict.items():
|
||||
if value:
|
||||
pump_pattern = get_pattern(name_c, get_pump(name_c, fixed_pump_name)['pattern'])
|
||||
print(pump_pattern)
|
||||
pump_pattern['factors'][modify_index] = float(value)
|
||||
print(pump_pattern['factors'][modify_index])
|
||||
cs = ChangeSet()
|
||||
cs.append(pump_pattern)
|
||||
set_pattern(name_c, cs)
|
||||
@@ -887,14 +891,15 @@ def run_simulation(name: str, simulation_type: str, modify_pattern_start_time: s
|
||||
tmp_file = './temp/simulation.result.out'
|
||||
shutil.copy(f'./temp/{name_c}.db.opt', tmp_file)
|
||||
output = Output(tmp_file)
|
||||
|
||||
node_result = output.node_results()
|
||||
link_result = output.link_results()
|
||||
link_flow = []
|
||||
for link in link_result:
|
||||
link_flow.append(link['result'][-1]['flow'])
|
||||
print(link_flow)
|
||||
num_periods_result = output.times()['num_periods']
|
||||
|
||||
print("simulation_type", simulation_type)
|
||||
print("before store result")
|
||||
|
||||
# print(num_periods_result)
|
||||
# print(node_result)
|
||||
# 存储
|
||||
@@ -903,8 +908,7 @@ def run_simulation(name: str, simulation_type: str, modify_pattern_start_time: s
|
||||
elif simulation_type.upper() == 'EXTENDED':
|
||||
influxdb_api.store_scheme_simulation_result_to_influxdb(node_result, link_result, modify_pattern_start_time,
|
||||
num_periods_result, scheme_Type, scheme_Name)
|
||||
|
||||
print("after store result")
|
||||
influxdb_api.fill_scheme_simulation_result_to_SCADA(scheme_Type=scheme_Type, scheme_Name=scheme_Name)
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
@@ -921,28 +925,29 @@ if __name__ == "__main__":
|
||||
globals.source_outflow_region_patterns, globals.realtime_region_pipe_flow_and_demand_patterns = get_realtime_region_patterns('bb', globals.source_outflow_region_id, globals.realtime_region_pipe_flow_and_demand_id)
|
||||
|
||||
# 打印字典内容以验证
|
||||
print("Reservoirs ID:", globals.reservoirs_id)
|
||||
print("Tanks ID:", globals.tanks_id)
|
||||
print("Fixed Pumps ID:", globals.fixed_pumps_id)
|
||||
print("Variable Pumps ID:", globals.variable_pumps_id)
|
||||
print("Pressure ID:", globals.pressure_id)
|
||||
print("Demand ID:", globals.demand_id)
|
||||
print("Quality ID:", globals.quality_id)
|
||||
print("Source Outflow Pattern ID:", globals.source_outflow_pattern_id)
|
||||
print("Realtime Pipe Flow Pattern ID:", globals.realtime_pipe_flow_pattern_id)
|
||||
print("Pipe Flow Region Patterns:", globals.pipe_flow_region_patterns)
|
||||
print("Source Outflow Region:", region_result)
|
||||
print('Source Outflow Region ID:', globals.source_outflow_region_id)
|
||||
print('Source Outflow Region Patterns:', globals.source_outflow_region_patterns)
|
||||
print("Non Realtime Region Patterns:", globals.non_realtime_region_patterns)
|
||||
print("Realtime Region Pipe Flow And Demand ID:", globals.realtime_region_pipe_flow_and_demand_id)
|
||||
print("Realtime Region Pipe Flow And Demand Patterns:", globals.realtime_region_pipe_flow_and_demand_patterns)
|
||||
# print("Reservoirs ID:", globals.reservoirs_id)
|
||||
# print("Tanks ID:", globals.tanks_id)
|
||||
# print("Fixed Pumps ID:", globals.fixed_pumps_id)
|
||||
# print("Variable Pumps ID:", globals.variable_pumps_id)
|
||||
# print("Pressure ID:", globals.pressure_id)
|
||||
# print("Demand ID:", globals.demand_id)
|
||||
# print("Quality ID:", globals.quality_id)
|
||||
# print("Source Outflow Pattern ID:", globals.source_outflow_pattern_id)
|
||||
# print("Realtime Pipe Flow Pattern ID:", globals.realtime_pipe_flow_pattern_id)
|
||||
# print("Pipe Flow Region Patterns:", globals.pipe_flow_region_patterns)
|
||||
# print("Source Outflow Region:", region_result)
|
||||
# print('Source Outflow Region ID:', globals.source_outflow_region_id)
|
||||
# print('Source Outflow Region Patterns:', globals.source_outflow_region_patterns)
|
||||
# print("Non Realtime Region Patterns:", globals.non_realtime_region_patterns)
|
||||
# print("Realtime Region Pipe Flow And Demand ID:", globals.realtime_region_pipe_flow_and_demand_id)
|
||||
# print("Realtime Region Pipe Flow And Demand Patterns:", globals.realtime_region_pipe_flow_and_demand_patterns)
|
||||
|
||||
# dump_inp(name='bb', inp="sensor_placement.inp", version='2')
|
||||
# 模拟示例1
|
||||
run_simulation(name='bb', simulation_type="realtime", modify_pattern_start_time='2025-02-25T23:45:00+08:00')
|
||||
# 模拟示例2
|
||||
# run_simulation(name='bb', simulation_type="extended", modify_pattern_start_time='2025-02-14T10:30:00+08:00', modify_total_duration=900,
|
||||
# scheme_Type="burst_Analysis", scheme_Name="scheme1")
|
||||
# run_simulation(name='bb', simulation_type="extended", modify_pattern_start_time='2025-03-10T12:00:00+08:00',
|
||||
# modify_total_duration=1800, scheme_Type="burst_Analysis", scheme_Name="scheme1")
|
||||
|
||||
# 查询示例1：query_SCADA_ID_corresponding_info
|
||||
# result = query_SCADA_ID_corresponding_info(name='bb', SCADA_ID='P10755')
|
||||
|
||||
Reference in New Issue
Block a user