CN105631556A - 一种非恒定的条件非线性最优参数扰动计算方法 - Google Patents

Abstract

本发明公开了一种非恒定的条件非线性最优参数扰动计算方法，包括如下步骤：假设一个预报系统表达为：X(t)＝M(p)(X₀)；假设一个随时间、空间变化的参数扰动记作p’；添加参数扰动之后上式变成：X’(t)＝M(p+p’)(X₀)；利用上述两个算式构造一个能够衡量参数扰动对预报误差影响大小的目标函数J；按照梯度定义估算目标函数J关于p’的各时刻、空间上的分量的梯度；再利用优化算法，计算出能使目标函数J最大的参数扰动p’^*。本发明可以计算出随时间和空间变化的最优参数扰动，可以用于分析不同时段、不同区域的参数扰动对预报系统的影响。

Description

一种非恒定的条件非线性最优参数扰动计算方法

技术领域

本发明涉及参数扰动计算领域，具体涉及一种非恒定的条件非线性最优参数扰动计算方法。

背景技术

在预报领域中(如天气预报、海况预报)，由于存在初始误差和参数误差，所以会导致预报结果出现误差。要提高预报结果的准确性就有必要研究初始误差和参数误差，在各种类型的初始误差和参数误差中，最典型的就是能导致预报误差最大的一类，也就是最优的初始扰动和参数扰动。要找到这样的最优扰动，就需要一定的计算方法。

目前已公开的类似的方法是条件非线性最优扰动方法(ConditionalNonlinearOptimalPerturbation,简称CNOP)，该方法有两方面内容：一是用于计算最优初始扰动(ConditionalNonlinearOptimalPerturbation-Initial,简称CNOP-I)；另一个则是用于计算最优参数扰动(ConditionalNonlinearOptimalPerturbation-Parameter,简称CNOP-P),但是其将参数扰动视作恒定进行优化计算，所得的最优参数扰动也是恒定的(不随时间、空间变化)，显然这种方法在某些参数扰动变化的系统中是不够完善的。

发明内容

为解决上述问题，本发明提供了一种非恒定的条件非线性最优参数扰动计算方法，将参数扰动作为变化的量来进行优化计算，计算出的最优参数扰动也会随时间或空间变化。

为实现上述目的，本发明采取的技术方案为：

一种非恒定的条件非线性最优参数扰动计算方法，包括如下步骤：

S1、假设一个预报系统表达为：

X(t)＝M(p)(X₀)(1)

式中，X为预报变量，X＝(x₁,x₂,…,x_n),n表示变量个数，X₀为初始状态，X(t)为预报结果，p＝(p₁,p₂…..p_m)为系统中的参数(共m个)，M(p)为参数为p时的非线性传播算子；p_i(i＝1,2,….,m)可以记作p_i(j,k)，j表示不同空间位置，k表示不同的时间点；p_i(j,k)可以是一个常数，也可以随时间、空间变化；

S2、假设一个随时间、空间变化的参数扰动记作p’，p’＝(p₁’,p₂’,…..,p_m’)。p_i’(i＝1,2,….,m)是随时间和空间变化的量，那么p_i’可以记作p_i’(j,k)，j表示空间位置，k表示时间点；

S3、添加参数扰动之后(1)变成：

X’(t)＝M(p+p’)(X₀)(2)

S4、利用(1)和(2)构造一个能够衡量参数扰动对预报误差影响大小的目标函数J，那么J就是p’的函数，即:

J＝J(p’)(3)

S5、根据梯度定义，估算目标函数J对p_i’(j,k)各分量的梯度函数G；

S6、设置参数扰动约束条件；

S7、在优化算法中输入多组不同的参数扰动初猜值，扰动约束条件，目标函数，目标函数对扰动的梯度函数，结果输出条件，即可计算出能使目标函数最大的参数扰动p’^*,p’^*即是最优参数扰动。

优选地，通过如下步骤计算p’^*：

S11、同S1；

S12、同S2；

S13、同S3；

S14、同S4,具体目标函数J为：

J＝-(||X(t)’-X(t)||₂)²/2，(4)

其中||||₂表示L2范数；

S15、根据梯度定义，估算目标函数J对p_i’(j,k)各分量的梯度函数G；具体为，令hp_i’(j,k)＝p_i’(j,k)+h,h取一个适当的小值，令hp_i’的其它项全等于p_i’的对应项，目标函数对p_i’(j,k)的偏导数即为：

J|p_i’(j,k)＝(J(hp_i’)-J(p_i’))/h(5)计算出J对p_i’(j,k)的各分量的偏导数，即计算出J关于p_i’(j,k)各分量的梯度G；

S16、同S6

S17、同S7。

优选地，通过如下步骤计算p’^*：

S21、同权利要求1,具体将p设为0，公式(1)具体变为：

X(t)＝M(0)(X₀)(6)

S22、同S2；

S23、同S3，公式(2)具体变为：

X’(t)＝M(p’)(X₀)(7)

S24、同S14；

S25、同S15；

S26、同S6,具体约束条件设置为：||p_i’||₂<＝r_i,或p_i’(j,k)∈[a_i,b_i](i＝1,2,…,m)，其中r_i、a_i、b_i取适当的值；

S27、同S7，具体使用的优化算法为SpectralProjectedGrad(SPG)算法。

本发明具有以下有益效果：

可以计算出随时间和空间变化的最优参数扰动，比上之前的CNOP方法要更精细，可以用于分析不同时段、不同区域的参数扰动对预报系统的影响。

具体实施方式

为了使本发明的目的及优点更加清楚明白，以下结合实施例对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

本发明实施例提供了一种非恒定的条件非线性最优参数扰动计算方法，包括如下步骤：

S1、假设一个预报系统表达为：

X(t)＝M(p)(X₀)(1)

式中，X为预报变量，X＝(x₁,x₂,…,x_n),n表示变量个数，X₀为初始状态，X(t)为预报结果，p＝(p₁,p₂…..p_m)为系统中的参数(共m个)，M(p)为参数为p时的非线性传播算子；p_i(i＝1,2,….,m)也可以记作p_i(j,k)，j表示不同空间位置，k表示不同的时间点；p_i(j,k)可以是一个常数，也可以随时间、空间变化；

S2、假设一个随时间变化的参数扰动记作p’，p’＝(p₁’,p₂’,…..,p_m’)。p_i’(i＝1,2,….,m)是随时间和空间变化的量，那么p_i’可以记作p_i’(j,k)，j表示空间位置，k表示时间点；

S3、添加参数扰动之后(1)变成：

X’(t)＝M(p+p’)(X₀)(2)

S4、利用(1)和(2)构造一个能够衡量参数扰动对预报误差影响大小的目标函数J，那么J就是p’的函数，即:

J＝J(p’)(3)

S5、根据梯度定义，估算目标函数J对p_i’(j,k)各分量的梯度函数G；

S6、设置参数扰动约束条件；

S7、在优化算法中输入多组不同的参数扰动初猜值，扰动约束条件，目标函数，目标函数对扰动的梯度函数，结果输出条件，即可计算出能使目标函数最大的参数扰动p’^*,p’^*即是最优参数扰动。

通过如下步骤计算p’^*：向优化算法中输入多组不同的参数扰动初猜值，扰动约束条件，目标函数，目标函数对扰动的梯度函数，结果输出条件，即可计算出最优参数扰动。

实施例

以下以洛伦兹方程为例，示范本发明。

洛伦兹方程为：

$\left\{\begin{array}{c}\frac{dx}{dt}=x+\mathrm{\&sigma;}(y-x)\\ \frac{dy}{dt}=x(\mathrm{\&rho;}-z)-y\\ \frac{dz}{dt}=xy-\mathrm{\&beta;}z\end{array}\right.---\left(8\right)$

其中(x,y,z)为预报变量，σ、ρ、β为参数。为简化计算，本实施例只计算参数ρ的最优扰动。

公式(8)中，时间积分间隔dt设为0.01，积分步数n＝10，预报变量初始值(x0,y0,z0)＝(10,10,10),参数σ＝10，ρ(1,2…,n)＝28，β＝8/3。

S1、按以上设置编写洛伦兹方程的程序代码，得到t时刻预报变量(x(t),y(t),z(t))；

S2、假设一个随时间变化的参数扰动记作ρ(1,2…,n)。

S3、参数ρ(1,2，…,n)添加随时间变化的扰动ρ’(1,2,…,n),输入程序,得到(x(t)’,y(t)’,z(t)’)；

S4、设置目标函数J为：

J(ρ’)＝-((x(t)’-x(t))²+(y(t)’-y(t))²+(z(t)’-z(t))²)/2(9)

S5、令hρ’(i)＝ρ’(i)+0.04,hp’的其它9项全等于ρ’的对应项，则目标函数J对ρ’(i),i＝1,2,…,n，的偏导数为：

J|ρ’(i)＝(J(hρ’)-J(ρ’))/0.04(10)计算出J对ρ’各分量的偏导数即计算出J对ρ’各分量的梯度G；

S6、扰动约束条件设为：

||ρ’||₂<＝8.8544

(11)

S7、在SPG优化算法中输入多组不同的参数扰动初猜值，扰动约束条件，目标函数，目标函数对扰动的梯度函数，结果输出条件，输出结果。

最优参数扰动输出保存于“fort.1202”文件中，具体为：

3.134808502

3.000556463

2.919795719

2.870522049

2.839016447

2.802882675

2.751427995

2.664029646

2.529738207

2.332206653

具体程序代码参见附带文件。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以作出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (3)

1.一种非恒定的条件非线性最优参数扰动计算方法，其特征在于，包括如下步骤：

S1、假设一个预报系统表达为：

X(t)＝M(p)(X₀)(1)

式中，X为预报变量，X＝(x₁,x₂,…,x_n),n表示变量个数，X₀为初始状态，X(t)为预报结果，p＝(p₁,p₂.....p_m)为系统中的参数(共m个)，M(p)为参数为p时的非线性传播算子；p_i(i＝1,2,....,m)也可以记作p_i(j,k)，j表示不同的空间位置，k表示不同的时间点；p_i(j,k)可以是一个常数，也可以随时间、空间变化；

S2、假设一个随时间、空间变化的参数扰动记作p’，p’＝(p₁’,p₂’,.....,p_m’)，p_i’(i＝1,2,....,m)是随时间和空间变化的量，那么p_i’可以记作p_i’(j,k)，j表示空间位置，k表示时间点；

S3、添加参数扰动之后(1)变成：

X’(t)＝M(p+p’)(X₀)(2)

S4、利用(1)和(2)构造一个能够衡量参数扰动对预报误差影响大小的目标函数J，那么J就是p’的函数，即:

J＝J(p’)(3)

S5、根据梯度定义，估算目标函数J对p_i’(j,k)各分量的梯度函数G；

S6、设置参数扰动约束条件；

S7、在优化算法中输入多组不同的参数扰动初猜值，扰动约束条件，目标函数，目标函数对扰动的梯度函数，结果输出条件，即可计算出能使目标函数最大的参数扰动p’^*,p’^*即是最优参数扰动。

2.根据权利要求1所述的一种非恒定的条件非线性最优参数扰动计算方法，其特征在于，通过如下步骤计算p’^*：

S11、同S1；

S12、同S2；

S13、同S3；

S14、同S4,具体目标函数J为：

J＝-(||X(t)’-X(t)||₂)²/2，(4)

其中||||₂表示L2范数；

S15、根据梯度定义，估算目标函数J对p_i’(j,k)各分量的梯度函数G；具体为，令hp_i’(j,k)＝p_i’(j,k)+h,h取一个适当的小值，令hp_i’的其它项全等于p_i’的对应项，目标函数对p_i’(j,k)的偏导数即为：

J|p_i’(j,k)＝(J(hp_i’)-J(p_i’))/h(5)

计算出J对p_i’(j,k)的各分量的偏导数，即计算出J关于p_i’(j,k)各分量的梯度G；

S16、同S6

S17、同S7。

3.根据权利要求1所述的一种非恒定的条件非线性最优参数扰动计算方法，其特征在于，通过如下步骤计算p’^*：

S21、同权利要求1,具体将p设为0，公式(1)具体变为：

X(t)＝M(0)(X₀)(6)

S22、同S2；

S23、同S3，公式(2)具体变为：

X’(t)＝M(p’)(X₀)(7)

S24、同S14；

S25、同S15；

S26、同S6,具体约束条件设置为：||p_i’||₂<＝r_i,或p_i’(j,k)∈[a_i,b_i](i＝1,2,…m)，其中r_i、a_i、b_i取适当的值；

S27、同S7，具体使用的优化算法为SpectralProjectedGrad(SPG)算法。