高精度计算（共8篇）
高精度计算 篇1
摘要：提出一种电力谐波参数高精度测量的新方法。依据傅立叶变换的时域收缩频域延伸的性质精确计算出电力信号各次谐波的频率, 而后通过窗函数插值方法精确计算出电力信号的各次谐波参数。讨论了时域收缩比对测量精度的影响, 还分析了截断窗函数对电力谐波估计精度的影响。用仿真试验和模拟测量证明了该算法的有效性。与其他加窗FFT插值校正分析方法相比, 所提方法计算耗时上有明显优势, 很适合用于DSP数字信号处理器, 是一种很实用的算法。
关键词：FFT,谐波分析,窗函数,频域内插
0 引言
随着电力电子技术和器件的发展, 非线性负荷在电力系统中的应用越来越广泛, 电力系统谐波污染日益严重, 谐波已成为影响电能质量的主要问题[1]。对谐波分量参数的高精度估计将有利于电能质量的评估和采取相应的必要治理措施。
快速傅立叶变换 (FFT) 是谐波分析最快捷的工具[2]。但非同步采样时, 快速傅立叶变换应用于谐波分析容易造成频谱泄漏和栅栏效应, 影响谐波相量计算的准确度。针对FFT算法的不足, 国内外学者提出的一系列加窗插值FFT算法[3,4,5]、多项式拟合的双谱线插值方法[6,7,8]和其他加窗FFT校正算法[9,10,11,12,13]等都有效地提高了谐波分析精度。但随着插值修正曲线拟合函数的阶次增高及谐波含有次数的增多, 谐波估计精度提高的同时计算量大量增加。本文将提出谐波高精度估计的另一条思路, 基于傅立叶变换时域收缩频域延伸的性质计算出各谐波频率, 进而精确计算谐波各参数的方法。
1 基于傅立叶变换时频收伸性的测频原理
设被分析电力信号为
式中:p∈[0, P], p为谐波次数, 是正整数;P为最高谐波次数;Ap为第p次谐波幅值;f1为电力信号基波频率;φp为第p次谐波初相角。
采样被分析电力信号得
式中:Ts为采样周期;wp=2πpf1 Ts是p次谐波信号的实际数字角频率。
对x (n) 加对称余弦窗 (Hanning窗、Nattall窗等) 函数w (n) , 得序列
xw (n) 的连续谱函数为
用FFT可求出xw (n) 的离散谱Xw (k) 实质上就是连续谱Xw (ejw) 在区间[-π, π]上以等间隔Δw=2π/N抽样的结果, N为分析数据截断长度。即
考虑到Xw (k) 的对称性, 仅用正频段来分析Xw (k) , 将上式简化为
若从截断数据N的起始点向后延伸用相同窗函数截断M长度数据得xw (m) , [m=0, 1, …, M], M
设第kp条和第lp条谱线分别为两次FFT变换p次谐波的最大谱线, 由于xw (m) 与xw (n) 有相同的相位, 则两次FFT变换的相位差为零, 即
式中:φpk和φpl分别为Xw (k) 和Xw (l) 在p次谐波处的峰值谱线kp和lp的相位;Δφpkl为Xw (k) 和Xw (l) 在p次谐波处的峰值谱线的相位差。上式简化中考虑了N□1和M□1, 即认为N-1=N, M-1=M。
因此得
2 时域收缩比M/N对测频误差的影响
式 (9) 可求得各次谐波的频率。由于各次谐波的频率是基波频率的整数倍, 而一般由于基波幅值大, 谐波相对于基波小, 谐波对基波的谱间干扰相对小, 因此, 当只需计算各次谐波的参数时, 可只计算出基波的频率, 各次谐波的频率乘以一个整数即可。
本文所提的谐波参数的测量是通过在频域内对窗函数插值获得的, 各谐波的幅值是先计算出它们的频率, 而后在频域内插值求得的, 因此, 各谐波频率的测量精度是谐波参数测量精度的关键。设f1*为电网的基波实际频率, f1为式 (9) 电网基波频率的测量值, 则频率测量的绝对误差为:f1-f1*。图1是采样频率Ts=5 120 Hz, 用Nattall窗函数截断数据长度N=1024时, 不同M/N值条件下计算基波频率误差曲线。观察频率误差曲线, 总体来说这种测频方法在10×M/N为整数时的测量精度都比较高, 频率测量的绝对误差小于0.001 (10-3) Hz, 并且, 曲线呈现两边高中间低的趋势, 因此可推测, 从图1也可看出取10M/N=5, 即M/N=0.5是恰当的。
3 窗函数插值计算电力谐波参数
设δωp为p次谐波处的峰值谱线的数字角频率与p次谐波信号的实际数字角频率的差值 (如图2) 。
再设βp=Ap/Xw (kpΔ-p) 为基波和各次谐波的校正系数。
由式 (6) 可得
由此可得P次谐波的幅值为
p次谐波的相位为
从式 (6) 和式 (7) 来看, 相位谱与截断窗函数无关, 但实际上, 由于φpk和φpl是由Xw (kp) 和Xw (lp) 的虚部和实部分别计算得到, 而Xw (kp) 和Xw (lp) 的虚部和实部的计算精度都要受“谱间干扰”的影响, 因此, 窗函数对φpk和φpl的计算精度有直接的影响。选择合适的窗函数截断分析数据对提高电力谐波参数的计算精度有不小的影响。加窗的目的是降低旁瓣“谱间干扰”, 而不能减小主瓣“谱间干扰”。人们希望窗函数频谱的旁瓣最大幅值愈小、旁瓣幅值衰减愈快愈好, 这意味着主瓣功率的集中度愈高愈好。主瓣功率的集中度与主瓣宽度相关联, 一般说来, 主瓣功率愈集中, 主瓣宽度也愈宽, 测量精度愈高。如汉宁窗的主瓣宽度为4个谱线分辨率, 即窗宽为4Δw, 而纳托尔窗的主瓣宽度为8个谱线分别率, 即8Δw, 纳托尔窗的主瓣功率集中度比汉宁窗高, 图3误差曲线也明显显示采用纳托尔窗截断分析数据具有更高的计算精度。图3是不同窗函数截断时的测量误差曲线, 它是由式 (1) 模型, 基波频率50.3 Hz, 以频率5120 Hz采样, 截断数据长度为N=1024, 取M=N/2进行仿真计算得到, 模型参数见表1。由图3曲线可见, 窗函数对频率的测量精度的影响还是比较大的, 选择窗函数旁瓣幅值小并且衰减快的窗函数 (如Nattall) 将有利于提高测量精度。
从式 (9) 、式 (11) 和式 (13) 来看, 本文方法与其他校正FFT类方法一样, 也可用于间谐波的分析。然而, 就校正FFT这类方法而言, 能否精确地分析被分析信号取决于被分析信号中信号频率间隔是否大于截断窗函数主瓣宽度, 因为加窗只能降低旁瓣“谱间干扰”, 而不能减小主瓣“谱间干扰”。按IEC标准[14], 谐波和间谐波的分析数据长度为10周波 (50 Hz) , 因此, 相应的分析频率分辨率为5 Hz。因此, 如果仅测量电力谐波, 所选窗函数的主瓣宽度只要小于等于10个谱线分别率Δw, 就不会发生“主瓣的谱间干涉”, 计算精度一般都可接受。当间谐波的频率与频率最近的谐波频率间隔比较大 (大于窗函数的主瓣宽度) 时, 也可同时测量间谐波, 但当2者的频率间隔小于窗函数主瓣宽度时, 将出现不能接受的测量误差。减小窗函数主瓣宽度有2条路径:一是选择主瓣宽度相对窄的窗函数, 如汉宁窗只有4Δw宽。当然这会降低测量精度, 但实际上, 一方面, “国标”对谐波测量精度的最高要求0.05%UN (Up<1%UN) , 另一方面, 由于信号传感器精度的限制, 追求过高的理论分析精度也没有太大意义, 一般实际测量绝对误差能达到10-2足已, 可见, 适当考虑选择主瓣宽度相对窄的窗函数是现实可行的。二是延长数据的截断长度, 频率分辨率与被分析数据截断长度成反比, 截断分析数据长度越长Δw越小。例如用汉宁窗截断长度为20周波 (200 ms) 的被分析数据, 则4Δw=10 Hz, 也即只要被分析信号中的任何谐波与间谐波的频率间隔都大于10 Hz, 分析误差是可接受的。倘若被分析电力信号中有频率间隔小于10 Hz的两个谐波和间谐波信号, 则它们在分析时将发生“主瓣谱间干涉”, 就可能带来难以接受的测量误差。因此, 一般说来, 校正FFT这类算法不太适用于复杂谐波 (既含有谐波, 又含有间谐波) 的分析。
对照文献[3-13]可见, 本文所提的电力谐波测量计算方法与各种FFT插值校正法或双谱线拟合方法相比有基本相同的测量计算精度, 但在计算耗时上则大大减小, 可应用于嵌入式系统。并且, 本文所提算法的主要计算量在2次FFT上, 而DSP信号处理器内带有硬件FFT, 不需要占用计算时间, 因此, 该算法更适合于应用在DSP数字信号处理器上, 用单片DSP数字信号处理器即可实现电力谐波的在线测量。
4 实验验证
本文采用如图4方法来模拟测量。首先由Matlab按式 (1) 模型和表1参数建模, 然后用d SPACE的DS1104板进行D/A转换, 将所建模型转换成模拟信号, 最后用“合众达公司”的DSP2812板, 通过16位A/D采样转换为数字信号, 并用本文算法计算, 完成模拟测量。A/D采样频率为5 120Hz、分析数据长度为1 024点 (10个周波) 。采用该模拟测量方法对验证算法的有效性有两个优势:一是它很好地模拟了实际测量过程, 分析数据中包含了信号调理误差、A/D采样误差和采样噪声;二是可精确知道被分析信号的参数。若采用实测数据分析, 由于并不知道被测实际信号的参数, 需要有一测量仪测出实际信号的参数。而任何测量仪是有误差的, 我们能知道的也只能是实际谐波信号有误差的“测量参数”, 显然, 用这一有误差的参数为“标准”来评价某算法的有效性有偏颇。
图5是采用本文所提方法, 用纳托尔窗截断1 024点数据, 由图4模拟测量得到的误差曲线。由误差曲线可见, 本文所提出的方法有足够用于实际测量的精度, 满足“国标”对谐波测量精度的最高要求0.05%UN (Up<1%UN) , 而计算量则大大降低。
5 结论
本文提出了一种电力谐波分析的新方法, 即先高精度地估计出含有谐波的电力信号中的基波参数 (幅值、频率和相位) , 进而, 通过窗函数频域内插值的方法精确计算出各次谐波的参数。通过同一电力信号的仿真试验证明, 本文提出的电力谐波高精度估计方法与各种加窗插值估计电力谐波的分析方法有基本等同的估计精度, 而在计算速度上有明显的优势, 并且, 本文所提算法的主要计算量在2次FFT上, 而DSP信号处理器内带有硬件FFT, 不需要占用计算时间, 因此, 该算法更适合于应用在DSP数字信号处理器上, 用单片DSP数字信号处理器即可实现电力谐波的在线测量, 因此, 是一种很有实用价值的对平稳电力谐波进行高精度估计的方法。
高精度计算 篇2
通过分析微小型组合导航系统中导航计算机的.功能要求,确定了导航计算机的硬件设计指标.提出以高性能浮点型TMS320C6713数字信号处理器(DSP)为核心处理器芯片,设计了高速18位差分A/D转换电路、16位D/A转换电路等多种接口电路,运用现场可编程逻辑技术有效地降低了系统设计复杂度.设计并实现了可靠的程序实现流程,使捷联惯性航姿系统脱离通用计算机平台.实验证明:基于DSP的组合导航系统的精度和实时性能够达到设计要求,并且,便携、价廉,对导航系统在微小型领域的广泛使用具有实际意义.
作 者：孙永荣 刘建业 刘瑞华 杜亚玲 SUN Yong-rong LIU Jian-ye LIU Rui-hua DU Ya-ling 作者单位：孙永荣,刘建业,杜亚玲,SUN Yong-rong,LIU Jian-ye,DU Ya-ling(南京航空航天大学,自动化学院,江苏,南京,210016)
刘瑞华,LIU Rui-hua(中国民航学院,空中交通管理学院,天津,300300)
高精度计算 篇3
随着计算机技术和应用的迅速发展,计算机网络已经深入到各个领域并发挥着及其重要作用,而网络中计算机的时间同步精度直接影响着各个网络系统的安全工作,特别是对于金融、通信、交通和军事等重要领域,计算机的时间同步问题就更显其重要性,对计算机时间同步精度测试是计算机时间同步的重要保障。
然而,到目前为止,还没有一个很好的解决方案,来实现对计算机时间同步误差的测试,通常测试方法是利用软件来获取2台计算机之间的时间差来完成,这种方法不仅精度低、不直观,而且无法对通过物理信号(如B码)实现计算机同步的情况来测试。通过研究和大量实验,探索出一种使用物理方法来进行计算机时间同步精度的测试方法,方法简单,测试数据可在时间间隔计数器上直接读取,测试结果一目了然,测试精度高,稳定可靠,适宜广泛推广使用。
1 解决方案
1.1 测试方法
在典型时间频率系统中,经常利用外部高精度时间频率信号来同步系统内部的计算机网络,如图1所示,5 MHz频率信号经过分频钟产生B码及1pps信号,时间服务器通过B码终端来校准本机时间,并通过时间同步协议(如最常用的NTP协议)使网内的其他计算机与自身同步,从而实现整个网络内计算机的时间同步。
如果计算机能够直接产生与自身时间同步的秒脉冲信号,则可以如图1所示,直接通过时间间隔计数器,测量时间源、时间服务器及各计算机之间的时间偏差,从而完成时间同步误差的测量。因此,该测试方法的关键技术即是控制计算机产生与自身时间同步的秒脉冲信号。
1.2 测试原理
利用计算机产生秒脉冲信号,可以尝试利用硬件电路,通过对计算机时间中断信号的访问,实现秒脉冲信号的输出,但事实上这种方法实现难度很大,制作成本也较高,不适宜推广使用。因此,通过大量的试验与研究,探索出通过软件来实现对计算机时钟频率的精确读取,并利用RS-232口的管脚来输出秒脉冲信号。
1.2.1 软件实现
每台计算机内部均有高精度定时器,就跟硬件时钟的晶振一样精确,在Windows平台环境下,系统提供了2个API函数来实现对该定时器的访问,这2个函数的声明如下:
BOOL QueryPerformanceCounter(LARGEINTEGER *lpPerformanceCount );
BOOL QueryPerformanceFrequency( LARGEINTEGER *lpFrequency );
其中,QueryPerformanceFrequency函数用来获取定时器的时钟频率,QueryPerformanceCounter函数用来获取定时器当前的计数值,因此,可以利用先后2次获得的计数值之差及时钟频率,计算出精确时延,或者实现精确定时,精度可达到μs级,其示例代码如下:
QueryPerformanceFrequency(&freq); // 得到定时器时钟频率
QueryPerformanceCounter(&li1); // 得到定时器当前计数器值
do
{
QueryPerformanceCounter(&li2); // 得到定时器当前计数器值
delay=(double)((li2.QuadPart li1.QuadPart) /freq.QuadPart);
} while(delay < uSec/1000000); // uSec为延时量,单位μs
由于Windows是一个分时系统,不是实时系统,通过系统提供的时间获取函数(如GetSystemTime等)得到的时间很不可靠,分辨率在10～15 ms左右,因此,要获得计算机准确的秒时刻,通过时间函数无法完成,但可以通过时间函数来获取尽量靠近整秒的时刻,剩余的时间利用前面所述的高精度定时器来精确延时,从而到达尽可能准确地整秒时刻,其实现流程如图2所示。
1.2.2 物理输出
计算机秒脉冲物理信号的输出可通过计算机RS-232口来实现,根据RS-232串口通信协议,可知DTR(数据准备)与RTS(请求数据)是利用高低电平来进行握手的,因此可根据针脚上的电平控制来输出秒脉冲信号,其示例代码如下:
EscapeCommFunction(hCom,SETRTS); // 设置RTS信号,置高电平
DelayUsec(1000); // 设置1 ms的时延(脉宽)
EscapeCommFunction(hCom,CLRRTS); // 清除RTS信号,置低电平
因此,只需要将普通串口线的另一端改装成BNC头,然后在测试计算机上安装秒信号产生软件,即可实现计算机秒脉冲的输出,从而方便地实现了图1所示的计算机时间同步精度测量,由于测试程序相同,测试线路可控,可以认为此结果即为服务器与终端的秒信号时间同步误差。
2 测试效果验证
2.1 测试环境
测试时,原子钟5 MHz信号经分频钟产生B码信号和秒信号,时间同步软件实时接收B码信号对计算机时间进行校准,秒脉冲产生软件运行在同一台计算机上,实时输出同步秒脉冲,2路秒信号送时间间隔计数器进行时差测量。测试设备如表1所示。
2.2 测试结果
图3为时间同步服务器输出秒与分频钟输出秒经过电缆时延修正后的钟差比对图,图中的部分数据抖动是因为时间同步软件的同步精度影响,从长期趋势来看,测试结果稳定,表明计算机秒脉冲输出稳定,因此可以得出,该项计算机时间同步误差测试方法是有效、可靠的。
3 结束语
通过本文描述的方法实现了利用物理手段对计算机同步精度的测试,测试精度高,方法先进,设备量少,操作简单,是目前计算机时间同步精度测试的先进方法,但由于此方法是物理信号直接测量,因此远距离的计算机无法通过此方法进行时间同步测量。此方法已经在windows XP、windows 2000 操作系统平台,并具备RS-232接口的计算机上进行了多次验证,并取得很好的效果。
摘要：针对目前计算机时间同步误差测试方法精度低、不直观等问题,研究并完成了一种利用脉冲信号进行时间间隔测量实现计算机时间同步误差的测试方法。该方法可以通过程序控制计算机输出物理脉冲信号进行时间间隔测量。介绍了该测试方法的原理和实现方案,并通过实验数据验证了该测试方法的先进性。该方法已经在Windows2000、WindowsXP、Windows2003等操作系统中成功验证,测试结果直观、精度高,可以广泛推广使用。
关键词：计算机,时间同步,秒脉冲,高精度定时器
参考文献
[1]龚建伟.VC++/Turbo C串口通信编程实践[M].北京:电子工业出版社,2004.
[2]邱仲潘.VC++6从入门到精通[M].北京:电子工业出版社,2005.
高精度计算 篇4
关键词：叠前反演,Zoeppritz方程,偏导计算,反演精度,方法验证
叠前地震反演技术能够利用不同偏移距地震变化信息提取地层弹性参数, 对储层岩性、物性以及含流体性研究具有十分重要意义。不断发展的石油勘探对叠前反演储层预测精度要求越来越高, 尽管精度影响因素众多, 但叠前反演方法是影响叠前反演技术精度的重要原因之一。
叠前反演技术的理论基础平面波非垂直入射理论, 在同一界面中地震波从不同角度入射会产生地震波的反射和折射。并且反射系数和折射系数随着入射角度变化而变化。Zoeppritz方程给出了其解析数学表达式, 精确描述了反射系数与入射角度、速度以及密度间的关系。叠前反演的核心就是计算反射系数与弹性参数的梯度变化关系;但由于Zoeppritz方程解析表达式非常复杂, 在叠前反演技术中直接构建介质弹性参数与反射系数的梯度变化关系难度大, 地球物理学者基于弱反射界面以及小角度入射的近似假设推导出了该方程的简化表达式[1—3], 目前大多数叠前反演方法通常采用简化表达式推导构建反射系数梯度方程;但由于简化式的近似假设会存产生两方面的误差:一方面上下地层弹性参数变化大的情况下, 简化式反演会产生一定的计算误差;另一方面在大入射角度地震反演计算时, 误差也会比较大。随着油气勘探深入, 地质目标越来越复杂, 对勘探精细程度要求越来越高, 需要提高叠前地层弹性参数反演的精度;同时由于地震采集技术的发展, 大偏移距宽入射角度的地震资料也越来越多, 基于简化式的叠前反演方法在这些地质目标和地震资料应用中会产生较大的误差, 反演准确性差和精度低。
针对目前叠前反演方法存在的问题, 直接由Zoeppritz方程计算反射系数与弹性参数的关系, 推导构建了基于Zoeppritz方程的高法精度和高效率的叠前高精度反演新方法。该方法没有近似条件限制, 对复杂地质条件和大入射角度地震资料都有较好适应性, 可为储层预测研究提供更可靠的反演成果。利用二维模地震型数据和实际地震资料对研究的高精度叠前反演方法进行了验证, 并与目前常用的基于简化式叠前反演方法进行了对比, 分析认为新方法对地震资料适应性更强, 储层表征和分辨能力更好。
1 叠前反演目标函数构建
根据广义线性反演的思路需要建立以下叠前反演目标函数
式中Vp表示纵波速度, Vs表示横波速度, ρ为密度, D为实际角道集记录, S为期望地震模型响应。R (Vp, Vs, ρ) 为纵波反射系数, W为地震子波。
将R (Vp, Vs, ρ) iΔ在初始模型响应处R (Vp, Vs, ρ) i泰勒展开, 为使运算简单通常将泰勒展开式中的二阶和二阶以上的省略, 得到
对两端求取ΔVp、ΔVs和Δρ的一阶偏导并令其为零从而达到目标函数最小, 令Δdi=D (Vp, Vs, ρ) iS (Vp, Vs, ρ) i, 则推出
通过这样一个方程来求解三个弹性参数的摄动量ΔVp、ΔVs和Δρ, 该方程涉及三个未知量, 为欠定方程求解时必然有很强的多解性。为减少其反演产生的多解性, 引入多个角度的地震资料来联合反演这三个参数。
若有三个角度的角道集数据, 分别为小角度D1、中角度D2和大角度D3, 以G (vp) 代替, 以G (vp) 代替, 以G (ρ) 代替, 则式 (3) 可化为:
将以上三式写为矩阵形式有
根据式 (7) 采用最小二乘法可得摄动量的表达式为
由于矩阵中元素G (vp) 为为, 而, 所以通过求解便可求的G (vp) , 即可以得到矩阵元素。求解方程 (8) 便可得到方程中的参数即反演结果, 而求解方程的核心是地震反射系数对介质弹性参数偏导关系的构建。
2 基于Zoeppritz方程的精确偏导计算
Zoeppritz方程给出了反射系数和入射角度、纵横波速度及密度等的复杂关系, 是叠前反演理论基础, 如式 (9) 。
式 (9) 中:Rp、Rs、Tp和Ts分别是纵波反射系数、转换横波反射系数、透射纵波、透射横波反射系数, α为纵波反射角, β转换横波反射角, α'透射纵波折射角, β'透射横波折射角, vp1、vs1、ρ1、vp2、vs2、ρ2分别为波在上下介质纵横波速度以及密度。
可令
则Zoeppritz方程可简写化为
式中:
对Zoeppritz方程进行简化得出的纵波反射系数与纵横波速度及密度的近似关系式, 相比复杂解析式描述准确性差, 基于简化式偏导方程的叠前反演运算精度低, 因此开展了基于Zoeppritz方程直接进行偏导方程计算构建研究。
利用利用斯奈尔定律和三角关系式用入射角度替代折射角和反射角, 减少Zoeppritz方程中的未知量个数。再在Zoeppritz方程中分别对上下地层纵横波速度以及密度6个变量求偏导得到6个偏导数方程组,
对所有6个变量可以统一用下式来表示:
式 (10) 中:R为反射系数, m=[vp1, vp2, vs1, vs2, ρ1, ρ2]。
通过计算A的逆以及A和B对弹性参数的偏导数就可以得到纵横波反射系数以及纵横波透射系数对上下地层纵横波速度及密度的偏导数, 式 (11) 为A对上介质纵波速度偏导。
对偏导为零的值进行了单独的分离, 提高计算精度的同时也提高了计算效率。这样就导出Jacobi偏导矩阵方程, 建立了各种反射系数对纵、横波速度、介质密度关系, 其中纵波反射系数对各个弹性参数的偏导是叠前高精度反演主要应用方程, 实现了基于Zoeppritz方程反射系数梯度矩阵的精确计算[4—6]。该偏导矩阵直接从Zoeppritz方程精确计算所得, 相比用其近似式构建的偏导方程更准确, 能够更高精度实现地层弹性参数的反演迭代运算。
3 二维模型验证
模型正演具有弹性参数明确、正演地震噪音少、反演子波确定等优势, 利用模型验证能够排除以上非方法因素对反演结果的影响, 因此首先基于二维正演数据开展高精度反演方法验证。设计二维理论模型弹性参数如图1所示, 弹性参数包含纵波速度、横波速度, 考虑到实际的地质弹性参数变化情况, 地层速度向深部逐渐增大。二维理论模型包含断层、楔状体、透镜体、薄互层等常见地质构造, 这些特殊构造的大小规模不等, 厚度不一, 模拟了地下复杂的地下情况。
考虑实际地震采集情况, 二维模型设计最大偏移距6 000 m, 最大深度2 300 m左右, 炮间距10 m, 估算最大入射角度达50°左右, 正演地震大偏移距信息多, 入射角度宽。采用声学波动方程正演得到了三个角度范围 (0°~15°, 16°~25°, 26°~30°) 地震剖面, 图2, 地震正演子波主频40 Hz左右, 频宽5~90 Hz, 与实际地震资料频谱相当。
基于正演分角度叠加地震分别开展了近似式和高精度反演处理, 并对两种反演结果开展了对比。图3为两种反演方法纵波速度反演结果对比剖面, 高精度反演对于层间差异表现更好, 无论是对薄层还是厚层, 高精度反演结果精度更高, 地层参数差异更加明显, 与模型吻合程度更好, 同时高精度反演结果对于断层面具有很好的表现。
4 实际资料验证
在模型验证基础上又开展了实际资料的验证, 选择目标为胜利油田KD地区浅层河流相储层, 目的层埋藏较浅 (1 100~1 500) m, 地震资料品质相对较好, 频带宽, 分辨率能力高, 能够避免一些资料因素带来的影响。地质目标单砂体薄 (2~12) m, 横向沉积变化快, 砂体纵向上相互叠置, 砂体关系复杂, 同时由于埋藏浅大偏移距信息相对丰富, 因此储层预测对叠前反演方法提出了更高的要求, 以往开展了很多叠前反演方法的应用未能满足储层预测高精度要求。
根据叠前反演对地震资料需要, 首先开展叠前精细速度分析和叠前偏移处理, 同时开展了能量均衡以及去噪处理, 保证了地震资料具有高信躁比和高保真性, 最后对分角度叠加地震资料分别开展了叠前高精度方法以及简化式方法反演处理, 可以得到纵横波速度、密度, 进一步可以得到提取储层敏感参数纵横波速度比、泊松比、拉梅系数乘密度, 利用地层敏感弹性参数对储层岩性流体识别精度更高。
图5为过KD地区KD104井的纵横波速度比反演剖面, 高精度反演结果可以清楚的表现三个河道砂体的展布, 砂体的形态和叠置关系清晰, KD104井钻遇了左边和中间的两个砂体, 反演结果与测井吻合较好。而简化式反演结果对砂体形态和叠置关系反映不清, 不能很好表现地震中储层的变化特征。实际资料表明高精度反演比简化式反演结果对储层的纵横向分辨能力更高, 特别是对于薄储层砂体表征能力更强。
5 结论
反射系数偏导方程是进行反演计算的核心, 目前叠前方法中大都应用Zoeppritz方程简化式计算的偏导方程, 由于简化式基于近似假设不能准确代表Zoeppritz方程所描述的反射系数与弹性参数间的关系, 因此造成反演精度低, 适应性差。直接利用Zoeppritz方程精确计算了反射系数偏导方程, 并对其开展了优化, 形成了基于精确反射系数偏导方程的叠前反演方法, 利用模型和实际资料对方法进行了测试, 相比基于简化式的叠前反演方法, 基于精确偏导方程的高精度叠前反演成果精度更高, 对大角度资料以及强变化界面适应性更强。
参考文献
[1] Bortfeld R.Approximations to the reflection andtransmission coefficients of plane longitudinal and transverse wave.Geophys Prosp, 1961; (4) :485—502
[2] Aki K I, Richards P G.Quantitative seismology.W.H.Freeman and Co, 1980:153—154
[3] Shuey, R T.A simplification of the Zoeppritz equations.Geophysics, 1985;50 (5) :609—614
[4] 孟宪军, 刘福平.地震反射系数的相角变化.石油地球物理勘探, 2010;45 (Z) :121—124Men X J, Liu F P.The phase angle change for seismic reflection coefficient.Oil Geophysical Prospecting, 2010;45 (Z) :121—124
[5] 刘福平, 孟宪军.反演纵横波速度的Jacobian矩阵及精确计算方法.中国科学:地球科学, 2010;40 (11) :1608—1616Liu F P, Men X J.Jacobian matrix for the inversion of P-and S-wave velocities and its accurate computation method.Sci China Earhth Sci, 2010;40 (11) :1608—1616
高精度计算 篇5
在各种测量及控制中, 通常都要求仪器体积小、便携, 所以单片机在测量仪器中的应用日益广泛。现代工业仪器智能化要求越来越高, 对复杂三角函数的计算量越来越大, 以目前单片机所具有的资源和技术能力, 要想达到上述要求, 对技术人员的编程技巧要求是非常高的, 工作量也非常庞大。基于以上原因, 目前特别需要研究出一种对单片机资源要求低, 易编程, 又能保证很高计算精度要求的算法。
1 单片机处理数据的常见方法
单片机处理数据的方法有两种:浮点数和定点数。浮点数数据运算在汇编环境下编程复杂、运算量大, 对内存容量要求高, 且对单片机的运算速度有很高要求, 这样普通51系列单片机无法满足复杂函数的计算。定点数在数据运算中运算量较小, 编程简单, 速度快, 但是精度低, 如何做到既不降低计算精度又能保证数据处理速度, 成为复杂函数在51单片机使用技术中的一个难点。
2 基本原理
求解三角函数采用的是综合运用数学中的插值法和查表法。
单纯使用查表法求解三角函数值时, 如果把0~90°每隔0.01°的函数值都存储在列表中, 需要存9 000个函数值, 一个函数值占两个字节, 9 000个函数值就是18 000 B, 如果精度要求更高所需的字节更大, 程序非常庞大, 还要占用很大的内存空间, 一般51系列单片机的内存已不能满足要求, 必须进行存储器扩展, 这样势必加大成本, 增大体积和功耗。所以单纯查表法求解三角函数计算精度就不能够要求太高, 否则占用存储空间太大。
利用插值法加查表法则既扩大了查表范围, 又不占用很大内存, 在表格中只存放小部分函数值, 其他函数值利用其相邻两个函数值, 通过插值计算来得到, 插值计算采用定点数运算, 计算速度比较快, 编程简单。
2.1 插值法原理
插值法原理如图1所示。
最简单的插值法为线性插值, 它通过相邻两个函数值来计算得出中间任意一点的函数值。设x为任意值, 它所对应的函数值为y, x1、x2为离x最近的两点, 它们相对应的函数值分别为y1、y2, 有x1
只要表格间距充分小, y的精度可以相当高。
2.2 插值法运用到三函数运算
根据上面公式将x换θ, 那么式 (1) 变为:
式中:θ为任意角度;θ2、θ1为θ2的相邻两点;y是θ函数值;y1、y2分别是θ2、θ1对应的函数值。
根据计算精度, 为使运算简便, 得出Δθ=0.32比较合适, 那么0.32×256=81.92°, 因此表格只要列出0~81.92°之间以0.32等间隔的角度函数值即可。
在单片机汇编环境下, 相邻两点的确定, 要遵循两点角度变换后的差是2的幂次方, 51单片机是8位机, 所以两点角度差应该是28, 即0X0100。
将式 (2) 作如下变换:
为了适合单片机编程, 设z、z1、z2、Δz分别为经过角度转换后的值, z=kθ, z1=kθ1, z2=kθ2, Δz=kθ2-kθ1, 因要使Δz=256, 那么256=k×0.32, 所以k=800, 即Δz=800×Δθ, 将角度扩大了800倍, 这样可以避免小数出现, 也有利于等间隔表的建立。z1、z2即为转换后的相邻两点。
在单片机中转换后的z应为十六进制, 即Δz=0X0100, 对应的十进制为256, 2Δz=0X0200, 对应的十进制为512, 依次类推。
经过以上角度变换后先将等间隔角度的函数值做成表, 放入存储中, 计算任意角度θ (0~81.92°) 的函数值时, 先进行角度转换:z=800×θ, 再换算成十六进制, 即z=0XH1H2H3H4, 取z的高字节为:
z1、z2即为与z相邻的两个值, 在表中查z1、z2相对应的函数值y1、y2, 则θ对应的函数值为:
y即为θ的函数值。
3 结语
该算法大大节省了单片机的存储器空间, 提高了运算速度, 且精度相当高。这种算法也可用于余弦、正切、余切或其他函数计算, 只需将公式中的参数做相应变动即可。特别适用于在低成本单片机上使用, 该算法已应用到89C51单片机上, 该单片机只有2 KB的存储容量, 计算精度达到了1/10 000。使用这种单片机可以降低设备的成本, 提高产品竞争力。所以这种算法, 具有很好的推广前景。
参考文献
[1]陈章龙.使用单片机大全[M].哈尔滨:黑龙江科技出版社, 1988.
[2]张友德.单片机微型机原理应用与实验[M].上海:复旦大学出版社, 1992.
[3]张毅刚.MCS-51单片机应用设计[M].2版.哈尔滨:哈尔滨工业大学出版社, 2004.
[4]徐瑞华.单片机原理与接口技术[M].北京:人民邮电出版社, 2008.
[5]白其峥.数学建模案例分析[M].北京:海洋出版社, 2000.
[6]李火林.数学模型及方法[M].南昌:江西高校出版社, 1997.
[7]陈理荣.数学建模导论[M].北京:北京邮电大学出版社, 1999.
[8]丁丽娟.数值计算方法[M].北京:北京理工大学出版社, 1997.
[9]刘满仓, 雷卫宁, 王春成.基于ARM的高精度数据采集系统设计[J].现代电子技术, 2012, 35 (8) :12-14.
二维数控精密转台精度计算与分析 篇6
关键词：二维数控,转台精度,指向误差
1 概述
针对非规则 (二维) 光学转鼓在扫描探测器中的广泛应用和产量需求, 自主研发能够和Ultraform250超精密切削机床配套使用的二维数控精密转台夹具。
由于二维数控精密转台是针对非规则光学转鼓的高效、高质量飞切加工的, 在加工过程中, 产生加工误差的因素很多。如果对转鼓规定一系列技术指标, 那么转台的有关误差对转鼓的技术指标的影响是很大的。因此, 要尽可能设法减少这些误差, 所占的比重越大, 留给补偿其他各类误差的空间就越小。其结果不是降低零件的加工精度, 就是增加加工难度。根据所生产转鼓的技术指标:平行度误差15角秒, 倾角误差20角秒, 转角误差20角秒, 对转台结构的精度进行逆推, 经过计算得出转台各主要误差指标。
2 转台的系统误差
由于转台的结构是二维数控工艺装置夹具结构, 对静态条件下的转台性能指标要求比较高, 下面对转台精度起主要影响作用的各项误差给出定义:
回转误差是指各轴在进行转动时, 轴端的运动轨迹所形成的包络线摆动范围的最小圆锥角度。它能够直接影响轴的空间指向。
垂直度误差是指转台空间直角坐标系中两轴轴线之间的实际空间夹角与90°角的差值。
位置误差是指在静止状态下, 转台的实际位置和理想位置的差值。
3 二维数控精密转台的指向误差
3.1 指向误差的定义 (见图1)
转台的指向误差指的是在横滚轴上确定一个单位向量λ0, 当转台转动一定的角度后, 转台的理想指向λ1与实际指向λ2之间的角度偏差。指向误差实际上是一种空间角度误差, 能够直接反应出转台的定位精度。在三维空间坐标中的单位向量的指向误差如图1所示。
指向误差的空间几何意义可以这样描述:设与横滚轴轴线重合的一个单位向量λ0, 按照欧拉变换的顺序, 将横滚轴和方位轴依次旋转一定角度, 则单位向量λ0回转到一个新的位置, 在新的方向上得到新的单位向量λ1, 则有:
式中, R0 (Ω) 为欧拉变换矩阵。
可是, 二维数控精密转台在做回转运动时, 会存在一系列的误差E, 这会导致单位向量λ0经过两轴回转后, 达不到预期的位置, 也就得不到新的单位向量λ1, 而是得到了单位向量λ2, 正因为这些误差的存在, 可推出公式:
这里误差E代表各个误差, 它主要由三类误差组成:轴之间的垂直度误差
α, 轴的回转误差β, 以及轴的位置误差γ, 则有:
那么指向误差则可表示为由单位向量λ1和单位向量λ2两个向量差的模长△Φ:
其中,
3.2 指向误差的算法
下面通过矩阵形式, 把指向误差的公式求解出来。
转台两轴依次绕X, Z轴转动X1, X2角度, 在没有误差的情况下, 单位向量λ0将变成λ1:λ1=R0 (Ω) λ0, 其中
由于存在误差因素E, 单位向量λ0会变成λ2。设垂直度误差α1;回转误差β1, β2;位置误差γ1, γ2。在转动两轴时, 根据轴的三类误差, 得到指向转换所产生的中间变量的欧拉变换矩阵:
X轴转动X1角度时, 产生位置误差γ1和回转误差β1:
Z轴转动x2角度时, 产生位置误差γ2和回转误差β2:
两轴之间垂直度误差α1:
于是, 得到含有三类误差的变换矩阵为
进而可以用公式 (4) , 求出指向误差△Φ。
4 转台的指向误差分配
由于指向误差受正交性误差α, 回转误差β和定位误差γ三类误差的影响, 但三类误差对指向误差的影响程度不同。所以要对指向误差进行分配。
根据控制变量法, 可计算出三类误差对指向误差的影响程度。当垂直度误差变化时, 分别取值0″, 1″, 2″, 3″, 4″, 5″得到各自的实际指向, 再计算相邻两个实际指向的余弦值, 从而得到相邻两个实际指向的夹角。
根据公式:
得到相邻两个实际指向的夹角:
同理, 当回转误差和定位误差变化时, 利用同样的方法, 观察实际指向的变化情况。则计算得出以下结论:
回转误差每增加1角秒, 指向误差大约增加1.4角秒;而垂直度误差和定位误差每增加1角秒, 指向误差均大约增加1角秒左右。
根据三类误差对指向误差的影响, 采用加权方法对三类误差进行误差分配, 设
其中, M、N、L为加权系数。
为了让产生指向误差的三类误差对指向误差的影响具有相同的程度, 则可确定加权系数为M=1, N=0.71, L=1。
由于三类误差均对转角误差、倾角误差有影响, 那么就可以把加权系数代入, 从而计算出三类误差。因为要保证转鼓相邻的两个加工面的转角误差不大于20角秒, 需要把各个面的加工误差控制在10角秒以内。计算得到横滚轴和方位轴的三类误差指标分别为:回转误差2.6角秒, 定位误差3.7角秒, 垂直度误差3.7角秒。同样, 方位轴在转动过程中, 也要保证每次回转时产生的倾角误差在10角秒以内, 则前面计算得到的三类误差满足要求。根据加工转鼓的加工工序可知道, 加工完一个面后, 需要转动横滚轴180°, 所以影响转鼓平行度误差的是横滚轴的回转误差和定位误差。由于横滚轴每次转动时产生的误差要控制在7.5角秒以内, 而前面计算出来的回转误差和定位误差之和小于这一数值, 则满足转鼓的加工要求。
结语
本文通过对二维数控精密转台的精度计算和分析, 确定了对转台精度起主要影响作用的三类误差指标。计算结果表明:三类误差对指向误差的影响程度不同。根据这一结论, 对误差进行分配。本文的结论可作为转台精度的技术要求, 为转台的设计提供依据, 因此对转台的制造具有重要的指导意义。
参考文献
[1]张文涛.金刚石飞切二维转鼓加工精密转台定位精度分析.新技术新工艺, 2011, 11:40-42
[2]白雪峰.单轴速率三轴位置惯性测试转台误差及传递分析.航天控制, 2006, 24 (02) :26-29
[3]Louis A DeMore.Design Study for a High-Accuracy Three-Axis Test Table.AIAA Guidance and Control Conference, 1985, 10 (01) :104-114
[4]周锐等.飞行仿真转台指向误差分析及误差分配.航天控制.1996, 14 (03) :66-72
提高工程量计算速度和精度的措施 篇7
在市场经济条件下, 我国建筑工程项目已实行公开招投标制度。在工程量清单计价模式下, 正确计算工程量是编制标书和控制工程造价的重要依据, 也是标书编制过程中最为重要和最为繁重的一项工作。只有工程量计算准确无误, 才能确保预算编制的质量。虽然工程计价软件解决了清单定额套用、工料分析、材料汇总、材料价差、费率等方面的计算问题, 但由于建筑工程图纸设计的类型多、变化大, 个别项目仍需要按图纸要求逐项计算工程量, 因此, 工程量计算仍然是概预算编制过程中一项非常繁重的工作任务。为了更好地开展工程量计算工作, 如何选择最佳项目计算线路组合, 如何优化计算程序、提高工作效率、减少重复劳动, 是每一位预算编审人员十分关注的问题。
1 工程量计算的方法
1.1 手工算量
自我国实行工程量计算方法以来, 随之就出现了手工算量, 目前这种算量方法仍然是我国工程量的算量主体, 对我国的基建工作做出了极其重要的贡献。在手工算量的长期应用和发展过程中, 许多熟练用户积累了丰富的工程量计算经验, 并总结形成了许多速算方法和速算表格, 在很大程度上提高了算量速度。
1.2 软件表格法算量
这种方法一般需要用户在软件中输入算量表达式, 程序自动进行汇总计算, 形成报表并打印。软件表格法算量实际上是用户手工算量方法的一种改进和延伸, 显著提高了用户的算量效率, 但并没有从根本上解脱算量人员的繁琐劳动。
1.3 软件自动算量
软件自动算量是目前算量方法中最具发展潜力的方法。该方法以计算规则为依据, 预算人员通过画图确定构件实体的位置, 并输入与算量有关的构件属性, 软件通过默认的计算规则, 自动计算得到构件实体的工程量, 自动进行汇总统计, 得到工程量清单。该算量方法简化了算量输入, 可以大幅度提高算量效率, 正在被越来越多的预算人员采用。目前, 工程量计算软件按照支持的图形维数不同, 又分为二维软件算量和三维软件算量两类。
2 提高工程量计算速度和精度的方法
2.1 确立基数主线
工程项目设计图纸是计算工程量的基础和依据, 预算编制人员首先必须熟悉图纸, 了解工程项目总体的平面布局、结构形式、设计标准和各分部的技术要求, 并确定计算过程中需要补充的相关资料, 估算出计算所需的工作日。在工程量计算过程中, 要统筹安排好计算程序, 尽量使前项计算结果成为以后各项目计算时共用的基数, 以减少重复劳动。根据共用基数这根主线编制的程序越科学, 计算就越简便、越快捷、越准确。
在计算前后有关联的各项工程量时, 应合理安排项目的计算顺序, 依照“由浅入深、先易后难、循序渐进”的原则进行。例如建筑面积是一个极其重要的多用途数据, 各项技术指标均要用此数据进行预算, 因而可放在较前的位置进行。在开始阶段, 编制人员对图纸尚不熟悉, 计算时可能有一定的难度, 但随着建筑面积计算工作的结束, 编制人员对施工图纸已比较熟悉, 对整个工程的建筑规模、设计标准、总体概况及各项技术要求已经明确, 计算基数主线也已确定, 此时再根据图纸计算整个工程项目的工程量就比较容易。
在工程量计算过程中, 确定基数主线是整个工程量计算的基础。例如在计算出底层建筑面积或楼 (地) 面面积后, 利用这一基数主线即可很快计算出楼面找平层、天棚抹灰、粉刷层的面积;用底层面积与垫层厚度相乘, 即可立即得到垫层的工程量、楼板的工程量等。应用基数主线的方法统筹进行计算, 不仅可以快速准确地计算出整个工程项目的工程量, 而且还可避免重复计算及漏算等差错。
2.2 优化计算程序
在建筑工程项目计算中, 墙体分部工程中含有较多的分项工程, 是整个项目预算中较大的项目。其中与墙体工程计算联系最多的是门和窗, 如果逐个计算较为麻烦, 而且也容易出错。为了提高墙体分部工程计算的速度和准确性, 应当优化计算程序, 将门窗的规格型号和分布状况及与门窗有关联的分项 (如窗栅等) 一并计算并统计出来, 将直观易算的门窗分部项目的面积和所占空间先行计算出来, 制成简表, 以备各类墙体计算时随时调用。这样就使墙体的工程量计算简单化, 在计算出各类墙体面积后扣除门窗所占部分, 即可完成整个墙体的计算。优化计算程序的方法既简便易行, 其计算精度又很高。
2.3 明确顺序分部计算
建筑工程项目预决算的编制方法多种多样, 每位编制人员所采用的方法可能不尽相同, 但不论采用哪种计算方法, 都应事先编制出计算项目的先后顺序, 并且条理清楚 (例如先按施工顺序、分部顺序、图纸顺序, 后按条、块顺序) , 对各分部项的基数主线要明确。
当编制人员采用某种顺序进行计算时, 应在明确顺序的基础上尽可能做到各分部的计算上下贯通, 每一序号上的所有项目一次计算完毕, 不应留有零散的尾项, 以防造成工程项目计算上的错误。在计算住宅楼工程量时, 因各分部、分项工程的图纸内容集中, 界限明确, 且各分部工程相互牵连的部分极少, 图纸基本上是独立的 (如楼梯、平台、阳台、雨篷、厨厕、屋面等) , 这样就可以将每一张图纸中的内容一次计算完成, 使得整个工程项目的计算范围逐渐缩小, 以减少错漏, 提高计算的准确性。
2.4 疑难部位后置
在预 (决) 算人员刚接到编制建筑工程项目图纸时, 可能会对其中的某些部位理解不深而感到难以下手, 此时最好先选择较容易的项目, 将尚未理解的部位做好标记暂缓一时。当较容易的项目全部计算完成后, 对图纸中各部位的内容已比较清楚, 原来做好标记的较难部位也已熟悉, 再对其进行计算就容易得多。这样做既可以提高计算工作的效率, 也可避免计算过程中的错漏。例如, 无论是项目实际施工还是预算分部的项目排列, 通常是将基础分部排在第一的位置, 但在预决算编制过程中, 由于基础分部为独立分部, 其端面变化较大, 计算较为复杂, 若放在前期计算难度较大, 工作效率相对较低;而如果将其放在后面计算, 因人员对整个工程项目的图纸、平面尺寸、结构形式、技术标准等均已比较熟悉, 这时再进行计算就比较容易, 其工作效率和计算的准确性也会相应提高。
3 结语
建筑工程项目的预决算及工程量的计算是一项十分复杂和繁琐的工作。为了提高工程量计算的工作效率, 减少重复劳动, 避免工作差错, 编制人员应在实践中不断总结经验, 找出工程量计算过程中较为便捷、准确的方法, 逐步掌握工程量计算的技巧, 充分应用现代计算工具和实践经验, 注重对可重复利用项目数据的收集和归类, 只有这样, 才能不断提高预算编制水平, 获得满意的计算速度和工作效果。
参考文献
[1]周和生, 尹贻林.建设项目全过程造价管理[M].天津:天津大学出版社, 2008.
高精度计算 篇8
由于直流潮流法模型简单、方程完全线性化、求解速度很快,因此国内外很多学者研究直流潮流在一些特殊场合的应用。在静态安全分析中,直流潮流可以快速筛选重载或过载支路[1,2];在调度计划中,直流潮流可以对机组组合给出的优化结果进行快速校核[3,4];在电力市场方面,直流潮流作为线性规划法的基础,在各种安全约束条件下,可以求得输电系统的可用输电能力(ATC)[5]。不仅如此,基于直流潮流计算的电网灵敏度模型,可以求得支路开断分布因子,服务于静态安全分析的辅助决策系统[6]。文献[7]在基于节点边际电价(LMP)的电力市场计算中考虑了安全约束优化潮流的目标,比较过全网交流潮流和直流潮流的差异性,结果表明直流模型非常接近于交流潮流模型。因而,对于LMP的电力市场计算,完全可以使用直流潮流来解决问题。
尽管在静态安全分析和电力市场方面大都采用直流潮流,由于忽略了电压和无功功率,有时计算结果并不理想,可能会带来校核结果的漏判、误判。这也引发了一些学者对直流潮流假设条件的详细分析,试图量化引起直流潮流的误差[8]。得出的结论是:实际运行电压应当尽可能地接近额定值。X/R的值应尽可能地高,至少在4以上,否则,标准直流潮流模型带来的误差是不可忽略的。
本文从标准直流潮流推导过程中方程的假设条件和忽略项出发,探讨了基于已知的初始稳态断面且适用于拓扑不变而仅有注入量变化的提高计算精度的类直流潮流算法。基于此方法,在算例中给出了多个IEEE标准系统和波兰实际电网的线性和随机模拟与仿真,结果表明:当网络注入量在-30%~60%范围内变化时,类直流潮流法比标准直流法在计算精度上有70%以上的改进。在一些对计算精度要求比较高的直流潮流应用场合,当满足拓扑不变且在注入量的变化区间时,类直流潮流法可以有效地提高潮流计算精度,使得直流潮流更加逼近交流潮流的精确解。
1 标准直流潮流算法简介
支路i-j的交流潮流方程表达式为:
式中:Vi和Vj分别为节点i,j的电压;θij为支路i-j两端的电压相角差;gij为支路电导;bij为支路电纳。
标准直流潮流算法与交流潮流算法相比,在计算过程中作了简化和近似[9]。简化过程如下。
1)正常运行的电网,各节点电压通常在额定电压附近,近似处理为Vi=Vj=1。
2)支路两端的电压相角差很小,有θij≈0°,则sinθij=θij,cosθij=1。
3)对超高压网络,线路电阻比电抗小得多,电阻可以忽略,即rij=0。
式(1)作以上简化后,得到标准直流潮流表达式为:
式中:θi和θj分别为节点i,j的电压相角;xij为支路电抗。
节点i的功率平衡方程为:
选取节点N为参考节点,相角θN=0°,参考节点的有功注入功率可由其他节点的注入功率唯一确定。令除去参考节点后的总节点数为n,则N=n+1。
得到实际应用中的直流潮流方程为:
式中:P和θ分别为除去参考节点后的n维有功功率和相角的列向量;B0为以1/xij为支路导纳的n×n阶节点导纳矩阵。
B0中的各元素为:
2 提高精度的类直流潮流法
首先对功率方程(式(1))进行恒等变换,支路有功功率和无功功率之间的关系可以表示为:
节点注入功率由支路功率之和表示,即
2.1 功率注入量P的修正
类直流潮流法中,用已知的交流潮流断面的支路无功功率值修正需要计算断面支路的有功功率。
由式(7)定义支路有功修正因子为Pij-cor;由式(8)定义节点注入功率修正因子为Pi-cor;定义节点注入功率向量P的修正因子为Pcor。
各修正因子分别为:
可以得到类直流潮流计算模型中,支路潮流P*ij、节点注入功率Pi*、节点注入功率向量P*分别为:
其中,有功修正因子Pij-cor,Pi-cor,Pcor都可以由已知的交流潮流断面获取,为已知常量。
2.2 矩阵参数B0的修正
由式(7)可知,在类直流潮流法中,用已知交流潮流断面的各个节点电压和相角值进行修正需要计算断面的支路导纳。
定义支路导纳比例修正因子为:
可以得到类直流潮流模型中支路导纳为:
则
进一步可得,类直流潮流方程为:
其中,类直流潮流模型中导纳矩阵B*的元素为:
由式(5)和式(19)可以看出,2种直流潮流的求解表达式在形式上达到了统一。
值得注意的是,Pij-cor,Pi-cor,Pcor,Bij-cor等修正因子在已知的一个交流潮流断面基础上只需求解一次并作为常量保存,后续基于这个相同网络拓扑而注入功率变化的断面直流潮流计算都复用这些已保存的修正因子。
3 类直流潮流法的计算流程
类直流法在计及标准直流潮流简化和忽略掉的计算项的基础上,借助一个初始断面的交流潮流精确解进行注入有功功率和支路电抗的修正,使得通过直流法求解的支路潮流值更精确。
下面以交流潮流求解的稳态断面case0为例来说明当网络拓扑不发生变化,而注入量变成P形成新的断面case1时,对应的类直流潮流法的求解过程。具体求解步骤如下。
步骤1:定义断面case0为初始稳态断面,并获取其交流潮流解。
步骤2:利用断面case0的交流潮流解来计算各个修正参数P0cor=P01-corP02-cor…P0n[]-corT,B0ij-cor=V0iV0j·(sinθ0ij)/θ0ij,并存储这些修正因子项。
步骤3:对case0拓扑下新的注入量P和导纳矩阵进行修正:Pi*=Pi-P0i-cor,B*ij=B0ij-cor/xij,求解P*,B0*。
步骤4:由P*=B0*θ求解θ。
步骤5:由Pij=P*ij+P0ij-cor=B0ij-corθij/xij+P0ij-cor求解支路潮流。
附录A图A1为类直流潮流算法的求解流程图。下面进一步从机理上对该方法进行分析。
在断面case0的基础上,P0-Pcor=B0*θ0;当注入功率P0有ΔP的功率偏差量变化时,P0+ΔP-Pcor=B0*(θ0+Δθ),由此得到:ΔP=B0*Δθ。可知,ΔP会带来Δθ的变化,由Δθ可计算支路功率的变化量。当ΔP=0时,Δθ=0,所求类直流潮流解与交流潮流解的结果相同。
4 算例分析
上面介绍的类直流潮流法基于一个已知的潮流断面,适用于网络拓扑不变而功率注入量发生变化的断面,这是类直流潮流法与标准直流潮流法的不同。
ΔP作为初始已知潮流断面P0基础上的功率注入偏量,在Δθ=(B0*)-1ΔP线性规律基础上有一定的线性求解范围。
对算例进行仿真,通过对类直流潮流法和标准直流法计算结果各个指标的比较,来说明这个可取范围的大小。
选取电网测试系统:IEEE-14,IEEE-39,IEEE-118,Polish system case2383。定义绝对误差Ea=|PDC-PAC|(PDC为直流潮流计算值;PAC为交流潮流计算值);Eamax为所有支路的绝对误差最大值;Ema为所有支路的平均绝对误差值。各个测试系统的节点、支路、拓扑信息见附录A表A1。
4.1 注入功率P的线性变化仿真
当节点注入功率P线性增长时,给出标准直流潮流法和类直流潮流法求解结果的差异性比较。
定义P=P0+(k-16)λP0,其中P0为初始断面有功功率向量,步长λ=4%,k=1,2,…,31,对应测试系统的节点注入功率P的增量ΔP在-60%~60%区间范围内线性变化。无功功率随有功功率同比例变化。
附录A图A2和A3分别为标准系统(IEEE-14,IEEE-39,IEEE-118)和Polish system采用标准直流法和类直流法求解的支路平均绝对误差比较曲线。可以看出:ΔP在-35%~60%范围波动时,类直流潮流算法的平均绝对误差都比标准直流潮流算法的小;对标准直流潮流算法而言,由于忽略了线路损耗和无功功率,当全网节点功率P增加时,支路潮流平均绝对误差也呈增长趋势,对应全网网络损耗的增加。对于类直流潮流法,由于初始交流潮流断面的使用,平均绝对误差在交流潮流断面平衡点(ΔP=0对应的断面)左右随着偏离平衡点距离的增加呈增大趋势。而且,在平衡点右侧,类直流法的误差相比标准直流法至少有50%的削减。
由附录A图A2和图A3中各个系统的标准直流潮流与类直流潮流计算的绝对误差比较可以看出,类直流潮流法的可适用求解区间为-30%~60%。
4.2 随机模拟仿真
在-30%~60%的注入量变化区间内对Polish system进行随机模拟。取随机模拟次数为1 000次。每次模拟对应生成一组个数为2 383的随机数,取值范围为-30%~60%,以表征Polish system的2 383个节点注入功率随机变化量。由于节点注入功率的随机性变化,必然有可能造成全网注入功率变化量的不平衡,在实际算例分析时,需对功率的不平衡量进行二次分配。不平衡功率以总量的形式按照全网各个可调发电机的容量按比例实施分配[10]。
附录A图A4为Polish system一次模拟对应的各个节点功率增量比例变化曲线。可以看出,节点功率增量变化在-30%~60%范围内呈随机分布,能够有效模拟各节点功率注入量的随机变化。
1)绝对误差比较
图1和图2分别为Polish system进行1 000次随机模拟对应的标准直流潮流与类直流潮流支路功率最大绝对误差和平均绝对误差变化曲线。图中黑线为线性趋势(下同)。由图1和图2可以看出:标准直流潮流的最大绝对误差区间范围为1.8~2.2,平均绝对误差区间范围为0.031~0.037;类直流潮流的最大绝对误差区间范围为0.3~0.7,平均绝对误差区间范围为0.008~0.013。类直流潮流法相比标准直流潮流法,在支路最大绝对误差和平均绝对误差上分别有75%和73%的改进。
2)相对误差比较
选用的Polish system中电压等级有400,220,110,15kV。
设重载线路功率限值为该支路电压等级基准功率的50%,考核支路功率误差限值为相应基准功率的5%。参考中国电网系统指标,Polish system各电压等级支路的基准功率、重载限值及支路功率误差限值见附录A表A2。
不同的电压等级对应不同的支路基准功率,为了方便比较类直流潮流法与标准直流法计算结果中全网支路的相对误差指标,定义相对误差Er=|PDC-PAC|/Sbranch(Sbranch为支路基准功率);Ermax为所有支路的相对误差最大值;Emr为所有支路的平均相对误差值。
图3和图4分别为Polish system进行1 000次随机模拟对应的类直流潮流与标准直流潮流支路功率最大相对误差和平均相对误差变化曲线。
由图3和图4可以看出:标准直流潮流的最大相对误差区间为23%~28%,平均相对误差区间为1.5%~1.8%;类直流潮流的最大相对误差区间为5%~10%,平均相对误差区间为0.4%~0.6%。类直流潮流法相比于标准直流潮流,在支路最大相对误差和平均相对误差上分别有了70.8%和73.5%的改进。
3)重载误差支路数比较
在静态安全分析中更关注的是重载支路的情况,下面对随机模拟的结果进行重载支路功率绝对误差的考核。定义潮流求解结果中,支路潮流大于附录A表A2中重载限值的支路为重载支路;直流潮流法与交流潮流法支路绝对误差大于附录A表A2中支路误差限值的支路为绝对误差比较大的支路,简称误差支路;既是重载支路又是误差支路的支路为重载误差支路。通过对随机模拟计算中重载误差支路数的统计,来比较类直流潮流法相比标准直流潮流法在重载支路计算准确性方面的优越性。
图5为Polish system进行1 000次随机模拟对应的类直流潮流与标准直流潮流重载误差支路数的统计比较曲线。可看出,标准直流潮流法重载误差支路数分布为38~52,而类直流潮流法重载误差支路数分布为0~8。类直流潮流法相比于标准直流潮流法,在重载误差支路数上有了95%的减少。
综上,对Polish system进行注入功率变化在区间-30%~60%范围的随机模拟,类直流潮流法相比标准潮流直流法,各种误差指标改进情况如表1所示。误差指标中的m代表重载误差支路数。
5 结语
由不同电网系统的线性功率增量模拟和对Polish system的随机模拟,确定当节点注入功率波动在-30%~60%时,采用类直流潮流法求解的结果比标准直流法在各个指标上整体有70%以上的改进。1 000次的随机模拟也表明该方法是可靠稳定的。其中,类直流潮流法和标准直流潮流法在表达形式上达到了统一,计算步骤相同,各个修正因子只需计算一次,因而在时间效率上保证了与标准直流潮流法的同步。
类直流法是基于一个已知的初始稳态潮流断面进行的计算,在对于拓扑不变仅有注入量变化的情况下是适用的。虽然类直流潮流法的使用有各种约束限制,由于具有高计算精度和快速求解速度,更好的应用场景值得探讨和挖掘。
附录见本刊网络版(http://aeps.sgepri.sgcc com.cn/aeps/ch/index.aspx)。
参考文献
[1]鲁宗相,郭勇基.配电系统可靠性及静态安全性评估中的潮流计算[J].清华大学学报:自然科学版,2000,40(3):32-35.LU Zongxiang,GUO Yongji.Flow analysis of distributionsystem used for reliability and static security evaluation[J].Journal of Tsinghua University:Science and Technology,2000,40(3):32-35.
[2]沈瑜,夏清,康重庆.发电联合转移分布因子及快速静态安全校核算法[J].电力系统自动化,2003,27(18):13-17.SHEN Yu,XIA Qing,KANG Chongqing.Generation jointlyshift distribution factor(GJSDF)and fast security analysisalgorithm in economic dispatch[J].Automation of ElectricPower Systems,2003,27(18):13-17.
[3]朱足君,李扬,陈文浩,等.节能发电调度中的有功安全校正策略[J].电力系统及其自动化学报,2008,20(6):47-50.ZHU Zujun,LI Yang,CHEN Wenhao,et al.Correctionstrategy of active power security in energy conservationgeneration dispatching[J].Proceedings of the CSU-EPSA,2008,20(6):47-50.
[4]汪洋,夏清,康重庆.考虑电网N-1闭环安全校核的最优安全发电计划[J].中国电机工程学报,2011,31(10):39-45.WANG Yang,XIA Qing,KANG Chongqing.Optimal securityconstrained generation scheduling considering closed-loop N-1security correction[J].Proceedings of the CSEE,2011,31(10):39-45.
[5]蔡晓波,吴政球,匡文凯,等.基于直流分布因子的可用输电能力灵敏度分析[J].电网技术,2006,30(18):45-48.CAI Xiaobo,WU Zhengqiu,KUANG Wenkai,et al.Study onsensitivity of available transfer capability based on DCdistribution factors[J].Power System Technology,2006,30(18):45-48.
[6]徐田,於益军,钱玉妹.能量管理系统中发电计划安全校核功能的设计[J].电力系统自动化,2006,30(10):88-92.XU Tian,YU Yijun,QIAN Yumei.Security correction designon generation plan in energy management system[J].Automation of Electric Power Systems,2006,30(10):88-92.
[7]OVERBYE T J,XU Cheng,YAN Sun.A comparison of theAC and DC power flow models for LMP calculations[C]//Proceedings of the 37th Hawaii International Conference onSystem Sciences,January 5-8,2004,Hawaii,USA.
[8]PURCHALA K,MEEUS L,VAN DOMMELEN D,et al.Usefulness of DC power flow for active power flow analysis[C]//IEEE Power Engineering Society General Meeting,June12-16,2005,Montreal,Canada:454-459.
[9]张伯明,陈寿孙,严正.高等电力网络分析[M].2版.北京:清华大学出版社,2007.