[发明专利]基于Johnson变换的非高斯风压模拟方法

权利要求书

1.一种基于Johnson变换的非高斯风压模拟方法，其特征在于，

步骤1：非高斯的Johnson变换模型

Suppose Y(t)是均值为μ_Y和标准差为σ_Y的非高斯风压过程；标准非高斯过程即标准非高斯序列，X(t)＝[Y(t)-μ_Y]/σ_Y转换为平稳标准高斯序列Z(t)，然后再表示为(Grigoriu1998)

x(t)＝g[z(t)] (1)

在这里F_X是X(t)的累积分布函数；是F_X的反函数；Φ_Z(·)为Z(t)的累积分布函数；

基于中心极限定理并参照Perason的四参数系统，将平移标准高斯序列Z(t)转换为非高斯序列X(t)的四参数转换模型，该模型被称为Johnson转化模型；Johnson转换模型表示为：

其中，J(·)分别为：

(a)无界转换模型，S_U

(b)有界转换模型，S_B

(c)对数正态转换模型，S_L

其中ε和γ是控制Johnson曲线位置的参数；λ和η是控制Johnson曲线尺度的参数，其值总是大于零；

为了下面分析的简化，定义需要注意的是线性变换不改变偏度和峰度；

式(2a)对应的逆转化模型(I-Johnson转换模型)为：

其中，J^-1(·)分别为：

(a)无界转换模型，S_U

(b)有界转换模型，S_B

(c)对数正态转换模型，S_L

同时，基于(3a)式，求出非高斯变量q的概率分布函数为：

(a)无界转化模型，S_U：

(b)有界转化模型，S_B：

(c)对数正态转换模型，S_L:

式中，是标准高斯过程Z的概率密度函数；

Johnson转换的适用范围中的S_L转换适用范围是偏度-峰度的一条曲线，该曲线的闭合表达式为(Hill et al.1976)：

α₄＝w⁴+2w³+3w²-3 (5b)

其中，α₃表示偏度；α₄表示峰度；w＝exp(η^-2)，

对于上述三种转换的选取，取决于指定的偏度和峰度；对于指定的偏度α₃，公式(5a)在复数域有三个根；由于w>0，故w取实数正根值；分别按下式计算：

式中，Q₁和Q₂表示为：

随后将求得的w带入公式(5b)得到α'₄；若α₄＜α'₄，选取S_B转换；若α₄>α'₄，选取S_U转换；当α₄＝α'₄时，即该点落在S_L曲线上，选取S_L转换；

步骤2：矩估计法

Hill et al(1976)基于矩估计法对Johnson转化模型的参数进行了估计；非高斯变量q的前四阶中心矩r_n(n＝1,2,3,4)表示为：

(1),无界转换模型，S_U:

r₁(q)＝-w^0.5sinh(Ω) (8a)

r₂(q)＝0.5(w-1)[wcosh(2Ω)+1] (8b)

r₃(q)＝-0.25w^0.5(w-1)²[w(w+2)sinh(3Ω)+3sinh(Ω)] (8c)

(2),有界转换模型，S_B

(3),有界对数正态转换模型，S_L:

r_n(q)＝exp(0.5n²η^-2-nγη^-1) (10)

其中，Ω＝γ/η；非高斯变量q的偏度α₃和峰度α₄表示为：偏度和峰度仅与γ和η有关；基于样本数据的偏度和峰度，即先求出γ和η；再根据关系和ε＝r₁(x)-λr₁(q)求出λ和ε；

步骤3：高斯过程的模拟

谱表示法是一种利用谱分解和三角函数技术叠加来模拟随机过程样本的传统方法；以一维单变量零均值平稳高斯随机过程为例，该方法表示为：

A_l＝(2S_Z(ω_l)△ω)^1/2,l＝0,1,2,....,N-1 (11b)

ω_l＝l△ω,l＝0,1,2,....,N-1 (11c)

△ω＝ω_u/N (12)

其中，为模拟的平稳高斯过程样本；N为频率分段数，一般取N＝2^μ，μ为正整数，以便于进行快速傅里叶变换；S_Z为高斯过程的目标功率谱；ω_l为频率；△ω为频率间隔；ω_u是上截频率，其值的确定见Shinozuka和Deodatis(1991)；φ_l为随机相位，取0～2π之间的随机数；

步骤4：功率谱密度匹配关系

在模拟得到了高斯平稳随机样本后，将这些样本通过Johnson转换模型即得到相应的平稳非高斯随机样本；但模拟前往往只知道非高斯过程X(t)的功率谱或相关函数；高斯与非高斯过程的功率谱密度函数存在偏差；因此，需要先基于X(t)的功率谱或相关函数得到相应的高斯过程的功率谱或相关函数；

Grigoriu(1998)给出了非高斯过程X(t)的自相关系数函数ρ_NG(τ)为：

其中，z₁＝z(t)，z₂＝z(t+τ)，表示二元标准高斯向量的联合概率密度函数，其表达式为：

当采用Johnson转换模型时，对于S_L和S_U转换模型，由式(13)获得求解ρ_NG(τ)的解析式：

(1)，对数高斯转换，即S_L转换

(2),无界转换，即S_U转换

对于S_B转换模型，式(13)很难获得ρ_NG(τ)的解析式，但能通过近似得到其解为

式中，该公式推导如下：

基于Hermite多项式，公式(14)展开为

式中：He_s(·)是s阶Hermite多项式，由下式计算

当|ρ_G|＜1时，式(16)会收敛，否则会发散；

将式(16)代入式(13),ρ_NG(τ)表示为

式中，由于Hermite多项式的正交性，I_1,0＝I_2,0＝0；

忽略s>M的项，式(18)近似为

当引入Gauss-Hermite求积公式，I_k,s近似表示为

式中，d是样本数目且建议d＝9～11；z_gh,a是Hermite多项式的解；w_gh,a是相应的权重，表示为

已知ρ_NG(τ)，通过牛顿迭代法求解上述非线性方程(15)从而求出ρ_G(τ)；其中ρ_NG(τ)和ρ_G(τ)应满足以下关系：

|ρ_NG(τ)|≤|ρ_G(τ)| (22)

上述公式的解不一定存在，当ρ_G(τ)＝-1和1时，相应的ρ_NG(τ)值分别为和当ρ_NG(τ)的取值区间在内时，上述公式有解；假如有解，还需要验证ρ_G是非负定矩阵(Gioffre et al.2000)；其中，ρ_G为

式中，τ_ij＝t_i-t_j,i,j＝1,2,...k，k为正整数；

假如上述矩阵不满足非负定，这就意味着指定的非高斯过程X(t)的功率谱密度函数S_X(ω)和概率分布F_X不匹配；此时，通过Shields et al.(2011)提出的方法来解决；

步骤5：整体模拟流程

基于Johnson转换模型的非高斯平稳随机过程整体模拟流程如下：

1，根据非高斯平稳随机过程Y(t)的目标功率谱密度S_Y(ω)和前四阶统计矩：μY,σ_Y,α₃,α₄，

获得标准化非高斯平稳随机过程X(t)的目标功率谱密度S_NG(ω)和前四阶统计矩：μ_X＝0,σ_X＝1,α₃,α₄；

2，根据非高斯过程X(t)的前四阶矩，确定Johnson转换类型及相应Johnson模型的参数；

3，基于功率谱密度和概率分布是否匹配，分为下面两种方式进行模拟：

3.1当功率谱密度和概率分布是匹配的时候,整个模拟过程如下：先基于逆傅里叶变换将S_NG(ω)转换为ρ_NG(τ)，将ρ_NG(τ)代入式(13)或(15a)、(15b)、(15c)进行迭代计算得到ρ_G(τ)，将ρ_G(τ)进行傅里叶变换得到S_Z(ω)，基于谱表示法模拟得到高斯过程的样本再将其代入(3a)式得到非高斯的模拟样本

3.2当功率谱密度和概率分布不匹配的时候，则模拟按照下面进行

3.2.1假定S_NG(ω)为S_G(ω)的初值

3.2.2，将进行IFFT得到通过(14)式计算求得将进行傅里叶变换(FFT)得到使用作为下一个高斯分布的估计值:

式中，表示作为第J_k+1次迭代高斯过程的谱密度；κ是优化迭代收敛速率的指数，当1.3≤κ≤1.5，迭代速度和精度都能够得到较好的保证(e.g.,Shields et al.2011)；

3.2.3，用更新的高斯过程的谱密度替换重复3.2.2，直到迭代产生的相对误差稳定在一常数值(或进一步迭代不能够得到更小的相对误差)，迭代结束；相对误差按下式计算：

4，根据得到的高斯过程的功率谱密度，采用谱表示法模拟并得到标准高斯过程样本再将通过Johnson转换模型得到归一化的非高斯过程模拟样本再由关系得到非高斯过程Y(t)的模拟样本。