CN108051676A - 一种雷电流幅值累积概率分布曲线拟合计算方法 - Google Patents

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

Info

Publication number
CN108051676A
CN108051676A CN201711328027.3A CN201711328027A CN108051676A CN 108051676 A CN108051676 A CN 108051676A CN 201711328027 A CN201711328027 A CN 201711328027A CN 108051676 A CN108051676 A CN 108051676A
Authority
CN
China
Prior art keywords
fitting
vector
value
formula
lightning
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201711328027.3A
Other languages
English (en)
Other versions
CN108051676B (zh
Inventor
雷梦飞
赵淳
王少华
李特
王剑
姜文东
王佩
王宇
胡洁
苏杰
李涛
任华
吴军
万家伟
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
State Grid Corp of China SGCC
Wuhan NARI Ltd
State Grid Zhejiang Electric Power Co Ltd
Electric Power Research Institute of State Grid Zhejiang Electric Power Co Ltd
Original Assignee
State Grid Corp of China SGCC
Wuhan NARI Ltd
State Grid Zhejiang Electric Power Co Ltd
Electric Power Research Institute of State Grid Zhejiang Electric Power Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by State Grid Corp of China SGCC, Wuhan NARI Ltd, State Grid Zhejiang Electric Power Co Ltd, Electric Power Research Institute of State Grid Zhejiang Electric Power Co Ltd filed Critical State Grid Corp of China SGCC
Priority to CN201711328027.3A priority Critical patent/CN108051676B/zh
Publication of CN108051676A publication Critical patent/CN108051676A/zh
Application granted granted Critical
Publication of CN108051676B publication Critical patent/CN108051676B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R31/00—Arrangements for testing electric properties; Arrangements for locating electric faults; Arrangements for electrical testing characterised by what is being tested not provided for elsewhere

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Steroid Compounds (AREA)
  • Complex Calculations (AREA)

Abstract

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

Description

一种雷电流幅值累积概率分布曲线拟合计算方法
技术领域
本发明涉及雷电防护技术领域，更具体地说，涉及一种雷电流幅值累积概率分布曲线拟合计算方法。
背景技术
雷电流幅值累积概率分布在宏观上体现雷电活动规律，也是开展雷击计算分析的重要参数之一。获取准确的雷电流幅值累积概率分布，对于认知雷电活动规律和特征、准确评价输电线路和杆塔的防雷性能、指导防雷工程设计至关重要。
早期由于雷电数据匮乏，我国利用新杭线磁钢棒记录的雷电流幅值数据，反演得出对数形式的雷电流幅值概率方程参数，并列入我国电力规程。目前在雷击计算分析中，常采用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为采用本发明描述方法和中值电流法、线性变换法对2005～2015年全国地闪雷电流幅值累积概率分布散点进行拟合的误差对比图。
具体实施方式
为了对本发明的技术特征、目的和效果有更加清楚的理解，现对照附图详细说明本发明的具体实施方式。
如图2所示，本发明提供一种雷电流幅值累积概率分布曲线拟合计算方法，包括以下步骤：
步骤1)：确定原型函数如式(16)，确定最小二乘准则下拟合误差计算式如式(17)；其中，R为N个散点的拟合残差平方和，j为散点的序号，(Xj,Yj)为第j个散点，y(Xj)为函数y(x)在Xj处的值，
。
步骤2)：在N个散点中，通过线性插值寻找中值电流拟合值a1*，基于中值电流拟合值a1*反算陡度参数b的拟合值b1*。如图1所示，在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*，由公式(18)计算得出，
；
然后基于中值电流拟合值a1*反算陡度参数b的拟合值b1*，对公式(16)作变换，按公式(19)计算b的拟合值b1*，
。
步骤3)：构造如公式(20)所示的Levenberg-Marquardt残差向量r(p)与如公式(21)所示的目标函数F(p)，
；
针对公式(16)的二元参数，对应拟合向量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。向量计算精度ε、阻尼因子μ、阻尼因子变化倍率ν的取值分别为ε＝1.0×10-5、μ＝0.001、ν＝10。
步骤5)：残差向量r(p)及Jacobi矩阵J(p)迭代，代入p值计算残差向量r(p)、Jacobi矩阵J(p)，Jacobi矩阵J(p)如公式(22)所示，其中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个元素。
以2005～2015年全国范围内广域雷电监测***探测的雷电地闪为例，使用本发明方法进行逐年雷电流幅值累积概率分布曲线拟合。在拟合之前，统计得到以5kA为间隔、雷电流幅值I取0、5、...、300kA处P(>I)的值，每年的雷电流幅值累积概率分布各61个散点(即N＝61)，并按雷电流幅值由小至大排序。本实施例参照本发明方法如图2所示的流程进行2015年全国地闪雷电流幅值累积概率分布散点的拟合计算，包括以下步骤：
步骤1)：散点个数N＝61，确定如公式(16)所示的原型函数，确定如公式(17)所示的拟合残差平方和作为最小二乘准则下拟合误差，用于评判拟合结果优劣。
步骤2)：在61个散点中，通过线性插值得出中值电流a1*＝23.52；基于中值电流a1*和公式(19)，计算陡度参数b的拟合值b1*＝2.81。
步骤3)：构造Levenberg-Marquardt残差向量与目标函数，针对公式(16)所示的二元参数，构造拟合向量p＝[a,b]T、元数M＝2，如公式(20)所示的残差函数向量r(p)，如公式(21)所示的目标函数F(p)；在最小二乘准则下，使F(p)值最小的向量p*即为最优解，拟合误差R亦最小。
步骤4)：Levenberg-Marquardt迭代初始化，设置初始点p＝p0＝[23.52,2.81]T、向量计算精度ε＝1.0×10-5，阻尼因子μ＝0.001、阻尼因子变化倍率ν＝10。
步骤5)：残差函数向量r(p)及Jacobi矩阵迭代，将p值代入公式(20)计算残差向量r(p)、代入公式(22)计算Jacobi矩阵J(p)，其中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、令μ＝μ/ν减小阻尼因子，并跳转至步骤6；若F(p+Δp)≥F(p)，则令μ＝μν放大阻尼因子，并跳转至步骤6)。
步骤9)：迭代计算结束，p*＝[23.11,2.42]T，即a*＝23.11、b*＝2.42。
根据计算结果，可以绘制出2015年全国地闪雷电流幅值累积概率分布曲线，并与散点放在同一个坐标系中对比，如图3所示，同时可以按公式(17)算出拟合误差R＝0.0018。同理，可全部得出2005～2015年逐年全国地闪雷电流幅值累积概率分布结果，并计算出拟合误差。图3为采用本发明描述方法拟合2015年全国地闪雷电流幅值累积概率分布散点得到的曲线图；图4为采用本发明描述方法和中值电流法、线性变换法对2005～2015年全国地闪雷电流幅值累积概率分布散点进行拟合的误差对比图。从图3和图4可以看出，本发明方法的拟合曲线与统计散点极为贴近，拟合误差R最小，且较另外两种拟合计算方法的误差低1～3个数量级。
上面结合附图对本发明的实施例进行了描述，但是本发明并不局限于上述的具体实施方式，上述的具体实施方式仅仅是示意性的，而不是限制性的，本领域的普通技术人员在本发明的启示下，在不脱离本发明宗旨和权利要求所保护的范围情况下，还可做出很多形式，这些均属于本发明的保护之内。

Claims (3)

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个散点的拟合残差平方和，
<mrow> <mi>y</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfrac> <mn>1</mn> <mrow> <mn>1</mn> <mo>+</mo> <msup> <mrow> <mo>(</mo> <mi>x</mi> <mo>/</mo> <mi>a</mi> <mo>)</mo> </mrow> <mi>b</mi> </msup> </mrow> </mfrac> <mo>,</mo> <mrow> <mo>(</mo> <mi>x</mi> <mo>&amp;GreaterEqual;</mo> <mn>0</mn> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <mi>R</mi> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>j</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </munderover> <msup> <mrow> <mo>(</mo> <msub> <mi>Y</mi> <mi>j</mi> </msub> <mo>-</mo> <mi>y</mi> <mo>(</mo> <msub> <mi>X</mi> <mi>j</mi> </msub> <mo>)</mo> <mo>)</mo> </mrow> <mn>2</mn> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> <mo>;</mo> </mrow>
步骤2)：采用中值电流法得到中值电流拟合值a1*，基于a1*反算出陡度参数的拟合值b1*；
步骤3)：构造如公式(3)所示的Levenberg-Marquardt残差向量r(p)与如公式(4)所示的目标函数F(p)，
<mrow> <mtable> <mtr> <mtd> <mrow> <mi>r</mi> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>=</mo> <msup> <mrow> <mo>&amp;lsqb;</mo> <msub> <mi>r</mi> <mn>1</mn> </msub> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>,</mo> <msub> <mi>r</mi> <mn>2</mn> </msub> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>,</mo> <mo>...</mo> <mo>,</mo> <msub> <mi>r</mi> <mi>N</mi> </msub> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>&amp;rsqb;</mo> </mrow> <mi>T</mi> </msup> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mi>r</mi> <mi>j</mi> </msub> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>=</mo> <msub> <mi>Y</mi> <mi>j</mi> </msub> <mo>-</mo> <mi>y</mi> <mrow> <mo>(</mo> <mi>p</mi> <mo>,</mo> <msub> <mi>X</mi> <mi>j</mi> </msub> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </mrow>
<mrow> <mi>F</mi> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfrac> <mn>1</mn> <mn>2</mn> </mfrac> <mi>r</mi> <msup> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mi>r</mi> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfrac> <mn>1</mn> <mn>2</mn> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>j</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </munderover> <msup> <mrow> <mo>(</mo> <msub> <mi>Y</mi> <mi>j</mi> </msub> <mo>-</mo> <mi>y</mi> <mo>(</mo> <mrow> <mi>p</mi> <mo>,</mo> <msub> <mi>X</mi> <mi>j</mi> </msub> </mrow> <mo>)</mo> <mo>)</mo> </mrow> <mn>2</mn> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </mrow> ；
针对公式(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个元素，
<mrow> <mi>J</mi> <mrow> <mo>(</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mtable> <mtr> <mtd> <mrow> <mo>-</mo> <mfrac> <mn>1</mn> <msubsup> <mi>p</mi> <mn>1</mn> <mn>2</mn> </msubsup> </mfrac> <msub> <mi>p</mi> <mn>2</mn> </msub> <msub> <mi>X</mi> <mn>1</mn> </msub> <msup> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mn>1</mn> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <mrow> <msub> <mi>p</mi> <mn>2</mn> </msub> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mi>y</mi> <msup> <mrow> <mo>(</mo> <mi>p</mi> <mo>,</mo> <msub> <mi>X</mi> <mn>1</mn> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> <mtd> <mrow> <mi>ln</mi> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mn>1</mn> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <msup> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mn>1</mn> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <msub> <mi>p</mi> <mn>2</mn> </msub> </msup> <mi>y</mi> <msup> <mrow> <mo>(</mo> <mi>p</mi> <mo>,</mo> <msub> <mi>X</mi> <mn>1</mn> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> </mtr> </mtable> </mtd> </mtr> <mtr> <mtd> <mtable> <mtr> <mtd> <mrow> <mo>-</mo> <mfrac> <mn>1</mn> <msubsup> <mi>p</mi> <mn>1</mn> <mn>2</mn> </msubsup> </mfrac> <msub> <mi>p</mi> <mn>2</mn> </msub> <msub> <mi>X</mi> <mn>2</mn> </msub> <msup> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mn>2</mn> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <mrow> <msub> <mi>p</mi> <mn>2</mn> </msub> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mi>y</mi> <msup> <mrow> <mo>(</mo> <mi>p</mi> <mo>,</mo> <msub> <mi>X</mi> <mn>2</mn> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> <mtd> <mrow> <mi>ln</mi> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mn>2</mn> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <msup> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mn>2</mn> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <msub> <mi>p</mi> <mn>2</mn> </msub> </msup> <mi>y</mi> <msup> <mrow> <mo>(</mo> <mi>p</mi> <mo>,</mo> <msub> <mi>X</mi> <mn>2</mn> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> </mtr> </mtable> </mtd> </mtr> <mtr> <mtd> <mn>...</mn> </mtd> </mtr> <mtr> <mtd> <mtable> <mtr> <mtd> <mrow> <mo>-</mo> <mfrac> <mn>1</mn> <msubsup> <mi>p</mi> <mn>1</mn> <mn>2</mn> </msubsup> </mfrac> <msub> <mi>p</mi> <mn>2</mn> </msub> <msub> <mi>X</mi> <mi>N</mi> </msub> <msup> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mi>N</mi> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <mrow> <msub> <mi>p</mi> <mn>2</mn> </msub> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mi>y</mi> <msup> <mrow> <mo>(</mo> <mi>p</mi> <mo>,</mo> <msub> <mi>X</mi> <mi>N</mi> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> <mtd> <mrow> <mi>ln</mi> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mi>N</mi> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <msup> <mrow> <mo>(</mo> <mfrac> <msub> <mi>X</mi> <mi>N</mi> </msub> <msub> <mi>p</mi> <mn>1</mn> </msub> </mfrac> <mo>)</mo> </mrow> <msub> <mi>p</mi> <mn>2</mn> </msub> </msup> <mi>y</mi> <msup> <mrow> <mo>(</mo> <mi>p</mi> <mo>,</mo> <msub> <mi>X</mi> <mi>N</mi> </msub> <mo>)</mo> </mrow> <mn>2</mn> </msup> </mrow> </mtd> </mtr> </mtable> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>5</mn> <mo>)</mo> </mrow> <mo>;</mo> </mrow>
步骤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.根据权利要求1所述的一种雷电流幅值累积概率分布曲线拟合计算方法，其特征在于，在所述步骤2)中，在N个散点中，通过线性插值寻找中值电流拟合值a1*，基于中值电流拟合值a1*反算陡度参数b的拟合值b1*。
3.根据权利要求2所述的一种雷电流幅值累积概率分布曲线拟合计算方法，其特征在于，在所述步骤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)计算得出，
<mrow> <msub> <mi>a</mi> <mn>1</mn> </msub> <mo>*</mo> <mo>=</mo> <mfrac> <mrow> <mo>(</mo> <mn>0.5</mn> <mo>-</mo> <msub> <mi>Y</mi> <mi>m</mi> </msub> <mo>)</mo> <mo>(</mo> <msub> <mi>X</mi> <mrow> <mi>m</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>-</mo> <msub> <mi>X</mi> <mi>m</mi> </msub> <mo>)</mo> </mrow> <mrow> <msub> <mi>Y</mi> <mrow> <mi>m</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>-</mo> <msub> <mi>Y</mi> <mi>m</mi> </msub> </mrow> </mfrac> <mo>+</mo> <msub> <mi>X</mi> <mi>m</mi> </msub> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>6</mn> <mo>)</mo> </mrow> </mrow> ；
然后基于中值电流拟合值a1*反算陡度参数b的拟合值b1*，对公式(2)作变换，按公式(7)计算b的拟合值b1*，
<mrow> <msub> <mi>b</mi> <mn>1</mn> </msub> <mo>*</mo> <mo>=</mo> <mfrac> <mn>1</mn> <mi>N</mi> </mfrac> <munderover> <mi>&amp;Sigma;</mi> <mrow> <mi>j</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </munderover> <mfrac> <mrow> <mi>ln</mi> <mrow> <mo>(</mo> <msubsup> <mi>Y</mi> <mi>j</mi> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msubsup> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow> <mrow> <mi>ln</mi> <mrow> <mo>(</mo> <msub> <mi>X</mi> <mi>j</mi> </msub> <mo>/</mo> <msub> <mi>a</mi> <mn>1</mn> </msub> <mo>*</mo> <mo>)</mo> </mrow> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>7</mn> <mo>)</mo> </mrow> <mo>.</mo> </mrow>
CN201711328027.3A 2017-12-13 2017-12-13 一种雷电流幅值累积概率分布曲线拟合计算方法 Active CN108051676B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201711328027.3A CN108051676B (zh) 2017-12-13 2017-12-13 一种雷电流幅值累积概率分布曲线拟合计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201711328027.3A CN108051676B (zh) 2017-12-13 2017-12-13 一种雷电流幅值累积概率分布曲线拟合计算方法

Publications (2)

Publication Number Publication Date
CN108051676A true CN108051676A (zh) 2018-05-18
CN108051676B CN108051676B (zh) 2020-04-21

Family

ID=62132137

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201711328027.3A Active CN108051676B (zh) 2017-12-13 2017-12-13 一种雷电流幅值累积概率分布曲线拟合计算方法

Country Status (1)

Country Link
CN (1) CN108051676B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110031676A (zh) * 2019-03-21 2019-07-19 国网浙江省电力有限公司电力科学研究院 基于傅里叶变换的电晕电流脉冲波形拟合方法
CN110363345A (zh) * 2019-07-12 2019-10-22 国网浙江省电力有限公司电力科学研究院 雷电流幅值概率预测方法及***

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102854364A (zh) * 2012-08-06 2013-01-02 广东电网公司东莞供电局 基于spss分段拟合地区性雷电流幅值概率分布的方法
CN104239739A (zh) * 2014-09-26 2014-12-24 华南理工大学 一种雷电流幅值概率分布函数的分段拟合方法与***
CN105912509A (zh) * 2016-04-28 2016-08-31 国网电力科学研究院武汉南瑞有限责任公司 一种多区段雷电流幅值累积概率分布统计方法
DE102016000930A1 (de) * 2015-09-04 2017-03-09 DEHN + SÖHNE GmbH + Co. KG. Verfahren zur Erfassung von Blitzstromparametern an Anlagen mit einer oder mehreren Fangeinrichtungen und Blitzstromableitpfaden

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102854364A (zh) * 2012-08-06 2013-01-02 广东电网公司东莞供电局 基于spss分段拟合地区性雷电流幅值概率分布的方法
CN104239739A (zh) * 2014-09-26 2014-12-24 华南理工大学 一种雷电流幅值概率分布函数的分段拟合方法与***
DE102016000930A1 (de) * 2015-09-04 2017-03-09 DEHN + SÖHNE GmbH + Co. KG. Verfahren zur Erfassung von Blitzstromparametern an Anlagen mit einer oder mehreren Fangeinrichtungen und Blitzstromableitpfaden
CN105912509A (zh) * 2016-04-28 2016-08-31 国网电力科学研究院武汉南瑞有限责任公司 一种多区段雷电流幅值累积概率分布统计方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
KARL LUNDENGARD等: "Estimation of Pulse Function Parameters for Approximating Measured Lightning Currents Using the Marquardt Least-Squares Method", 《2014 INTERNATIONAL SYMPOSIUM ON ELECTROMAGNETIC COMPATIBILITY》 *
李瑞芳等: "雷电流幅值概率计算公式", 《电工技术学报》 *
赵淳等: "输电线路走廊雷电流幅值分布统计方法", 《高电压技术》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110031676A (zh) * 2019-03-21 2019-07-19 国网浙江省电力有限公司电力科学研究院 基于傅里叶变换的电晕电流脉冲波形拟合方法
CN110363345A (zh) * 2019-07-12 2019-10-22 国网浙江省电力有限公司电力科学研究院 雷电流幅值概率预测方法及***

Also Published As

Publication number Publication date
CN108051676B (zh) 2020-04-21

Similar Documents

Publication Publication Date Title
CN112819312B (zh) 气候变化情景下干旱社会经济暴露度评估方法和***
CN105320809B (zh) 一种针对风电场空间相关性的风速预测方法
CN106600037B (zh) 一种基于主成分分析的多参量辅助负荷预测方法
CN107609774B (zh) 一种基于思维进化算法优化小波神经网络的光伏功率预测方法
CN106908668B (zh) 一种实测地面合成电场数据的处理方法及***
CN107220907B (zh) 一种采用秩和比综合评价的谐波污染用户分级方法
CN108051676B (zh) 一种雷电流幅值累积概率分布曲线拟合计算方法
CN106202707B (zh) 一种基于灰色置信区间的结构应力-强度干涉模型集合可靠性分析方法
CN106940830B (zh) 未来气候变化对生物多样性影响与风险综合评估技术
CN111414703A (zh) 一种滚动轴承剩余寿命预测方法及装置
CN102682348A (zh) 复杂装备部件维修级别优化***及其建立方法
CN111931395A (zh) 一种降低应变场重构误差的传感器测点优化方法
Jiang et al. Adaptive gaussian process for short-term wind speed forecasting
CN108830405B (zh) 基于多指标动态匹配的实时电力负荷预测***及其方法
CN109827662B (zh) 基于逆高斯分布低值绝缘子红外检测温度阈值的判定方法
Zulfi et al. The development rainfall forecasting using Kalman filter
CN111597766A (zh) 基于粒子滤波采样的全过程动态仿真长期孤网稳定性预测方法及***
CN111104640A (zh) 一种基于层次分析法的降雨观测评价方法及***
CN109934394A (zh) 一种基于灰色和马尔科夫理论的需求侧响应预测方法
CN113610436B (zh) 一种承灾体动态脆弱性评估方法及***
CN107590346B (zh) 基于空间多重相关解集算法的降尺度校正模型
CN110598914A (zh) 一种多因素影响下矿井灾害气体浓度区间预测方法及***
CN114092545A (zh) 一种适用于球形靶标拟合的自适应栅格搜索方法
CN110472818B (zh) 一种快速评价扰动湿地恢复力的方法
CN110221147B (zh) 基于多复合优化算法的电能质量检测分析方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant