CN113591020B - 一种基于轴对称盒空间滤波的非线性系统状态估计方法 - Google Patents
一种基于轴对称盒空间滤波的非线性系统状态估计方法 Download PDFInfo
- Publication number
- CN113591020B CN113591020B CN202110848659.2A CN202110848659A CN113591020B CN 113591020 B CN113591020 B CN 113591020B CN 202110848659 A CN202110848659 A CN 202110848659A CN 113591020 B CN113591020 B CN 113591020B
- Authority
- CN
- China
- Prior art keywords
- axisymmetric
- state
- box space
- box
- space
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 105
- 238000001914 filtration Methods 0.000 title claims abstract description 18
- 238000005259 measurement Methods 0.000 claims abstract description 85
- 230000008569 process Effects 0.000 claims abstract description 65
- 238000004364 calculation method Methods 0.000 claims description 16
- 238000013016 damping Methods 0.000 claims description 9
- 239000011159 matrix material Substances 0.000 claims description 7
- 238000006243 chemical reaction Methods 0.000 claims description 6
- 238000005070 sampling Methods 0.000 claims description 5
- PXFBZOLANLWPMH-UHFFFAOYSA-N 16-Epiaffinine Natural products C1C(C2=CC=CC=C2N2)=C2C(=O)CC2C(=CC)CN(C)C1C2CO PXFBZOLANLWPMH-UHFFFAOYSA-N 0.000 claims description 3
- 238000005520 cutting process Methods 0.000 claims description 3
- 238000005312 nonlinear dynamic Methods 0.000 claims description 3
- 230000011218 segmentation Effects 0.000 claims description 3
- 230000009466 transformation Effects 0.000 claims description 3
- 230000000452 restraining effect Effects 0.000 abstract 1
- 238000004422 calculation algorithm Methods 0.000 description 5
- 238000009826 distribution Methods 0.000 description 3
- 238000001514 detection method Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008602 contraction Effects 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000003860 storage Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
- G06F17/12—Simultaneous equations, e.g. systems of linear equations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Data Mining & Analysis (AREA)
- Computational Mathematics (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Theoretical Computer Science (AREA)
- General Engineering & Computer Science (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Algebra (AREA)
- Operations Research (AREA)
- Life Sciences & Earth Sciences (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Evolutionary Biology (AREA)
- Probability & Statistics with Applications (AREA)
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
- Complex Calculations (AREA)
Abstract
本发明公开了一种基于轴对称盒空间滤波的非线性系统状态估计方法，属于状态估计技术领域。所述方法包括获取工业控制系统中的非线性弹簧‑质量‑阻尼器系统模型，获取初始区间盒子的顶点坐标，将非线性模型线性化，利用区间盒子包裹线性化过程误差，根据闵可夫斯基和运算得到k+1时刻预测区间盒子状态集，通过测量得到区间盒子测量集，并利用区间盒子测量集对k+1时刻的状态预测集进行更新，将区间盒子拆分成超平面进行约束得到两个区间盒子的交集，完成对k+1时刻系统状态的估计；解决了在未知噪声干扰下的写线性系统状态估计问题；利用区间盒子包裹误差使得估计状态集更加紧致、准确；区间盒子的拆分在保持准确性的同时减少了运算复杂度。
Description
技术领域
本发明涉及一种基于轴对称盒空间滤波的非线性系统状态估计方法，属于状态估计技术领域。
背景技术
质量-弹簧-阻尼系统是一种比较常见的机械震动系统，在人们的日常生活中也随处可见这种系统，例如汽车缓冲器就是一种可以减耗运动能量的装置，是保证驾驶员行车安全的必备装置，再者在建筑抗震加固措施中引入阻尼器，改变结构的自振特性，增加结构阻尼，吸收地震能量，降低地震作用对建筑物的影响，因此研究质量-弹簧-阻尼结构是具有现实意义的。
传统的状态估计方法在进行状态估计时，会假定噪声和扰动是已知或满足一定概率分布的变量，在此基础上进行状态估计，然而实际系统的运行环境比较复杂，扰动和噪声并非是已知或满足一定概率分布的变量，而是很难用特定的统计特性描述，所以会导致检测结果不准确。
为克服因为实际运行环境中扰动和噪声并非是已知或满足一定概率分布的变量而导致的检测结果不准确的问题，现有状态估计方法中通常会采用集员滤波的方法，即采用空间几何体，比如区间、椭球、多胞体等描述测量数据、扰动和噪声，通过弹簧系统的理论输出对噪声干扰下的系统实际弹簧位置和速度进行估计。同时在对非线性系统进行处理时，为了尽量降低计算复杂度，会进行线性化处理，而线性化处理过程往往会带来误差，利用集员滤波的方法可以很好地包裹处理该误差，实现非线性系统地状态估计。
而上述集员滤波的方法中，采用不同空间图形进行滤波时，其后续计算过程具有较大的差异，而且计算量差异也较大。比如，采用椭圆几何体，后续状态集的迭代过程中是根据最小迹椭球原理求取相关参数，并以此参数进行椭球的外包裹估计，椭圆几何体具有形状规则且计算简单的特点，但也有很明显的保守性，存在较大的空间冗余；全对称多胞体则是转换成多面体的形式，保守性虽好但是随着维度的增加，伴随着矩阵的计算，计算复杂程度也将呈指数型增加，高维度的运算时复杂度过大；超平行体虽然保守性很好，但求取最合适的超平面过程计算复杂，需要转化为多层最优化问题。因此，亟需一种新的状态估计方法以解决现有方法对非线性系统状态预测结果保守性较大且计算量大的问题。
发明内容
为了进一步提高对于非线性系统受到外界未知干扰时的状态估计的准确度，本发明提供了一种基于轴对称盒空间滤波的非线性系统状态估计方法，进行滤波扩展，同时创新性的提出了基于空间维度拆分的迭代方式，提高估计准确度的前提下也降低了计算复杂度。
一种基于轴对称盒空间滤波的非线性系统状态估计方法，所述方法包括：
步骤一：建立非线性系统模型；
步骤二：对非线性系统的状态函数和测量函数分别进行线性化处理，得到状态函数和测量函数的线性化系统模型；
步骤三：利用轴对称盒空间对状态函数线性化过程和测量函数线性化过程产生的误差分别进行包裹，获得k时刻状态函数线性化误差轴对称盒空间Bef(aef,def)和测量函数线性化误差轴对称盒空间Beh(aeh,deh)；
步骤四：获取k时刻的预测过程噪声的轴对称盒空间Bw(aw,dw)，并与状态函数线性化误差轴对称盒空间Bef(aef,def)、预测状态轴对称盒空间Bx(ax,dx)进行闵可夫斯基和运算，获得集员表示下的k+1时刻预测状态集Pk+1；
步骤五：针对k+1时刻的测量过程噪声的轴对称盒空间Bv(av,dv)、测量函数线性化误差轴对称盒空间Beh(aeh,deh)两个测量过程产生的轴对称盒空间进行闵可夫斯基和运算，根据空间状态方程得到测量状态集Sk+1，再计算得到测量状态集的轴对称盒空间
步骤六：利用测量状态集Sk+1更新预测状态集Pk+1，得到k+1时刻的状态估计集，完成对不确定噪声干扰下的非线性系统状态预测，得到k+1时刻的状态可行集Xk+1。
可选的，所述步骤六进行更新预测状态集Pk+1时，利用轴对称盒空间边界正交的性质，在进行轴对称盒空间的交集运算时，对每个轴对称盒空间进行拆分，将其拆分成n组超平面进行分组迭代更新，其中n为轴对称盒空间的空间维度。
可选的，所述非线性系统为非线性质量-弹簧-阻尼系统时，所述步骤一包括：
针对非线性质量-弹簧-阻尼系统，根据弹簧的运动特性，获取连续时间动态运动方程：
将连续时间动态运动方程离散化为空间状态表达式，得到非线性系统的非线性动态模型：
其中，函数f(xk)为非线性系统的状态函数，函数h(xk)为非线性系统的测量函数；wk为预测过程噪声，vk为测量过程噪声；
其中，xk为系统状态变量，x1,k和x2,k分别为物理量弹簧位置和弹簧速度；ΔT为系统采样时间；k0、kd、c为弹簧动态性能参数。
可选的，所述步骤二中，
状态函数的线性化系统模型为：
其中，J1,k为状态函数f的雅可比矩阵，xk为k时刻状态变量，表示当前状态估计点，/>为泰勒展开高阶项，用区间进行包裹表示为e1,k＝[e1,k,min,e1,k,max]，即为预测过程线性化误差；
测量函数的线性化系统模型为：
其中，J2,k为测量函数h的雅可比矩阵，xk+1为k+1时刻状态变量，为泰勒展开高阶项，用区间进行包裹表示为e2,k＝[e2,k,min,e2,k,max]，即为测量过程线性化误差。
可选的，所述步骤三包括：
假设非线性状态函数f(xk)是一个DC函数，则f＝[f1,f2,f3…,fi]T中每一个部分fi都是DC函数，i＝1,2,3…,j；j为任意自然数，表示可能的维度；
以凸函数gi(x)和ri(x)的差近似表示fi(x)：fi(x)＝gi(x)-ri(x)，其中，
其中α>0；
基于以上假设，根据凸函数特点构建函数：
其中，
u1和u2分别为gi(x)和ri(x)在当前状态估计点的次梯度，表达式如下：
根据上述表达式，带入到fi(x)＝gi(x)-ri(x)中得到：
fLi为fi的线性化部分；
用轴对称盒空间表示为：
其中，fL(xk)为将状态变量带入到线性化方程后的结果；
进一步将误差转换成轴对称盒空间的表达形式，得到线性化误差轴对称盒空间B(ae,de)；
其中，ae为线性化误差轴对称盒空间的中心，de为线性化误差轴对称盒空间的半径；
根据公式(13)计算k时刻状态函数线性化误差轴对称盒空间Bef(aef,def)和测量函数线性化误差轴对称盒空间Beh(aeh,deh)。
可选的，所述步骤四包括：
利用传感器测量得到预测过程噪声的边界后，根据公式(13)计算得到k时刻的预测过程噪声的轴对称盒空间Bw(aw,dw)；
根据公式(14)求得当前系统的预测状态集Pk+1：
依据线性化后得到的线性系统进行仿射变换，将当前k时刻的轴对称盒空间状态集Bx(ax,dx)转换成k+1时刻的轴对称盒空间状态集即
利用公式(15)计算得到每一时刻的线性化误差轴对称盒空间，再与测量过程噪声轴对称盒空间Bv(av,dv)进行闵可夫斯基和，得到误差轴对称盒空间，其中，初始0时刻的Bx(ax,dx)为传感器测量得到。
可选的，所述步骤五包括：
对非线性的测量函数进行线性化，得到：
依据公式(13)得到测量函数线性化误差轴对称盒空间Beh(aeh,deh)，再通过闵可夫斯基和运算，与测量过程噪声轴对称盒空间Bv(av,dv)得到测量过程误差轴对称盒空间Bh(ah,dh)，其中，测量过程噪声轴对称盒空间Bv(av,dv)由传感器实际测量得到的边界计算得到；
得到测量过程的误差轴对称盒空间Bh(ah,dh)后，即可根据公式(2)的测量函数计算得到测量状态集Sk+1，再对测量状态进行轴对称盒空间表达；k+1时刻测量集轴对称盒空间表达式为：
其中，Rn表示n维实数集；yk+1表示输出变量；表示线性化的测量函数；表示测量函数线性化一阶误差；/>表示测量函数线性化高阶误差；/>表示轴对称盒空间中心；/>表示轴对称盒空间对角阵；m表示轴对称盒空间维度。
可选的，所述步骤六包括：
针对轴对称盒空间
B(a，d)＝{x∈Rn|x＝a+diag{d}m，||m||∞≤1} (18)
利用轴对称盒空间边界正交的性质，在进行轴对称盒空间的交集运算时，对每个轴对称盒空间进行拆分，将其拆分成n组超平面，其中n为轴对称盒空间的空间维度；
超平面的表达式如下：
其中S′k+1为轴对称盒空间测量集，为超平面测量集；
将一个n维空间盒子Sk+1转换成n组超平面，利用每组超平面对另一个相交的轴对称盒空间Pk+1进行切分；n组超平面在数值层面上等效为2n个约束条件，再将超平面转换2n个约束条件，进行线性规划；
通过线性规划，在利用所有超平面完成切分后，即得到轴对称盒空间交集的新区间边界[βi,min,βi,max]，重新计算得到这一时刻可行集轴对称盒空间即两个轴对称盒空间取交集运算后的新的轴对称盒空间，其中，
最终得到k+1时刻的状态可行集Xk+1＝Pk+1∩Sk+1，其中，Pk+1为预测状态集，Sk+1为测量状态集。
可选的，所述非线性系统为非线性质量-弹簧-阻尼系统时，系统采样时间ΔT＝0.1s；弹簧动态性能参数k0＝1.5，kd＝3，c＝1.24。
可选的，所述预测过程噪声wk和测量过程噪声vk分别为满足条件wk∈[-0.002,0.002]和vk∈[-0.001,0.001]的均匀噪声。
本发明有益效果是：
通过利用轴对称盒空间包裹线性化误差的方法，解决了利用扩展集员算法进行非线性状态估计时保守性过差，计算复杂度过高的问题，达到了利用轴对称盒空间包裹线性化误差进行状态估计时提高算法保守性的效果，同时根据空间维度拆分盒空间，转换成超平面进行分组迭代，每组的迭代就转换成了两组约束条件的线性规划，保证准确性的同时简化了运算。在保证实现非线性系统状态估计的功能的同时，减少了保守性，降低了计算复杂度。
附图说明
为了更清楚地说明本发明实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1是本发明一个实施例中公开的一种基于轴对称盒空间滤波的非线性系统状态估计方法的流程图。
图2是本发明一个实施例中针对具体的非线性质量-弹簧-阻尼系统进行状态估计的弹簧位置状态变量的估计区间变化图。
图3是本发明一个实施例中针对具体的非线性质量-弹簧-阻尼系统进行状态估计的弹簧速度状态变量的估计区间变化图。
图4是本发明一个实施例中通过滤波得到状态估计集员空间表示仿真图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明实施方式作进一步地详细描述。
实施例一：
本实施例提供一种基于轴对称盒空间滤波的非线性系统状态估计方法，质量-弹簧-阻尼系统作为一种典型的非线性系统，本实施例以应用于质量-弹簧-阻尼系统(后简称为弹簧系统)中进行状态估计为例进行说明，参见图1，所述方法包括：
步骤一：建立非线性系统模型；
针对非线性质量-弹簧-阻尼系统，根据弹簧的运动特性，获取连续时间动态运动方程：
将连续时间动态运动方程离散化为空间状态表达式，得到非线性系统的非线性动态模型：
其中，函数f(xk)为非线性系统的状态函数，函数h(xk)为非线性系统的测量函数；wk为预测过程噪声，vk为测量过程噪声；
弹簧系统相关参数如下：x1,k和x2,k分别为物理量弹簧位置和弹簧速度。ΔT＝0.1s，为系统采样时间；k0＝1.5，kd＝3，c＝1.24，为弹簧动态性能参数。其他非线性系统中，包含对应的系统动态性能参数。
本实施例中，弹簧系统初始状态设为x0＝[0.2,0.2]T，同时设初始状态可行集为其中，/>d0＝diag{0.01,0.01}；预测过程噪声wk和测量过程噪声vk分别为有界噪声，本实施例中二者分别满足条件wk∈[-0.002,0.002]和vk∈[-0.001,0.001]的均匀噪声，xk为系统状态变量，其中，k表示时刻。
步骤二：对非线性系统的状态函数进行线性化处理，得到状态函数的线性化系统模型：
其中，J1,k为状态函数f的雅可比矩阵，为泰勒展开高阶项，用区间进行包裹表示为e1,k＝[e1,k,min,e1,k,max]，即为预测过程线性化误差。
同理，对非线性系统的测量函数进行线性化处理，得到测量函数的线性化系统模型：
其中，J2,k为测量函数h的雅可比矩阵，为泰勒展开高阶项，用区间进行包裹表示为e2,k＝[e2,k,min,e2,k,max]，即为测量过程线性化误差。
本实施例中，由于弹簧系统的非线性动态模型中测量函数h(xk)本身即为线性化函数，因此直接给出了测量函数的线性化系统模型；对于其他的非线性系统的测量函数，直接进行线性化处理得到对应的测量函数的线性化系统模型。
步骤三：利用轴对称盒空间对状态函数线性化过程和测量函数线性化过程产生的误差分别进行包裹，获得k时刻状态函数线性化误差轴对称盒空间Bef(aef,def)和测量函数线性化误差轴对称盒空间Beh(aeh,deh)。
假设非线性状态函数f(xk)是一个DC函数(每个连续函数都可以用两个凸函数的差来近似)，则认为f＝[f1,f2,f3…,fi]T函数中的每一个部分fi(i＝1,2,3…,j)都是DC函数，其中，j为任意自然数，表示可能的维度。即存在凸函数gi(x)和ri(x)满足关系
fi(x)＝gi(x)-ri(x)，其中，
其中α>0，为任意参数。
基于以上假设，根据凸函数特点可以构建函数：
其中，u1和u2分别为gi(x)和ri(x)在当前状态估计点/>的次梯度，表达式如下：
根据上述表达式，带入到误差的关系式即fi(x)＝gi(x)-ri(x)中得到：
其中，fLi为fi的线性化部分；
在集员滤波算法中，状态量为轴对称盒空间表示，例如在k时刻，xk∈VB，VB为状态的轴对称盒空间。因此最终推导得到：
fL(xk)为将状态变量带入到线性化方程后的结果。
进一步将误差转换成轴对称盒空间的表达形式，得到线性化误差轴对称盒空间B(ae,de)。其中：
ae为线性化误差轴对称盒空间的中心，de为线性化误差轴对称盒空间的半径；公式(13)为轴对称空间盒子求取的通用公式，分别根据公式(13)求得非线性状态函数f(xk)和测量函数h(xk+1)的线性化过程产生的状态函数线性化误差轴对称盒空间Bef(aef,def)和测量函数线性化误差轴对称盒空间Beh(aeh,deh)。
步骤四：获得k时刻的预测过程噪声的轴对称盒空间Bw(aw,dw)，并与状态函数线性化误差轴对称盒空间Bef(aef,def)、预测状态轴对称盒空间Bx(ax,dx)进行闵可夫斯基和运算，获得集员表示下的k+1时刻预测状态集Pk+1；
将当前的预测状态轴对称盒空间集Bx(ax,dx)，与状态函数线性化过程所产生的状态函数线性化误差的轴对称盒空间Bef(aef,def)和预测过程噪声的轴对称盒空间Bw(aw,dw)进行闵可夫斯基和运算，可以求得当前系统的预测状态集Pk+1，其中，预测过程噪声的轴对称盒空间Bw(aw,dw)由传感器测量得到未知噪声边界后，利用公式(13)计算得到；
根据公式(14)求得当前系统的预测状态集Pk+1：
依据线性化后得到的线性系统进行仿射变换，将当前k时刻的轴对称盒空间状态集Bx(ax,dx)转换成k+1时刻的轴对称盒空间状态集即
利用公式(15)计算得到每一时刻的线性化误差轴对称盒空间集员，再与测量过程噪声轴对称盒空间Bv(av,dv)进行闵可夫斯基和，得到误差轴对称盒空间，其中，初始0时刻的Bx(ax,dx)为传感器测量得到。
采用上述状态预测过程，可以得到非线性系统每一时刻的系统状态变量的区间，并以轴对称盒空间的形式进行表示。
步骤五：针对k+1时刻的测量过程噪声轴对称盒空间Bv(av,dv)、测量函数线性化误差轴对称盒空间Beh(aeh,deh)等两个测量过程产生的轴对称盒空间进行闵可夫斯基和运算，根据空间状态方程得到测量状态集Sk+1，再计算得到测量集轴对称盒空间
借鉴卡尔曼滤波思想对非线性的测量函数进行线性化，得到
并依据式(13)可以得到测量函数线性化误差轴对称盒空间Beh(aeh,deh)，再通过闵可夫斯基和运算，与测量过程噪声轴对称盒空间Bv(av,dv)得到测量过程误差轴对称盒空间Bh(ah,dh)，其中，测量过程噪声轴对称盒空间Bv(av,dv)由传感器实际测量得到的边界计算得到。
得到测量过程的误差轴对称盒空间Bh(ah,dh)后，即可根据公式(2)计算得到测量状态集Sk+1，再对测量状态进行轴对称盒空间表达。k+1时刻测量集轴对称盒空间表达式为：
步骤六：利用测量状态集Sk+1更新预测状态集Pk+1，得到k+1时刻的状态估计集，完成对不确定噪声干扰下的非线性系统状态预测，得到状态可行集Xk+1。
在传统的集员滤波过程中，测量更新这一收缩步是通过不同集员间的取交集运算得到的，在不同的空间结构的表示下，存在着不同的交集计算方式。为了尽可能减少交集求取过程中误差的干扰，本申请所提出轴对称盒空间求交集策略，利用轴对称盒空间2n条边界在空间中彼此正交的性质，在进行闵可夫斯基和运算时不存在冗余项，可以得到准确的轴对称盒空间交集。
针对轴对称盒空间
B(a，d)＝{x∈Rn|x＝a+diag{d}m，||m||∞≤1} (18)
利用轴对称盒空间边界正交的性质，在进行轴对称盒空间的交集运算时，对每个轴对称盒空间进行拆分，将其拆分成n组超平面，其中n为轴对称盒空间的空间维度。此时,超平面的表达式如下：
其中S′k+1为轴对称盒空间测量集，为超平面测量集。
将一个n维空间盒子Sk+1转换成n组超平面，利用每组超平面对另一个相交的轴对称盒空间Pk+1进行切分。n组超平面在数值层面上可以理解为2n个约束条件，再将超平面转换2n个约束条件，进行线性规划。
通过线性规划，在利用所有超平面完成切分后，即得到轴对称盒空间交集的新区间边界[βi,min,βi,max]，重新计算得到这一时刻可行集轴对称盒空间即两个轴对称盒空间取交集运算后的新的轴对称盒空间，其中，
最终得到k+1时刻的状态可行集Xk+1＝Pk+1∩Sk+1，其中，Pk+1为预测状态集，Sk+1为测量状态集。
根据上述基于轴对称盒空间的状态估计方法，对构造的非线性质量-弹簧-阻尼系统模型进行状态估计。
在本实例中，在预定时间范围内，执行步骤一至步骤六后，得到状态量的估计盒空间，实现真实值的估计。图2和图3分别为采用不同方法对弹簧系统中弹簧位置状态变量的估计区间变化图和弹簧速度状态变量的估计区间变化图。
从图2和图3可以看出，现有的两种估计方法和本申请所提出的估计方法都能实现对非线性系统的状态估计，状态真值位于估计区间范围内，但是本申请所提出的估计算法估计的到状态值区间更小，即状态集员包裹的保守性更好，更加准确。
图4为多维状态量的空间表示，可以看到真实值是位于估计得到的状态集轴对称盒空间内的，但相比于现有的算法，轴对称盒空间的包裹性更好，更加接近于真值，同时在收缩速度上也更快，在较短时间内可以得到一个较好的估计状态集，而其他现有方法要进行多步迭代后才能取得大致相同的效果，所以本发明所提出的基于轴对称盒空间的集员估计方法在准确性，包裹性，和运算速度上都具备更好的使用性能。
本发明实施例中的部分步骤，可以利用软件实现，相应的软件程序可以存储在可读取的存储介质中，如光盘或硬盘等。
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
Claims (10)
1.一种基于轴对称盒空间滤波的非线性系统状态估计方法，其特征在于，所述方法包括：
步骤一：建立非线性系统模型；
步骤二：对非线性系统的状态函数和测量函数分别进行线性化处理，得到状态函数和测量函数的线性化系统模型；
步骤三：利用轴对称盒空间对状态函数线性化过程和测量函数线性化过程产生的误差分别进行包裹，获得k时刻状态函数线性化误差轴对称盒空间Bef(aef,def)和测量函数线性化误差轴对称盒空间Beh(aeh,deh)；
步骤四：获取k时刻的预测过程噪声的轴对称盒空间Bw(aw,dw)，并与状态函数线性化误差轴对称盒空间Bef(aef,def)、预测状态轴对称盒空间Bx(ax,dx)进行闵可夫斯基和运算，获得集员表示下的k+1时刻预测状态集Pk+1；
步骤五：针对k时刻的测量过程噪声的轴对称盒空间Bv(av,dv)、测量函数线性化误差轴对称盒空间Beh(aeh,deh)两个测量过程产生的轴对称盒空间进行闵可夫斯基和运算，根据空间状态方程得到测量状态集Sk+1，再计算得到测量状态集的轴对称盒空间
步骤六：利用测量状态集Sk+1更新预测状态集Pk+1，得到k+1时刻的状态估计集，完成对不确定噪声干扰下的非线性系统状态预测，得到k+1时刻的状态可行集Xk+1；
所述针对轴对称盒空间表示为：
B(a,d)＝{x∈Rn∣x＝a+diag{d}m,||m||∞≤1} (18)
其中，a是轴对称盒空间的中心，d表示轴对称盒空间的半径，利用轴对称盒空间边界正交的性质，在进行轴对称盒空间的交集运算时，对每个轴对称盒空间进行拆分，将其拆分成n组超平面，其中n为轴对称盒空间的空间维度。
2.根据权利要求1所述的方法，其特征在于，所述步骤六进行更新预测状态集Pk+1时，利用轴对称盒空间边界正交的性质，在进行轴对称盒空间的交集运算时，对每个轴对称盒空间进行拆分，将其拆分成n组超平面进行分组迭代更新，其中n为轴对称盒空间的空间维度。
3.根据权利要求2所述的方法，其特征在于，所述非线性系统为非线性质量-弹簧-阻尼系统时，所述步骤一包括：
针对非线性质量-弹簧-阻尼系统，根据弹簧的运动特性，获取连续时间动态运动方程：
将连续时间动态运动方程离散化为空间状态表达式，得到非线性系统的非线性动态模型：
其中，函数f(xk)为非线性系统的状态函数，函数h(xk)为非线性系统的测量函数；wk为预测过程噪声，vk为测量过程噪声；
其中，xk为系统状态变量，x1,k和x2,k分别为物理量弹簧位置和弹簧速度；ΔT为系统采样时间；k0、kd、c为弹簧动态性能参数。
4.根据权利要求3所述的方法，其特征在于，所述步骤二中，
状态函数的线性化系统模型为：
其中，J1,k为状态函数f的雅可比矩阵，xk为k时刻状态变量，表示当前状态估计点，为泰勒展开高阶项，用区间进行包裹表示为e1,k＝[e1,k,min,e1,k,max]，即为预测过程线性化误差；
测量函数的线性化系统模型为：
其中，J2,k为测量函数h的雅可比矩阵，xk+1为k+1时刻状态变量，为泰勒展开高阶项，用区间进行包裹表示为e2,k＝[e2,k,min,e2,k,max]，即为测量过程线性化误差。
5.根据权利要求4所述的方法，其特征在于，所述步骤三包括：
假设非线性状态函数f(xk)是一个DC函数，则f＝[f1,f2,f3...,fi]T中每一个部分fi都是DC函数，i＝1,2,3…,j；j为任意自然数，表示可能的维度；
以凸函数gi(x)和ri(x)的差近似表示fi(x)：fi(x)＝gi(x)-ri(x)，其中，
其中α>0，为任意参数；
基于以上假设，根据凸函数特点构建函数：
其中，
u1和u2分别为gi(x)和ri(x)在当前状态估计点的次梯度，表达式如下：
根据上述表达式，带入到fi(x)＝gi(x)-ri(x)中得到：
fLi为fi的线性化部分；
用轴对称盒空间表示为：
其中，fL(xk)为将状态变量带入到线性化方程后的结果；
进一步将误差转换成轴对称盒空间的表达形式，得到线性化误差轴对称盒空间B(ae,de)；
其中，ae为线性化误差轴对称盒空间的中心，de为线性化误差轴对称盒空间的半径；
根据公式(13)计算k时刻状态函数线性化误差轴对称盒空间Bef(aef,def)和测量函数线性化误差轴对称盒空间Beh(aeh,deh)。
6.根据权利要求5所述的方法，其特征在于，所述步骤四包括：
利用传感器测量得到预测过程噪声的边界后，根据公式(13)计算得到k时刻的预测过程噪声的轴对称盒空间Bw(aw,dw)；
根据公式(14)求得当前系统的预测状态集Pk+1：
依据线性化后得到的线性系统进行仿射变换，将当前k时刻的轴对称盒空间状态集Bx(ax,dx)转换成k+1时刻的轴对称盒空间状态集即
利用公式(15)计算得到每一时刻的线性化误差轴对称盒空间，再与测量过程噪声轴对称盒空间Bv(av,dv)进行闵可夫斯基和，得到误差轴对称盒空间，其中，初始0时刻的Bx(ax,dx)为传感器测量得到。
7.根据权利要求6所述的方法，其特征在于，所述步骤五包括：
对非线性的测量函数进行线性化，得到：
依据公式(13)得到测量函数线性化误差轴对称盒空间Beh(aeh,deh)，再通过闵可夫斯基和运算，与测量过程噪声轴对称盒空间Bv(av,dv)得到测量过程误差轴对称盒空间Bh(ah,dh)，其中，测量过程噪声轴对称盒空间Bv(av,dv)由传感器实际测量得到的边界计算得到；
得到测量过程的误差轴对称盒空间Bh(ah,dh)后，即可根据公式(2)计算得到测量状态集Sk+1，再对测量状态进行轴对称盒空间表达；k+1时刻测量集轴对称盒空间表达式为：
其中，Rn表示n维实数集；yk+1表示输出变量；表示线性化的测量函数；/>表示测量函数线性化一阶误差；/>表示测量函数线性化高阶误差；/>表示轴对称盒空间中心；/>表示轴对称盒空间对角阵；m表示轴对称盒空间维度。
8.根据权利要求7所述的方法，其特征在于，所述步骤六包括：
针对轴对称盒空间
B(a,d)＝{x∈Rn∣x＝a+diag{d}m,||m||∞≤1} (18)
利用轴对称盒空间边界正交的性质，在进行轴对称盒空间的交集运算时，对每个轴对称盒空间进行拆分，将其拆分成n组超平面，其中n为轴对称盒空间的空间维度；
超平面的表达式如下：
其中S′k+1为轴对称盒空间测量集，为超平面测量集；
将一个n维空间盒子Sk+1转换成n组超平面，利用每组超平面对另一个相交的轴对称盒空间Pk+1进行切分；n组超平面在数值层面上等效为2n个约束条件，再将超平面转换2n个约束条件，进行线性规划；
通过线性规划，在利用所有超平面完成切分后，即得到轴对称盒空间交集的新区间边界[βi,min,βi,max]，重新计算得到这一时刻可行集轴对称盒空间即两个轴对称盒空间取交集运算后的新的轴对称盒空间，其中，
最终得到k+1时刻的状态可行集Xk+1＝Pk+1∩Sk+1，其中，Pk+1为预测状态集，Sk+1为测量状态集。
9.根据权利要求8所述的方法，其特征在于，所述非线性系统为非线性质量-弹簧-阻尼系统时，系统采样时间ΔT＝0.1s；弹簧动态性能参数k0＝1.5，kd＝3，c＝1.24。
10.根据权利要求9所述的方法，其特征在于，所述预测过程噪声wk和测量过程噪声vk分别为满足条件wk∈[-0.002,0.002]和vk∈[-0.001,0.001]的均匀噪声。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110848659.2A CN113591020B (zh) | 2021-07-23 | 2021-07-23 | 一种基于轴对称盒空间滤波的非线性系统状态估计方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110848659.2A CN113591020B (zh) | 2021-07-23 | 2021-07-23 | 一种基于轴对称盒空间滤波的非线性系统状态估计方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113591020A CN113591020A (zh) | 2021-11-02 |
| CN113591020B true CN113591020B (zh) | 2024-03-01 |
Family
ID=78250371
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110848659.2A Active CN113591020B (zh) | 2021-07-23 | 2021-07-23 | 一种基于轴对称盒空间滤波的非线性系统状态估计方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113591020B (zh) |
Families Citing this family (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114563104B (zh) * | 2022-02-14 | 2022-11-08 | 江南大学 | 一种基于非线性电热模型的动力电池双温度估计方法 |
| CN117194866B (zh) * | 2023-09-11 | 2024-03-26 | 哈尔滨理工大学 | 一种基于质量弹簧阻尼系统的分布式滤波方法 |
| CN117910201B (zh) * | 2023-10-31 | 2024-07-16 | 湖南盛鼎科技发展有限责任公司 | 用于雷达数据的多胞体双滤波状态估计算法 |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103294931A (zh) * | 2013-06-28 | 2013-09-11 | 上海交通大学 | 基于改进的非线性鲁棒滤波算法的系统状态估计方法 |
| US9002678B1 (en) * | 2014-01-10 | 2015-04-07 | King Fahd University Of Petroleum And Minerals | Unified approach to detection and isolation of parametric faults using a kalman filter residual-based approach |
| CN104795811A (zh) * | 2015-03-17 | 2015-07-22 | 河海大学 | 一种电力系统区间状态估计方法 |
| CN109597864A (zh) * | 2018-11-13 | 2019-04-09 | 华中科技大学 | 椭球边界卡尔曼滤波的即时定位与地图构建方法及系统 |
| CN109612470A (zh) * | 2019-01-14 | 2019-04-12 | 广东工业大学 | 一种基于模糊容积卡尔曼滤波的单站无源导航方法 |
| CN111985093A (zh) * | 2020-08-03 | 2020-11-24 | 哈尔滨工程大学 | 一种带噪声估计器的自适应无迹卡尔曼滤波状态估计方法 |
| CN112305418A (zh) * | 2020-10-13 | 2021-02-02 | 江南大学 | 一种基于混合噪声双重滤波的电机系统故障诊断方法 |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US11501193B2 (en) * | 2017-12-18 | 2022-11-15 | Mitsubishi Electric Research Laboratories, Inc. | Model-based control under uncertainty |
-
2021
- 2021-07-23 CN CN202110848659.2A patent/CN113591020B/zh active Active
Patent Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103294931A (zh) * | 2013-06-28 | 2013-09-11 | 上海交通大学 | 基于改进的非线性鲁棒滤波算法的系统状态估计方法 |
| US9002678B1 (en) * | 2014-01-10 | 2015-04-07 | King Fahd University Of Petroleum And Minerals | Unified approach to detection and isolation of parametric faults using a kalman filter residual-based approach |
| CN104795811A (zh) * | 2015-03-17 | 2015-07-22 | 河海大学 | 一种电力系统区间状态估计方法 |
| CN109597864A (zh) * | 2018-11-13 | 2019-04-09 | 华中科技大学 | 椭球边界卡尔曼滤波的即时定位与地图构建方法及系统 |
| CN109612470A (zh) * | 2019-01-14 | 2019-04-12 | 广东工业大学 | 一种基于模糊容积卡尔曼滤波的单站无源导航方法 |
| CN111985093A (zh) * | 2020-08-03 | 2020-11-24 | 哈尔滨工程大学 | 一种带噪声估计器的自适应无迹卡尔曼滤波状态估计方法 |
| CN112305418A (zh) * | 2020-10-13 | 2021-02-02 | 江南大学 | 一种基于混合噪声双重滤波的电机系统故障诊断方法 |
Non-Patent Citations (4)
| Title |
|---|
| Efficient extended cubature Kalman filtering for nonlinear target tracking;HE R K等;International Journal of Systems Science;第52卷(第2期);392 – 406 * |
| 一种基于凸多面体的集员滤波故障诊断方法;沈艳霞等;控制与决策;第33卷(第1期);150 – 156 * |
| 基于集员估计的线性MIMO系统故障检测方法;张峻峰等;上海应用技术学院学报(自然科学版);第15卷(第2期);149 – 152 * |
| 非线性系统中心差分集员估计方法;沈强等;控制理论与应用;第36卷(第8期);1239 – 1249 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113591020A (zh) | 2021-11-02 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN113591020B (zh) | 一种基于轴对称盒空间滤波的非线性系统状态估计方法 | |
| Ragwitz et al. | Markov models from data by simple nonlinear time series predictors in delay embedding spaces | |
| CN108304612B (zh) | 基于噪声补偿的迭代平方根ckf的汽车雷达目标跟踪方法 | |
| Tzeng | Design of fuzzy wavelet neural networks using the GA approach for function approximation and system identification | |
| KR20160021209A (ko) | 동적 모델 식별을 모니터링 및 가변 구조 또는 가변 동작 조건을 가진 동적 기기의 제어 방법 및 시스템 | |
| CN106767780B (zh) | 基于Chebyshev多项式插值逼近的扩展椭球集员滤波方法 | |
| CN114428803B (zh) | 空压站运行优化方法、系统、存储介质及终端 | |
| CN101763087B (zh) | 一种基于非线性共轭梯度法的工业过程动态优化系统及方法 | |
| CN111191186B (zh) | 一种生产车间移动机器人位姿定位的多胞体滤波方法 | |
| CN113946990A (zh) | Sto-btcn热误差预测模型建模方法及其迁移学习方法 | |
| CN112883508A (zh) | 一种基于平行空间滤波的弹簧阻尼系统状态估计方法 | |
| CN113052265A (zh) | 基于特征选择的移动对象轨迹简化算法 | |
| Miao et al. | Fault diagnosis of wheeled robot based on prior knowledge and spatial-temporal difference graph convolutional network | |
| CN111761583B (zh) | 一种智能机器人运动定位方法及系统 | |
| Sun et al. | Belief space planning for mobile robots with range sensors using iLQG | |
| CN115114985A (zh) | 一种基于集合理论的传感器系统分布式融合方法 | |
| CN116067600A (zh) | 幕墙工作模态参数识别、故障诊断与健康状态监测方法 | |
| CN113486523B (zh) | 一种线性变参数振动系统全局辨识方法 | |
| CN113848876B (zh) | 一种低通讯和计算消耗的多agv定距巡航系统故障检测方法 | |
| Tian et al. | An improved method for NURBS surface based on particle swarm optimization BP neural network | |
| CN113031510B (zh) | 面向复杂多轴系统的b样条曲线计算加速方法 | |
| CN102801426B (zh) | 一种时序数据拟合及压缩方法 | |
| CN113608442A (zh) | 基于特征函数的非线性状态模型系统的状态估计方法 | |
| Ma et al. | Research on robot SLAM of RBPF improved with weight accumulation | |
| Zheng et al. | Neurodynamics Adaptive Reward and Action for Hand-to-Eye Calibration With Deep Reinforcement Learning |
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 |