基于状态空间模型的浮式结构时域响应分析方法

摘要

本发明公开一种基于状态空间模型的浮式结构时域响应分析方法，包括步骤：(1)用状态空间模型转换Cummins方程；(2)基于极值和留数构造状态空间模型；(3)极值和留数的求解及浮式结构动力响应分析；本发明提出了一种求解浮式结构时域动力响应的新思路——用极值和留数来构造延迟函数的状态空间模型，然后用状态空间模型替换Cummins方程中的卷积项，再用四阶Runge‑Kutta法求解浮体结构的动力响应。在工程领域，能够为浮式海洋结构物的设计及检测工作提供新的技术手段，在浮式结构动力响应分析中具有较高的工程应用前景及推广价值。

说明书

技术领域

本发明涉及一种基于状态空间模型的浮式结构时域响应分析方法。

背景技术

随着海洋油气资源向深海的开发和生产活动，传统的固定导管架平台已不能满足于钻井和生产的需要，这促使许多新型浮式结构和半潜式平台发展起来，动力响应分析是浮式海洋结构在设计阶段的重要内容，对浮式海洋平台的安全运行至关重要。

浮式结构时域运动方程即著名的Cummins方程由Cummins于1962年提出，该方程假设浮式结构的运动是一系列脉冲响应的线性组合，然后将速度势分解为瞬时效应和记忆效应并分别求解。记忆效应的引入提高了动力响应分析的准确性，但Cummins方程中的卷积项不利于动力控制系统的分析和设计，考虑到Cummins方程的时域模拟效率，1964年，Ogilvie将Cummins时域运动方程转到频域来估算海洋结构的动力响应。

为了提高浮体结构时域响应估计的分析效率，在Cummins方程中，有多种方法可以近似替换卷积项，其中一种就是通过状态空间模型替换卷积项，此方法因其高效率和高精度而被许多研究者广泛讨论。但是，传统的方法是直接求解微分方程，比如传统Newmark-β法，存在求解高阶多项式根的病态问题，而且，在计算状态模型矩阵的过程中，需要设置初值，导致计算效率降低、精度差。

发明内容

本发明提出一种基于状态空间模型的浮式结构时域响应分析方法，该方法用极值和留数表示延迟函数，通过构造延迟函数的状态空间模型来替换Cummins方程中的卷积项，最后对浮式结构的动力响应进行分析，为浮式海洋结构物的设计检测工作提供了新的、有效的技术手段，可有效提高了分析计算精度及效率。

本发明是采用以下的技术方案实现的：一种基于状态空间模型的浮式结构时域响应分析方法，包括以下步骤：

步骤A、为提高海洋浮式结构动力响应的计算效率，用状态空间模型转换Cummins方程；

步骤B、基于极值和留数构造状态空间模型：用极值和留数描述延迟函数，然后通过极值留数与传递函数之间的关系，构造延迟函数的状态空间模型矩阵；

步骤C、极值和留数的求解及浮式结构动力响应分析。

进一步的，所述步骤A具体包括以下步骤：

(1)建立时域上的Cummins方程：

式中，K为延迟函数，A(∞)是无限频率的附加质量，M是质量矩阵，C是静水恢复力系数矩阵，x(t)是六自由度浮动结构的响应向量，t为时间，f^exc(t)是激励力；

(2)用状态空间模型代替卷积项：

式中，是系统输出，是相应的输入，i＝1,2,…,6，k＝1,2,…,6，M_ik表示6×6的质量矩阵，x_k(t)表示第k个自由度的位移。

(3)用状态空间模型表示系统的输入输出关系：

式中，[A_ik>ik>ik>ik]是与K_ik(t)相关的状态空间矩阵，K_ik(t)表示6×6的延迟函数矩阵，每个延迟函数K_ik(t)对应一个状态向量Z_ik(t)，该步骤中用状态空间模型替换Cummins方程中的卷积项，并用延迟函数的状态空间模型矩阵表示系统的输入输出关系，大大提高了浮式结构时域响应计算的效率。

进一步的，所述步骤B具体包括以下步骤：

(1)建立浮式结构动力系统的传递函数：

Y(s)是的拉普拉斯变换，U(s)是的拉普拉斯变换，q_n-1,…,q₀及p_m,p_m-1,…,p₀是系数；

(2)令λ_l(l＝1,2,…,n)作为分母sⁿ+q_n-1s^n-1+…+q₁s+q₀＝0的根，则得到：(s-λ₁)(s-λ₂)…(s-λ_n)＝0相应的留数：

则将传递函数写为：

对和U(s)＝sZ_l(s)-λ_lZ_l(s)进行拉普拉斯逆变换，得到时域表达式：

和

式中，是状态向量的元素，z_l(t)是状态变量向量元素，也是Z_l(t)的逆变换，λ_l是极值，l＝1,2,3……；

(3)将y(t)和u(t)转换为状态空间模型的矩阵形式：

则延迟函数的状态空间模型表示为：

式中，

C＝[γ₁γ₂…γ_N],D＝0。

进一步的，所述步骤C包括：

(1)采用等间距采样得到离散延迟函数，采样间隔Δt，则离散延迟函数表示为：

K_ik,l＝K_ik(lΔt)(l＝0，1，2，...，N-1)。

采用复指数分解，用极值和留数描述延迟函数：

(2)构造Hankel矩阵：

式中，ξ和η表示H(l)的行和列的数目；

(3)令l＝0及l＝1，对Hankel矩阵H(0)和H(1)应用奇异值分解技术，得到

的特征值是z_ik,n(n＝1,2,...,p)，由λ_ik,n＝ln(z_ik,n)/Δt和求出延迟函数的极值和留数；

(4)构建A_ik，B_ik和C_ik，并将其组装到一个多自由度状态空间模型里，将多自由度状态空间模型引入Cummins方程，如下面公式所示：

最后采用四阶Runge-Kutta方法求解微分方程，即实现对浮式结构时域响应的分析计算，本步骤中在计算状态空间模型时，采用一阶微分方程。以有效避免分析求解时可能出现的方程病态问题，使得分析计算结果更加精确和稳定。

与现有技术相比，本发明的优点和积极效果在于：

1)本发明创造性的提出通过建立延迟函数的状态空间模型矩阵来替换Cummins方程中的卷积项，是一种新的浮式结构延迟函数的表示方法，该新的时域分析方法为分析计算浮式海洋结构物的动力响应提供了新的技术手段；

2)在估算状态空间模型的时候，用一阶微分方程计算极值和留数，有效避免传统求解方式可能出现的病态问题，使得计算分析结果更加精确和稳定；

3)同时，此方法在分析计算状态空间模型矩阵的过程中不需要任何初值，使得分析计算效率得到了很大的提升，在工程领域，为浮式海洋结构物的设计检测工作提供了新的、有效的技术手段，具有较高的工程应用前景。

附图说明

图1为本发明实施例所述方法原理示意图；

图2为本发明实施例中延迟函数极值、留数分解重构与原始延迟函数的对比示意图；

图3为时域上本发明实施例方法得到的浮式结构运动响应与传统Newmark-β法结果对比示意图；

图4为频域上本发明实施例方法得到的浮式结构运动响应与传统Newmark-β法结果对比示意图。

具体实施方式

为了能够更加清楚地理解本发明的上述目的、特征和优点，下面结合附图及实施例对本发明做进一步说明。需要说明的是，在不冲突的情况下，本申请的实施例及实施例中的特征可以相互组合。

本实施例公开一种基于状态空间模型的浮式结构响应频域分析方法，如图1所示，利用极值和相应的留数来获得延迟函数的状态空间模型，用计算出的极值和留数表示延迟函数，构造延迟函数的状态空间模型来替换Cummins方程中的卷积项，最后用四阶Runge-Kutta法求解浮体结构的动力响应。本发明方法将高阶微分方程转换成一阶微分方程，而不是直接求解微分方程，大大减少了传统方法中需要确定高阶多项式根的病态问题；同时，该方法在计算状态模型矩阵的过程中不需要任何初值，大大提升了计算效率，具体分析如下：

步骤S1、为提高海洋浮式结构动力响应的计算效率，用状态空间模型转换Cummins方程：

(1)对于时域上的Cummins方程：

式中，K为延迟函数，A(∞)是无限频率的附加质量，M是质量矩阵，C是静水恢复力系数矩阵，x(t)是六自由度浮动结构的响应向量，f^exc(t)是激励力；

(2)用状态空间模型代替卷积项：

式中，是系统输出，是相应的输入，i＝1,2,…,6，k＝1,2,…,6，M_ik表示6×6的质量矩阵，x_k(t)表示第k个自由度的位移；

(3)用状态空间模型表示系统的输入输出关系：

式中，[A_ik>ik>ik>ik]是与K_ik(t)相关的状态空间矩阵，K_ik(t)表示6×6的延迟函数矩阵，每个延迟函数K_ik(t)对应一个状态向量Z_ik(t)，该步骤中用状态空间模型替换Cummins方程中的卷积项，并用延迟函数的状态空间模型矩阵表示系统的输入输出关系，大大提高了浮式结构时域响应计算的效率。

步骤S2、基于极值和留数构造状态空间模型：用极值和留数描述延迟函数，然后通过极值留数与传递函数之间的关系，构造延迟函数的状态空间模型矩阵；

(1)浮式结构动力系统的传递函数：

Y(s)是的拉普拉斯变换，U(s)是的拉普拉斯变换，q_n-1,…,q₀及p_m,p_m-1,…,p₀是系数；

(2)令λ_l(l＝1,2,…,n)作为分母sⁿ+q_n-1s^n-1+…+q₁s+q₀＝0的根，则得到：(s-λ₁)(s-λ₂)…(s-λ_n)＝0相应的留数：则传递函数可写为：

通过对和U(s)＝sZ_l(s)-λ_lZ_l(s)进行拉普拉斯逆变换，得到时域的表达式：

和

式中，是状态向量的元素，z_l(t)是状态变量向量元素，也是Z_l(t)的逆变换。λ_l是极值；

(3)将y(t)和u(t)转换为状态空间模型的矩阵形式：

则延迟函数的状态空间模型可以表示为：

式中，

C＝[γ₁>2 …>N],D＝0。

步骤S3、极值和留数的求解及浮式结构的动力响应分析：

(1)采用等间距采样得到离散延迟函数，采样间隔Δt，则离散延迟函数表示为：

K_ik,l＝K_ik(lΔt)(l＝0，1，2，...，N-1)。

应用复指数分解的方法，用极值和留数描述延迟函数：

(2)构造Hankel矩阵：

式中，ξ和η表示H(l)的行和列的数目；

(3)令l＝0及l＝1，对Hankel矩阵H(0)和H(1)应用奇异值分解技术，得到

的特征值是z_ik,n(n＝1,2,...,p)，由λ_ik,n＝ln(z_ik,n)/Δt和求出延迟函数的极值和留数；

(4)构建A_ik，B_ik和C_ik，并将其组装到一个多自由度状态空间模型里，将多自由度状态空间模型引入Cummins方程，如下面公式所示：

最后采用四阶Runge-Kutta方法求解微分方程，即实现对浮式结构时域响应的分析；本步骤中在计算状态空间模型时，采用一阶微分方程。以有效避免分析求解时可能出现的方程病态问题，使得分析计算结果更加精确和稳定。

下面为进一步说明本方法的有效性，本实施例以SESAM软件中自带的半潜式平台模型为例进行介绍，从SESAM软件中提取出该模型相关的水动力参数以及模型结构数据，利用所提取的附加阻尼矩阵，通过离散的数值积分方法得到延迟函数K(t)。图2给出的信号对比图说明本发明方法所求得的极值与留数满足了原始信号的所有成分，重构之后的延迟函数能够代表原延迟函数，证明了本方法的正确性。

根据JONSWAP波谱，利用SESAM的DeepC模块计算了半潜式平台的波浪力，在此波浪力作用下，求解结构动力响应：一种是用基于状态空间模型的四阶Runge-Kutta法，另一种是用传统的Newmark-β法。在求解微分方程之前，将单自由度的估算状态空间模型组装到一个多自由度模型中，由SESAM推导出质量矩阵、附加质量矩阵和静水恢复力系数矩阵。本方法中设微分方程的初始条件为零，两种方法求解结构动力响应的结果如图3所示，相应的振幅与频率如图4所示。从图3和图4中可以看出响应吻合良好，说明了本发明方法的正确性和准确性，为计算分析浮式海洋结构物的动力响应提供了新的技术手段。

以上所述，仅是本发明的较佳实施例而已，并非是对本发明作其它形式的限制，任何熟悉本专业的技术人员可能利用上述揭示的技术内容加以变更或改型为等同变化的等效实施例应用于其它领域，但是凡是未脱离本发明技术方案内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与改型，仍属于本发明技术方案的保护范围。