CN105894130A - 一种用于城市供水管网监测点的优化布置方法 - Google Patents

Abstract

本发明公开了一种用于城市供水管网监测点的优化布置方法。本发明分别改进了传统的灵敏度分析和聚类分析方法，并结合该两种方法得到优化的测点分布，使测点布置不仅具有灵敏度好，相对分散，有利于划定异常区域进行精确侦测的特点，同时较现有测点分布更广、更均匀，且关键位置处布置合理，有利于模型校核。

Description

一种用于城市供水管网监测点的优化布置方法

技术领域

本发明属于城市供水管网领域，具体是一种用于城市供水管网监测的测点优化布置方法。

背景技术

供水管网的测点布置是指在供水管网中分布安装一定数量的流量、压力等传感器，用于监测管网运行工况，为管网调度分析提供数据支持。由于我国城市供水管网规模庞大，管线复杂，测流、测压点在管网初建、扩建过程中部署不多，且缺乏系统考虑，导致测量偏颇，供水管网运行状态不能全面、准确的反映，进而影响安全供水调度的科学性分析。

目前，关于供水管网压力测点优化布置研究主要基于灵敏度分析[1]和聚类分析[2]两类。在传统的灵敏度分析中，通过对每个异常节点寻找灵敏度最大的侦测点构建目标函数，无法保证平均灵敏度足够大；而以平均灵敏度作为目标函数，无法保证能寻优到灵敏度最大的节点。且灵敏度标准差方法会使测点分布不均匀，因为节点灵敏度的相对离散，不能保证实际节点在地理上相对分散，不能良好地反映整个管网状况。在聚类分析中，编网法在对大量节点进行聚类时存在效率过慢、分类结果可靠性低等问题。

有鉴于此，本发明针对管网异常事件侦测和实时模型校正的需要，提出了一种城市供水管网监测点优化布置方法。

文献[1]王伟哲，郄志红，刘美侠，吴鑫淼，基于改进遗传算法的供水管网故障监测点布置优化[J]，水力发电学报，2012,31(1):15-19.

文献[2]周书葵，许仕荣，城市供水管网水压监测点优化布置的研究[J]，南华大学学报(自然科学版)，2005(01).

发明内容

本发明针对现有技术的不足，提出了一种最优灵敏度分析和模糊聚类分析相结合的方法，分别改进了传统的灵敏度分析和聚类分析方法，并结合该两种方法得到优化的测点分布，使测点布置不仅具有灵敏度好，相对分散，有利于划定异常区域进行精确侦测的特点，同时较现有测点分布更广、更均匀，且关键位置处布置合理，有利于模型校核。

为实现以上目的，本发明采取以下步骤：

1.获取模糊相似矩阵

(1)灵敏度系数矩阵

令管网模型节点k的流量(也可以是如管道摩阻系数等其它参量)变化，使所有节点i水压发生变化，将所有节点i的水压变化值ΔH_i与节点自身水压的变化值ΔH_k求比值，获得节点的灵敏度值λ_ik，即：

${\mathrm{\λ}}_{ik}=\frac{{\mathrm{\ΔE}}_i}{{\mathrm{\ΔE}}_k}$

据此，可求得所有节点的灵敏度系数，用矩阵形式表示为：

其中，第i行代表节点i发生变化时，对管网所有节点的影响系数向量，第j列代表管网中所有节点发生变化时，对节点j的影响系数向量。

(2)标准化

对灵敏度系数矩阵进行平移标准差变换处理，公式如下：

${R^{\′}}_{ij}=\frac{R_{ij}-\stackrel{\‾}{R_j}}{S_j}(i=1,2,\mathrm{...}n,j=1,2,\mathrm{...}n)$

${\stackrel{\‾}{R}}_j=\frac{1}{n}{\Σ}_{i=1}^n{R}_{ij},S_j=\sqrt{\frac{1}{n}{\Σ}_{i=1}^n{(R_{ij}-\stackrel{\‾}{R_j})}^2}$

其中R_ij是R的第i行第j列的值，是R的第j列的平均值，S_j是R的第j列的标准差，R′_ij为所求标准化矩阵第i行第j列的值。

(3)模糊相似化

采用相似系数法中的相关系数法求解标准化矩阵，公式如下：

${R^{\′\′}}_{ij}=\underset{i=1}{\overset{m}{\Σ}}\underset{j=1}{\overset{m}{\Σ}}\frac{{\Σ}_{k=1}^m\left(\right|{R^{\′}}_{ik}-\stackrel{\‾}{{R^{\′}}_i}\left|\right|{R^{\′}}_{kj}-\stackrel{\‾}{{R^{\′}}_j}\left|\right)}{\sqrt{{\Σ}_{k=1}^m{({R^{\′}}_{ik}-\stackrel{\‾}{{R^{\′}}_i})}^2}\sqrt{{\Σ}_{k=1}^m{({R^{\′}}_{kj}-\stackrel{\‾}{{R^{\′}}_j})}^2}}$

$\begin{array}{cc}{\stackrel{\‾}{R^{\′}}}_i=\frac{1}{m}{\Σ}_{k=1}^m{R^{\′}}_{ik}& {\stackrel{\‾}{R^{\′}}}_j=\frac{1}{m}{\Σ}_{k=1}^m{R^{\′}}_{kj}\end{array}$

其中，是R′的第i行的平均值，是R′的第j列的平均值。

2.基于最优灵敏度确定测点

确定每个异常节点j对管网所有节点感应灵敏度的平均值和最大值，先分别平方再求和，最后开根号，增加测点之间的平均欧氏距离之和来控制测点之间的分散程度。构建目标函数如下：

$\max f=\sqrt{{\mathrm{\ϵ}}_1\underset{j=1}{\overset{M}{\Σ}}({\left(\frac{1}{N}\underset{i=1}{\overset{M}{\Σ}}\left|\frac{{\mathrm{\ΔH}}_i}{{\mathrm{\ΔH}}_j}\right|\right)}^2+{\mathrm{\ϵ}}_2\max {\left(\frac{{\mathrm{\ΔH}}_i}{{\mathrm{\ΔH}}_j}\right)}^2)}+{\mathrm{\ϵ}}_3\frac{1}{M(M-1)/2}\underset{i=1}{\overset{M}{\Σ}}\underset{j=1}{\overset{M}{\Σ}}||{\mathrm{cor}}_i-{\mathrm{cor}}_j||$

上式中，N代表测点的数量，M代表可能的爆管点数,i＝1,2,3,...,N。ε₁、ε₂、ε₃是引入的权重系数，代表该部分在整个目标函数中的重要程度，由反复实验获得较好的给定值。ΔHi/ΔH_j是测点i对可能爆管点j的灵敏度系数，这里使用模糊相似系数矩阵R″，cor_i、cor_j是测点i和测点j的二维地理坐标，可以通过模型文件得到。第一部分不仅保证了平均灵敏度足够大同时保证对最大灵敏度节点的寻优，第二部分可在实际地理空间上控制测点之间的分散程度。通过matlab使用实数编码的遗传算法求解该目标函数，便得到基于最优灵敏度分析的测点分布。

3.基于模糊聚类确定测点

使用Kmeans聚类分析方法实现快速聚类，以计算节点灵敏度到原型的海明距离为目标，通过求取极值得到迭代运算的调整规则。算法通过随机初始聚类中心向量V进行最优分类，使得下式所示的目标函数J最小。

$J={\Σ}_{n=1}^N{\Σ}_{k=1}^K{r}_{nk}{||x_n-u_k||}^2$

u_k是分类中心向量，x_n是观测坐标向量，r_nk是0或1，若n被归类至k组，则r_nk＝1，否则r_nk＝0；J的目标是最小化观测向量与聚类中心之间的距离。r_nk与u_k寻优时，一般是通过迭代运算逐步求出，固定u_k，选择最优的r_nk，再固定r_nk，求最优的u_k，在选择r_nk时，将数据归类到距其最近的中心点可使得目标函数最小，求u_k时，令J对u_k求导等于0，可得到：

$u_k=\frac{{\Σ}_n{r}_{nk}{x}_n}{{\Σ}_n{r}_{nk}}$

u_k是每个分类块数据的平均值，因为每次都能够使J取最小，因此J的最小值可以求得，步骤如下：

(1)设定k个节点的灵敏度作为分类中心的初始点，即对应k个测点中心；

(2)计算所有节点的灵敏度到k个分类中心的距离，按照最近邻法归类；

(3)计算每一类的灵敏度距离平均值，距离最小的点作为新的分类中心；

(4)重复步骤(2)、(3)，在达到迭代次数上限或满足要求时停止，便得到基于kmeans模糊聚类分析的测点分布。

4综合确定测点

根据步骤3得出的测点分布，确定步骤2得出的测点分布在关键出入水口处或重要管段处布置缺失的测点个数K以及对应的区域坐标Ф_(xi,yi)(i＝1，2，...,k)。比较步骤2的测点中距离接近的测点灵敏度大小，把灵敏度相对较小的K的测点移动到Ф_(xi,yi)(i＝1，2，...,k)处，便得到最终测点集。

本发明的有益效果如下：

通过标准-模糊相似化方法改进灵敏度系数矩阵数据之间的相关性，提高灵敏度寻优和模糊聚类的准确度。

为更好地进行管网异常事件侦测，对目标函数灵敏度算法进行改进，即保证灵敏度总和足够大又保证对最大灵敏度节点的寻优，同时引入空间距离离散度因子，在地理上保证优化后的测点相对分散。

使用了Kmeans聚类分析方法实现快速聚类，虽然结果与初始的分类中心有关，但在足够多的迭代次数下，解具有相对稳定性，在关键出入水口处皆有分布。与现有测点分布相比，能为供水管网提供更全面的监测数据。

分析以上两种方法得到的测点分布的优点，绘制出既有利于异常事件侦测又有利于模型校正的最终结果图。

附图说明

图1用于异常事件侦测的测点优化示意图；

图2用于模型校正的测点优化示意图；

图3最终测点分布图。

具体实施方式

为使本发明实现的技术创新点易于理解，下面结合附图和实例，对本发明的实现方式进一步详细叙述，具体步骤如下：

1.获得模糊相似矩阵

(1)灵敏度系数矩阵

S市DMA共有节点5377个(不计入水源)，为选取较好的工况，采用最大工作时，中午12点的DMA的模型数据作为试验数据。节点的需水量示例如表格1所示：

表1节点流量

其中，节点编号按照索引进行排序。对于某些节点的需水量为0，比如索引5375，令其需水量为很小值，0.0000001。运行EPANET模型，进行延时模拟，得到节点水量变化前的节点压力向量R⁰：

$R^0=\[R_{11}^0{R}_{12}^0{R}_{13}^0...R_{1n-k}^0\]$

k是节点的水源数目，代表节点1,2,3...n-k改变需水量前的水头或者是压力值。改变管网中每一个节点的需水量，使其变化20％，运行模型，得到变化后的管网节点压力向量

其中，表示节点i的需水量改变后，所有普通节点压力值大小。在每次模型运行结束后，节点i需水量恢复原值，依次改变每个普通节点需水量，得到的n-k个压力变化行向量就构成(n-k)×(n-k)维矩阵，为简洁表示，以下n-k统一表示成m。通过有限差分处理，获取灵敏度系数矩阵如下：

由于矩阵维数是5377×5377，显示内容太多，这里只给出10×10的矩阵示例：

$\begin{array}{c}R=\\ \begin{array}{cccccccccc}\[& & & & & & & & & \\ 1.00& 0.47& 0.41& 0.41& 0.37& 0.50& 0.42& 0.31& 0.31& 0.30\\ 0.86& 1.00& 0.85& 0.73& 0.64& 0.94& 0.75& 0.51& 0.51& 0.47\\ 0.84& 0.95& 1.00& 0.75& 0.70& 0.83& 0.71& 0.55& 0.55& 0.58\\ 0.58& 0.56& 0.52& 1.00& 0.79& 0.57& 0.49& 0.38& 0.38& 0.40\\ 0.60& 0.58& 0.55& 0.90& 1.00& 0.58& 0.51& 0.43& 0.43& 0.46\\ 0.69& 0.67& 0.60& 0.60& 0.55& 1.00& 0.73& 0.41& 0.41& 0.45\\ 0.47& 0.46& 0.42& 0.41& 0.38& 0.65& 1.00& 0.31& 0.31& 0.34\\ 0.89& 0.89& 0.92& 0.93& 0.95& 0.77& 0.78& 1.00& 1.00& 0.88\\ 0.73& 0.74& 0.75& 0.74& 0.74& 0.72& 0.70& 0.83& 1.00& 0.70\\ 0.68& 0.68& 0.68& 0.68& 0.68& 0.68& 0.68& 0.61& 0.61& 1.00\\ \]& & & & & & & & & \end{array}\end{array}$

(2)标准化

对灵敏度系数矩阵中每个元素进行平移标准差变换处理，公式如下：

${R^{\′}}_{ij}=\frac{R_{ij}-\stackrel{\‾}{R_j}}{S_j},(i=1,2,\mathrm{...}n,j=1,2,\mathrm{...}n)$

$\begin{array}{cc}{\stackrel{\‾}{R}}_j=\frac{1}{n}{\Σ}_{i=1}^n{R}_{ij},& S_j=\sqrt{\frac{1}{n}{\Σ}_{i=1}^n{(R_{ij}-\stackrel{\‾}{R_j})}^2}\end{array}$

其中R_ij是R的第i行第j列的值，是R的第j列的平均值，S_j是R的第j列的标准差，R′_ij为所求标准化矩阵第i行第j列的值。

(3)模糊相似化

针对标准化矩阵，采用相似系数法中的相关系数法计算各个分类对象之间的相似程度，从而建立模糊相似矩阵，公式如下：

${R^{\′\′}}_{ij}=\underset{i=1}{\overset{m}{\Σ}}\underset{j=1}{\overset{m}{\Σ}}\frac{{\Σ}_{k=1}^m\left(\right|{R^{\′}}_{ik}-\stackrel{\‾}{{R^{\′}}_i}\left|\right|{R^{\′}}_{kj}-\stackrel{\‾}{{R^{\′}}_j}\left|\right)}{\sqrt{{\Σ}_{k=1}^m{({R^{\′}}_{ik}-\stackrel{\‾}{{R^{\′}}_i})}^2}\sqrt{{\Σ}_{k=1}^m{({R^{\′}}_{kj}-\stackrel{\‾}{{R^{\′}}_j})}^2}}$

$\begin{array}{cc}{\stackrel{\‾}{R^{\′}}}_i=\frac{1}{m}{\Σ}_{k=1}^m{R^{\′}}_{ik}& {\stackrel{\‾}{R^{\′}}}_j=\frac{1}{m}{\Σ}_{k=1}^m{R^{\′}}_{kj}\end{array}$

其中，是R′的第i行的平均值，是R′的第j列的平均值，R″_ij为所求模糊相似化矩阵第i行第j列的值。

由于矩阵维数是5377×5377，显示内容太多，这里只给出10×10的矩阵示例：

$\begin{array}{c}{R}^{\′\′}=\\ \begin{array}{cccccccccc}\[& & & & & & & & & \\ 1.00& 0.27& 0.34& 0.40& 0.36& 0.34& 0.36& 0.09& 0.16& 0.11\\ 0.27& 1.00& 0.34& 0.39& 0.27& 0.49& 0.49& -0.08& 0.09& -0.24\\ 0.34& 0.34& 1.00& 0.84& 0.86& 0.76& 0.61& 0.47& 0.73& 0.59\\ 0.40& 0.39& 0.84& 1.00& 0.98& 0.81& 0.60& 0.33& 0.49& 0.35\\ 0.36& 0.27& 0.86& 0.98& 1.00& 0.76& 0.51& 0.44& 0.60& 0.45\\ 0.34& 0.49& 0.76& 0.81& 0.76& 1.00& 0.70& 0.21& 0.35& 0.21\\ 0.36& 0.49& 0.61& 0.60& 0.51& 0.70& 1.00& 0.03& 0.18& 0.17\\ 0.09& -0.08& 0.47& 0.33& 0.44& 1.21& 0.03& 1.00& 0.72& 0.29\\ 0.16& 0.09& 0.73& 0.49& 0.60& 0.35& 0.18& 0.72& 1.00& 0.64\\ 0.11& -0.24& 0.59& 0.35& 0.45& 0.21& 0.17& 0.29& 0.64& 1.00\\ \]& & & & & & & & & \end{array}\end{array}$

2.基于最优灵敏度确定测点

设定测点数目M为20个，N是DMA中所有普通节点的个数，即5377个。ΔHi/ΔH_j是测点i对可能爆管点j的灵敏度系数，这里使用模糊相似系数矩阵R″，ε₁、ε₂、ε₃∈(0，1)且ε₁+ε₂+ε₃＝1，具体值根据优化后的测点分布比较原有测点分布，通过动态调整最终确定。cor_i、cor_j是测点i和测点j的二维地理坐标，可以通过EPANET模型文件得到。

目标函数如下：

$\max f=\sqrt{{\mathrm{\ϵ}}_1\underset{j=1}{\overset{M}{\Σ}}({\left(\frac{1}{N}\underset{i=1}{\overset{M}{\Σ}}\left|\frac{{\mathrm{\ΔH}}_i}{{\mathrm{\ΔH}}_j}\right|\right)}^2+{\mathrm{\ϵ}}_2\max {\left(\frac{{\mathrm{\ΔH}}_i}{{\mathrm{\ΔH}}_j}\right)}^2)}+{\mathrm{\ϵ}}_3\frac{1}{M(M-1)/2}\underset{i=1}{\overset{M}{\Σ}}\underset{j=1}{\overset{M}{\Σ}}||{\mathrm{cor}}_i-{\mathrm{cor}}_j||$

通过matlab使用实数编码的遗传算法求解目标函数，设定参数为：种群规模50，遗传代数200，选择概率0.95，变异概率0.05。其中，选择、交叉、变异方法分别使用轮盘赌选择、实数交叉和单点变异方法。上述目标函数是求解最大值问题，通过进化迭代，得出用于异常事件侦测的测点优化布置结果，如图1右侧所示，其中，圆形标记代表优化测点，为方便比较，图1左侧给出了原有测点分布图。

3.基于模糊聚类确定测点

利用以上步骤中得到的模糊相似矩阵结合Kmeans方法进行分类。使用matlab自带的Kmeans工具箱进行求解。

现有的DMA内部测点有14个，但依据现有测点分布，使实时模型定位难度较大，因此，考虑增加测压点，这里设定测点数目M为20，Kmeans函数参数中分类距离设为“cityblock”；设定重复次数为5；设定“emptyaction”为“drop”，即在出现上述情况时丢弃原有分类块，继续迭代求解，通过运算，得出分类结果。

通过分类得到20个分类块后，先统计每一个块中的分类点对应的节点索引，然后使用步骤1中模糊相似化的方法得到模糊相似矩阵。

计算每个节点除自身外的所有其余节点的灵敏度相似系数值之和，作为每个节点与其余节点的相似距离，公式如下所示：

$dis={\Σ}_{j=1,i\≠j}^m{R^{\′\′}}_{i,j}$

cent＝arc(min(dis))

最后取与其他节点相似距离和最小的节点为中心点，也就是优化布置后的测点。Matlab下利用节点坐标信息，绘出20个优化测点(用方形表示)在管网中的分布，如图2右侧所示，为便于与现有测点比较，图2左侧给出了原有测点分布图。

4.综合确定测点

获取图1、图2中的优化测点在EPANET模型中对应的节点坐标信息，以图中1的优化测点分布为蓝本，对比图2中的优化测点分布，确定需要移动到由图2中的优化测点确定的关键出入水口处或重要管段处(图3b已圈出)的测点个数K＝4和坐标索引Ф_(xi,yi)(i＝1，2，3，4)。通过matlab分别计算图1中优化测点间距离最接近的四个区域(图3a已圈出)的灵敏度大小，把灵敏度相对较小的4个测点分配到Ф_(xi,yi)(i＝1，2，3，4)处(分别是左下角、右上角、最右侧和中间靠下)，得到最终的优化测点在管网中的分布，如图3c所示(圆形代表图1中的优化测点，方形代表图2中的优化测点)，为了便于比较，图3a给出用于异常事件侦测的优化测点分布图，图3b给出用于模型校正的优化测点分布图，图3c右侧给出了最终测点分布图。

Claims (1)

1.一种用于城市供水管网监测点的优化布置方法，其特征在于该方法包括以下步骤：

步骤1.获取模糊相似矩阵

(1)建立灵敏度系数矩阵

令管网模型节点k的流量变化，使所有节点i水压发生变化，将所有节点i的水压变化值ΔH_i与节点自身水压的变化值ΔH_k求比值，获得节点的灵敏度值λ_ik，即：

据此，可求得所有节点的灵敏度系数，用矩阵形式R表示为：

其中，第i行代表节点i发生变化时，对管网所有节点的影响系数向量，第j列代表管网中所有节点发生变化时，对节点j的影响系数向量；

(2)标准化

对灵敏度系数矩阵进行平移标准差变换处理，公式如下：

其中R_ij是R的第i行第j列的值，是R的第j列的平均值，S_j是R的第j列的标准差，R′_ij为所求标准化矩阵第i行第j列的值；

(3)模糊相似化

采用相似系数法中的相关系数法求解标准化矩阵，公式如下：

其中，是R′的第i行的平均值，是R′的第j列的平均值；

步骤2.基于最优灵敏度确定测点

确定每个异常节点j对管网所有节点感应灵敏度的平均值和最大值，先分别平方再求和，最后开根号，增加测点之间的平均欧氏距离之和来控制测点之间的分散程度；构建目标函数如下：

上式中，N代表测点的数量，M代表可能的爆管点数；ε₁、ε₂、ε₃是引入的权重系数，代表该部分在整个目标函数中的重要程度,ΔHi/ΔH_j是测点i对可能爆管点j的灵敏度系数；cor_i、cor_j是测点i和测点j的二维地理坐标，通过模型文件得到；通过matlab使用实数编码的遗传算法求解该目标函数，便得到基于最优灵敏度分析的测点分布；

步骤3.基于模糊聚类确定测点

使用Kmeans聚类分析方法实现快速聚类，以计算节点灵敏度到原型的海明距离为目标，通过求取极值得到迭代运算的调整规则；算法通过随机初始聚类中心向量V进行最优分类，使得下式所示的目标函数J最小；

u_k是分类中心向量，x_n是观测坐标向量，r_nk是0或1，若n被归类至k组，则r_nk＝1，否则r_nk＝0；J的目标是最小化观测向量与聚类中心之间的距离；r_nk与u_k寻优时，通过迭代运算逐步求出，固定u_k，选择最优的r_nk，再固定r_nk，求最优的u_k，在选择r_nk时，将数据归类到距其最近的中心点可使得目标函数最小，求u_k时，令J对u_k求导等于0，可得到：

u_k是每个分类块数据的平均值，因为每次都能够使J取最小，因此J的最小值可以求得，步骤如下：

(1)设定k个节点的灵敏度作为分类中心的初始点，即对应k个测点中心；

(2)计算所有节点的灵敏度到k个分类中心的距离，按照最近邻法归类；

(3)计算每一类的灵敏度距离平均值，距离最小的点作为新的分类中心；

(4)重复步骤(2)、(3)，在达到迭代次数上限或满足要求时停止，便得到基于kmeans模糊聚类分析的测点分布；

步骤4.综合确定测点

根据步骤3得出的测点分布，确定步骤2得出的测点分布在关键出入水口处或重要管段处布置缺失的测点个数K以及对应的区域坐标Φ_(xi,yi)；比较步骤2的测点中距离接近的测点灵敏度大小，把灵敏度相对较小的K的测点移动到Φ_(xi,yi)处，便得到最终测点集。