一种雷电流幅值累积概率分布曲线拟合计算方法技术

本发明专利技术公开了一种雷电流幅值累积概率分布曲线拟合计算方法，在N个统计散点前提下，以IEEE推荐表达式为原型为原型函数、以最小二乘为准则，结合中值电流法与Levenberg‑Marquardt算法，获得一种雷电流幅值累积概率分布曲线拟合计算方法，使得拟合曲线与统计散点最为接近、残差平方和最小。利用本发明专利技术方法能够拟合得到准确的雷电流幅值累积概率分布曲线，客观反映雷电活动规律和特征，用于评价输电线路和杆塔的防雷性能，以及指导防雷工程设计。

【技术实现步骤摘要】
一种雷电流幅值累积概率分布曲线拟合计算方法
本专利技术涉及雷电防护
，更具体地说，涉及一种雷电流幅值累积概率分布曲线拟合计算方法。
技术介绍
雷电流幅值累积概率分布在宏观上体现雷电活动规律，也是开展雷击计算分析的重要参数之一。获取准确的雷电流幅值累积概率分布，对于认知雷电活动规律和特征、准确评价输电线路和杆塔的防雷性能、指导防雷工程设计至关重要。早期由于雷电数据匮乏，我国利用新杭线磁钢棒记录的雷电流幅值数据，反演得出对数形式的雷电流幅值概率方程参数，并列入我国电力规程。目前在雷击计算分析中，常采用IEEE推荐公式，如公式(8)。公式(8)中，I为雷电流幅值，P(>I)为雷电流幅值大于I的概率，a、b为常数参数，IEEE推荐值为a＝31、b＝2.6。a的物理意义为中值电流，即雷电流幅值大于a的概率为50％，b为表征雷电流幅值累积概率分布曲线陡度的参数，b值越大曲线越陡。IEEE推荐参数是综合上世纪全世界数百次自然雷电观测后拟合的总体结果，雷电活动随时间、空间变化的差异性非常大，不宜直接采用，实际应用中往往根据某个时间段内某个特定地区或线路走廊的雷电地闪监测数据统计拟合得到。随着广域雷电地闪监测系统覆盖全国，且绝大部分雷电探测站运行时间超过10年，已积累了海量雷电地闪监测数据，具备了开展雷电参数时空差异性统计研究的条件。广域雷电地闪监测系统针对特定地区或线路走廊记录每一次地闪雷电流幅值及地闪总次数，由此可按定义计算出大于任意幅值的雷电流占地闪总次数百分比，实际应用中通常每隔数kA选择一个计算点，结果作为雷电流幅值累积概率分布的统计散点。由统计散点得到分布曲线的拟合过程，即曲线拟合，要保证拟合误差尽可能小，依赖于拟合方法的准确有效性。I的取值，理论上是(0,∞)，但自然雷电地闪中200kA以上的雷电流已属罕见，因此雷电流幅值累积概率公式的拟合，本质是利用一个区间如[0,200]或[0,300]等之上的观测统计散点近似代表(0,∞)上的分布曲线。实际观测统计中，a的取值在20～40之间较常见；b的值大于1，较常见是在2～3之间。本专利技术旨在提出一种将雷电流幅值累积概率分布统计散点拟合为近似分布曲线的计算方法。
技术实现思路
本专利技术要解决的技术问题在于，在N个统计散点前提下，以IEEE推荐表达式为原型为原型函数、以最小二乘为准则，结合中值电流法与Levenberg-Marquardt(列文博格-马夸尔特)算法，获得一种雷电流幅值累积概率分布曲线拟合计算方法，使得拟合曲线与统计散点最为接近、残差平方和最小。本专利技术解决其技术问题所采用的技术方案是：构造一种雷电流幅值累积概率分布曲线拟合计算方法，包括以下步骤：步骤1)：将雷电地闪监测数据统计计算I取不同值时对应的P(>I)值所得N个散点按I从小至大排序，记此N个散点分别为(X1,Y1)、(X2,Y2)、....、(XN,YN)，(Xj,Yj)表示第j个散点，j为散点的序号，确定如公式(9)所示的原型函数和如公式(10)所示的最小二乘准则下拟合误差计算式，y(Xj)为函数y(x)在Xj处的值，R表示N个散点的拟合残差平方和，；步骤2)：采用中值电流法得到中值电流拟合值a1*，基于a1*反算出陡度参数的拟合值b1*；步骤3)：构造如公式(11)所示的Levenberg-Marquardt残差向量r(p)与如公式(12)所示的目标函数F(p)，；针对公式(9)的二元参数，对应拟合向量p＝[a,b]T、元数M＝2；残差函数向量r(p)表征在N个统计散点处的残差组成的向量，每个元素rj(p)为M元函数；在最小二乘准则下，使F(p)值最小的向量p*即为最优解，残差平方和亦最小；步骤4)：Levenberg-Marquardt迭代初始化，设置初始点p＝p0、向量计算精度ε、阻尼因子μ、阻尼因子变化倍率ν，其中，p0由中值电流法的结果代入，即p0＝[a1*,b1*]T；步骤5)：残差向量r(p)及Jacobi矩阵J(p)迭代，代入p值计算残差向量r(p)、Jacobi矩阵J(p)，Jacobi矩阵J(p)如公式(13)所示，其中p1表示向量p第1个元素、p2表示向量p第2个元素，；步骤6)：求解增量向量，代入p值计算矩阵S(p)＝JT(p)J(p)+μdiag(JT(p)J(p))，构造增量正规方程S(p)Δp＝-JT(p)r(p)，求解出Δp，其中diag(JT(p)J(p))表示由JT(p)J(p)矩阵主对角元素组成的对角矩阵；步骤7)：精度判断，若满足|Δp|<ε则停止迭代、计算结束，此时p即为p*，跳转至步骤9)，否则跳转至步骤8)；步骤8)：迭代判断，计算F(p+Δp)和F(p)的值并比较大小，若F(p+Δp)<F(p)，则令p＝Δp+p更新p、令μ＝μ/ν减小阻尼因子，并跳转至步骤5)；若F(p+Δp)≥F(p)，则令μ＝μν放大阻尼因子，并跳转至步骤6)；步骤9)：Levenberg-Marquardt法迭代计算结束时，a*＝p*1、b*＝p*2，其中p*1表示向量p*的第1个元素、p*2表示向量p*的第2个元素。优选地，在所述步骤2)中，在N个散点中，通过线性插值寻找中值电流拟合值a1*，基于中值电流拟合值a1*反算陡度参数b的拟合值b1*。优选地，在所述步骤2)中，在N个散点中，若存在Yk＝0.5则直接得出中值电流拟合值a1*＝Xk；否则，在0.5附近找出最近的两点(Xm,Ym)和(Xm+1,Ym+1)，且满足Ym>0.5>Ym+1，过(Xm,Ym)、(Xm+1,Ym+1)两点的直线与y＝0.5交点横坐标即为a1*，由公式(6)计算得出，；然后基于中值电流拟合值a1*反算陡度参数b的拟合值b1*，对公式(2)作变换，按公式(15)计算b的拟合值b1*，。实施本专利技术一种雷电流幅值累积概率分布曲线拟合计算方法，具有以下有益效果：1、在相同条件下，本专利技术的计算结果较线性变换拟合、中值电流法拟合等常见简单方法的拟合结果，残差平方和显著减小，约减少1～2个数量级，拟合曲线与统计散点最为接近，在常见拟合方法的结果中最优。2、本专利技术使用中值电流法的结果作为初始值开始Levenberg-Marquardt迭代，初始值计算简单易实现，且避免了选择初始值的盲目性。3、在相同条件下，本专利技术由于使用中值电流法的结果作为初始值开始Levenberg-Marquardt迭代，较选择IEEE推荐值作初始值开始Levenberg-Marquardt迭代，减少了迭代步数，减少约26％的运算时间。4、本专利技术使用Levenberg-Marquardt迭代过程中，迭代判断使用了简化的评价方法，仅根据迭代后目标函数值大小变化决定是否接受本轮迭代计算结果，这种评价方法简单易实现；5、本专利技术以IEEE推荐表达式为原型，保证计算结果的通用性。附图说明下面将结合附图及实施例对本专利技术作进一步说明，附图中：图1为本专利技术利用中值电流法求解其拟合值a1*的示意图；图2为一种雷电流幅值累积概率分布曲线拟合计算方法的流程框图；图3为采用本专利技术描述方法拟合2015年全国地闪雷电流幅值累积概率分布散点得到的曲线图；图4为采用本专利技术描述方本文档来自技高网...

【技术保护点】
一种雷电流幅值累积概率分布曲线拟合计算方法，其特征在于，包括以下步骤：步骤1)：将雷电地闪监测数据统计计算I取不同值时对应的P(>I)值所得N个散点按I从小至大排序，记此N个散点分别为(X1,Y1)、(X2,Y2)、....、(XN,YN)，(Xj,Yj)表示第j个散点，j为散点的序号，确定如公式(1)所示的原型函数和如公式(2)所示的最小二乘准则下拟合误差计算式，y(Xj)为函数y(x)在Xj处的值，R表示N个散点的拟合残差平方和，

【技术特征摘要】
1.一种雷电流幅值累积概率分布曲线拟合计算方法，其特征在于，包括以下步骤：步骤1)：将雷电地闪监测数据统计计算I取不同值时对应的P(>I)值所得N个散点按I从小至大排序，记此N个散点分别为(X1,Y1)、(X2,Y2)、....、(XN,YN)，(Xj,Yj)表示第j个散点，j为散点的序号，确定如公式(1)所示的原型函数和如公式(2)所示的最小二乘准则下拟合误差计算式，y(Xj)为函数y(x)在Xj处的值，R表示N个散点的拟合残差平方和，步骤2)：采用中值电流法得到中值电流拟合值a1*，基于a1*反算出陡度参数的拟合值b1*；步骤3)：构造如公式(3)所示的Levenberg-Marquardt残差向量r(p)与如公式(4)所示的目标函数F(p)，；针对公式(1)的二元参数，对应拟合向量p＝[a,b]T、元数M＝2；残差函数向量r(p)表征在N个统计散点处的残差组成的向量，每个元素rj(p)为M元函数；在最小二乘准则下，使F(p)值最小的向量p*即为最优解，残差平方和亦最小；步骤4)：Levenberg-Marquardt迭代初始化，设置初始点p＝p0、向量计算精度ε、阻尼因子μ、阻尼因子变化倍率ν，其中，p0由中值电流法的结果代入，即p0＝[a1*,b1*]T；步骤5)：残差向量r(p)及Jacobi矩阵J(p)迭代，代入p值计算残差向量r(p)、Jacobi矩阵J(p)，Jacobi矩阵J(p)如公式(5)所示，其中p1表示向量p第1个元素、p2表示向量p第2个元素，

【专利技术属性】
技术研发人员：雷梦飞，赵淳，王少华，李特，王剑，姜文东，王佩，王宇，胡洁，苏杰，李涛，任华，吴军，万家伟，
申请(专利权)人：国网浙江省电力有限公司电力科学研究院，国网浙江省电力有限公司，国网电力科学研究院武汉南瑞有限责任公司，国家电网公司，
类型：发明
国别省市：浙江,33