CN102252897B - 一种深水顶张式立管弯曲振动分析方法 - Google Patents

Abstract

本发明涉及海洋深水立管的研究方法，具体涉及一种深水顶张式立管弯曲振动分析方法。该方法考虑了大位移引起的截面转动对立管弯曲刚度的影响，将截面转动引起的张力作用方向变化转化为立管弯曲刚度的一部分，比现有方法的几何刚度更完整，且更符合深水立管的实际受力和变形状态；同时，该方法还考虑了剪切变形对立管弯曲刚度的影响，将张力与剪切刚度的比值引入弯曲刚度，修正了现有方法的不合理假定，使深水顶张式立管的弯曲振动分析更加符合实际。

Description

一种深水顶张式立管弯曲振动分析方法

技术领域

本发明涉及海洋深水立管的研究方法，具体涉及一种深水顶张式立管弯曲振动分析方法。

背景技术

深水顶张式立管是深水海洋立管的一种主要类型，立管的下端与海底井口连接，管壁张力直接作用在井口上，而管内的流体源源不断地从海底穿过立管流向浮式平台。弯曲振动是深水顶张式立管在海洋环境荷载作用下的主要振动形式，也是引起立管疲劳损伤的主要因素。因此，弯曲振动响应分析是深水顶张式立管设计的重要内容，也是疲劳分析的基础。目前，深水立管的弯曲振动分析方法是基于Eular-Bernoulli梁的纯弯曲理论，该理论建立在小变形假设的基础上，即梁的挠度远远小于其横截面尺寸。

基于该理论的深水顶张式立管弯曲振动方程可表示为：

$\mathrm{EI}\frac{{\∂}^{4}y}{\∂x^4}-\frac{\∂}{\∂x}\left(T\frac{\∂y}{\∂x}\right)+\stackrel{\‾}{m}\frac{{\∂}^{2}y}{\∂x^2}+c\frac{\∂y}{\∂t}=q(x,t)---\left(1\right)$

式中：y为立管横向弯曲位移；

x为立管的轴向坐标；

t为时间；

EI为立管横截面抗弯刚度；

T为立管张力，是时间和立管轴向坐标的函数，即：T＝T(x，t)；

为立管单位长度的质量；

c为阻尼系数；

q(x，t)为作用在立管上的流体荷载。

式(1)是小变形条件下的深水顶张式立管弯曲振动方程，没有考虑大位移引起的立管截面转动对横向弯曲振动的影响，即张力作用方向改变引起的弹性回复力，而仅仅考虑张力作为几何刚度对弯曲刚度的贡献：

$\left[K_G\right]=T{\∫}_0^L{\left[N^{\′}\right]}^T\left[N^{\′}\right]\mathrm{dx}---\left(2\right)$

式中：[K_G]——立管的几何刚度矩阵；

[N]——立管单元的插值函数矩阵。

现有的深水顶张式立管弯曲振动分析方法的主要缺陷如下：

1.没有考虑深水顶张式立管弯曲振动大位移引起的横截面转动对弯曲刚度的影响。

现有技术是建立在纯弯曲理论和小变形假设基础上的，因此，仅考虑张力的几何刚度部分，而没有考虑由于截面转动引起的张力作用方向变化对立管弯曲刚度的影响。

2.没有考虑剪切变形对深水立管弯曲振动的影响。

现有技术基于纯弯曲理论，因此，没有考虑深水顶张式立管剪切变形。如果仅仅从深水顶张式立管的外观几何形状来看，它属于细长梁，其长细比较大，可以采用纯弯曲模型。但立管是管状结构，其截面的剪切刚度较小，因此，其剪切变形的影响较大，特别是张力与截面剪切刚度的比值较大，从而对弯曲刚度影响较大。

发明内容

本发明的目的在于针对现有技术的缺陷，提供一种深水顶张式立管弯曲振动分析方法，使深水顶张式立管的弯曲振动分析更加符合实际。

本发明的技术方案如下：一种深水顶张式立管弯曲振动分析方法，该方法同时考虑大位移引起的截面转动和剪切变形，提供的深水顶张式立管弯曲振动分析模型的方程如下：

$\stackrel{\‾}{m}\frac{{\∂}^{2}y}{\∂t^2}+c\frac{\∂y}{\∂t}+(1+2.6k{\mathrm{\ϵ}}_T)\mathrm{EI}\frac{{\∂}^{4}y}{\∂x^4}-2.6k{\mathrm{\ϵ}}_{T}T\frac{{\∂}^{2}y}{\∂x^2}-(1+2.6k{\mathrm{\ϵ}}_T)\frac{\∂T}{\∂x}\frac{\∂y}{\∂x}=q(x,t)$

式中：y为立管横向弯曲位移；

x为立管的轴向坐标；

t为时间；

EI为立管横截面抗弯刚度；

T为立管张力，是时间和立管轴向坐标的函数，即：T＝T(x，t)；

为立管单位长度的质量；

c为阻尼系数；

k为剪切不均匀系数；

ε_T为张力引起的立管轴向应变；

q(x，t)为作用在立管上的流体荷载；

利用上述方程，计算立管弯曲振动的加速度、速度、位移、应力和应变随时间的变化。

进一步，如上所述的深水顶张式立管弯曲振动分析方法，该方法计算立管弯曲振动的加速度、速度、位移、应力和应变随时间变化的具体过程如下：

(1)将立管划分为若干个单元；

(2)将单元的位移函数表示为插值函数的形式：y＝[N]{a}，

其中，[N]是单元插值函数矩阵，{a}是单元节点位移向量；

(3)采用伽辽金方法，将深水顶张式立管弯曲振动分析模型的方程转换为矩阵方程如下：

$\left[M\right]\left\{\stackrel{\·\·}{a}\right\}+\left[C\right]\left\{\stackrel{\·}{a}\right\}+\left[K\right]\left\{a\right\}=\left\{F\right\}$

式中：

阻尼矩阵[C]＝α[M]+β[K]，α，β为瑞雷阻尼系数；

荷载向量

为加速度矢量；

为速度矢量；

{a}为位移矢量；

l为单元长度；

n为单元数量；

(4)基于立管的初始张力求出立管的初始轴向应变ε_T，并将其带入步骤(3)中刚度矩阵[K]的表达式，计算出立管的初始刚度矩阵；

(5)由步骤(3)中质量矩阵[M]、阻尼矩阵[C]和荷载向量{F}的表达式分别计算出立管的质量矩阵、阻尼矩阵和荷载向量；

(6)将立管的质量矩阵、阻尼矩阵、刚度矩阵和荷载向量代入步骤(3)中的矩阵方程，采用逐步积分法按下式求出立管第一个时间增量Δt后的加速度增量速度增量

和位移增量{Δa_i}，此时i＝0：

$\left(\frac{6}{\mathrm{\Δ}{t}^2}\right[M]+\frac{3}{\mathrm{\Δt}}[C]+[K\left]\right)\left\{\mathrm{\Δ}{a}_i\right\}=\left(\frac{6}{\mathrm{\Δt}}\right[M]+3[C\left]\right)\left\{{\stackrel{\·}{a}}_i\right\}+\left(3\right[M]+\frac{\mathrm{\Δt}}{2}[C\left]\right)\left\{{\stackrel{\·\·}{a}}_i\right\}+\left\{\mathrm{\Δ}{F}_i\right\};$

{ΔF_i}是与时间增量Δt对应的荷载增量；

(7)由加速度增量、速度增量和位移增量按下式计算第一个时间增量后的加速度、速度和位移，此时i＝0：

a_i+1＝a_i+Δa_i

${\stackrel{\·}{a}}_{i+1}={\stackrel{\·}{a}}_i+{\mathrm{\Δ}\stackrel{\·}{a}}_i;$

${\stackrel{\·\·}{a}}_{i+1}={\stackrel{\·\·}{a}}_i+{\mathrm{\Δ}\stackrel{\·\·}{a}}_i$

(8)根据步骤(7)中得到的位移计算出立管的应力和轴向应变；

(9)将步骤(8)中得到的轴向应变代入步骤(3)中刚度矩阵[K]的表达式，计算出第一个时间增量后的立管刚度矩阵；

(10)重复步骤(5)～(9)，直至时间t达到要求的时长，即可计算出立管弯曲振动的加速度、速度、位移、应力和轴向应变随时间的变化。

本发明的有益效果如下：本发明在深水顶张式立管的弯曲振动分析中，考虑了大位移引起的截面转动对立管弯曲刚度的影响，将截面转动引起的张力作用方向变化转化为立管弯曲刚度的一部分，比现有方法的几何刚度更完整，且更符合深水顶张式立管的实际受力和变形状态；同时，本发明还考虑了剪切变形对立管弯曲刚度的影响，将张力与剪切刚度的比值引入弯曲刚度，修正了现有方法的不合理假定，使深水顶张式立管的弯曲振动分析更加符合实际。

附图说明

图1为本发明的深水顶张式立管的弯曲振动分析方法流程图；

图2为立管小变形的微元体模型示意图；

图3为立管大变形的微元体模型示意图；

图4为纯弯曲梁的微元体模型示意图。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细描述。

本发明采用同时考虑大位移引起的截面转动和剪切变形的深水顶张式立管弯曲振动分析模型，方程如下：

$\stackrel{\‾}{m}\frac{{\∂}^{2}y}{\∂t^2}+c\frac{\∂y}{\∂t}+(1+2.6k{\mathrm{\ϵ}}_T)\mathrm{EI}\frac{{\∂}^{4}y}{\∂x^4}-2.6k{\mathrm{\ϵ}}_{T}T\frac{{\∂}^{2}y}{\∂x^2}-(1+2.6k{\mathrm{\ϵ}}_T)\frac{\∂T}{\∂x}\frac{\∂y}{\∂x}=q(x,t)---\left(3\right)$

式中：y为立管横向弯曲位移；

x为立管的轴向坐标；

t为时间；

EI为立管横截面抗弯刚度；

T为立管张力，是时间和立管轴向坐标的函数，即：T＝T(x，t)；

为立管单位长度的质量；

c为阻尼系数；

k为剪切不均匀系数；

ε_T为张力引起的立管轴向应变；

q(x，t)为作用在立管上的流体荷载。

上式中，2.6kε_T即为张力和剪切变形对立管弯曲刚度的影响。

深水顶张式立管弯曲振动分析方法的具体步骤如下：

(S1)将立管划分为若干个单元；

(S2)将单元的位移函数表示为插值函数的形式：

y＝[N]{a}                    (4)

其中，[N]是单元插值函数矩阵，{a}是单元节点位移向量；

(S3)采用伽辽金方法(此为本领域的公知常识)，

$\underset{0}{\overset{L}{\∫}}\mathrm{\δy}[\stackrel{\‾}{m}\frac{{\∂}^{2}y}{\∂t^2}+c\frac{\∂y}{\∂t}+(1+2.6k{\mathrm{\ϵ}}_T)\mathrm{EI}\frac{{\∂}^{4}y}{\∂x^4}-$

$2.6k{\mathrm{\ϵ}}_{T}T\frac{{\∂}^{2}y}{\∂x^2}-(1+2.6k{\mathrm{\ϵ}}_T)\frac{\∂T}{\∂x}\frac{\∂y}{\∂x}-q(x,t)]=0---\left(5\right)$

上式中，L是立管的总长度，δy是弯曲位移的变分，在公式(5)中的作用是权函数，此为公知常识；

将深水顶张式立管弯曲振动分析模型的方程(公式(3))转换为矩阵方程如下：

$\left[M\right]\left\{\stackrel{\·\·}{a}\right\}+\left\{C\right\}\left\{\stackrel{\·}{a}\right\}+\left[K\right]\left\{a\right\}=\left\{F\right\}---\left(6\right)$

式中：

阻尼矩阵[C]＝α[M]+β[K]，α，β为瑞雷阻尼系数；        (9)

荷载向量

为加速度矢量；

为速度矢量；

{a}为位移矢量；

l为单元长度；

n为单元数量；

(S4)基于立管的初始张力求出立管的初始轴向应变ε_T，并将其带入步骤(S3)中刚度矩阵[K]的表达式(公式(8))，计算出立管的初始刚度矩阵；

(S5)由步骤(S3)中质量矩阵[M]、阻尼矩阵[C]和荷载向量{F}的表达式(公式(7)、(9)、(10))分别计算出立管的质量矩阵、阻尼矩阵和荷载向量；

(S6)将立管的质量矩阵、阻尼矩阵、刚度矩阵和荷载向量代入步骤(S3)中的矩阵方程(公式(6))，采用逐步积分法按下式求出立管第一个时间增量Δt(此时t＝Δt)后的加速度增量

速度增量

和位移增量{Δa_i}，此时i＝0：

$\left(\frac{6}{\mathrm{\Δ}{t}^2}\right[M]+\frac{3}{\mathrm{\Δt}}[C]+[K\left]\right)\left\{\mathrm{\Δ}{a}_i\right\}=\left(\frac{6}{\mathrm{\Δt}}\right[M]+3[C\left]\right)\left\{{\stackrel{\·}{a}}_i\right\}+\left(3\right[M]+\frac{\mathrm{\Δt}}{2}[C\left]\right)\left\{{\stackrel{\·\·}{a}}_i\right\}+\left\{\mathrm{\Δ}{F}_i\right\},$

此处，{ΔF_i}是与时间增量Δt对应的荷载增量；

(S7)由加速度增量、速度增量和位移增量按下式计算第一个时间增量后的加速度、速度和位移，此时i＝0：

a_i+1＝a_i+Δa_i

${\stackrel{\·}{a}}_{i+1}={\stackrel{\·}{a}}_i+\mathrm{\Δ}{\stackrel{\·}{a}}_i;$

${\stackrel{\·\·}{a}}_{i+1}={\stackrel{\·\·}{a}}_i+\mathrm{\Δ}{\stackrel{\·\·}{a}}_i$

(S8)根据步骤(S7)中得到的位移计算出立管的应力和轴向应变；

(S9)将步骤(S8)中得到的轴向应变代入步骤(S3)中刚度矩阵[K]的表达式(公式(8))，计算出第一个时间增量后的立管刚度矩阵；

(S10)重复步骤(S5)～(S9)，直至时间t达到要求的时长，即可计算出立管弯曲振动的加速度、速度、位移、应力和轴向应变随时间的变化(即时程)。

本发明考虑了深水顶张式立管弯曲振动的横向大位移引起的横截面转动对立管弯曲刚度的影响，比现有方法的几何刚度更完整，且更符合深水顶张式立管的实际受力和变形状态。

现有方法建立在小变形假设的基础上，即假定立管变形后的微元段两截面仍保持平行，如图2所示。但深水立管的横向弯曲振动属于大变形(大位移)，其变形后的微元段两横截面有相对转动，如图3所示。因此，变形后的张力不再保持水平方向，其横向分量将提供部分回复力，从而增大了弯曲刚度。本发明将张力对回复力的贡献从仅仅考虑微元体两截面因剪切变形发生相对错动(见图2)的几何刚度发展成同时考虑截面相对错动和相对转动(见图3)耦合的几何刚度。

另外，本发明考虑了深水立管剪切变形对弯曲刚度的影响，修正了现有方法的不合理假定。

现有方法采用图2微元体模型，这与纯弯曲梁(不考虑剪切变形)不同，两横截面之间有上下的错动。而纯弯曲梁的微元体模型如图4所示。

如果按照图4的纯弯曲梁模型，则张力对弯曲刚度没有影响，即几何刚度等于零[K_G]＝0。现有技术为了能够考虑张力梁的几何刚度而采用了图2的模型，这意味着考虑了剪切变形，但却采用了纯弯曲的数学模型(此为公知技术)：

$\frac{\mathrm{dy}}{\mathrm{dx}}=\mathrm{\θ}$

式中：θ为弯曲引起的梁截面转角，即挠度的导数等于转角。这表明，现有技术的数学模型和立学模型不一致。本发明克服了现有技术的这一不足，采用了考虑剪切变形的剪弯梁数学模型：

$\frac{\mathrm{dy}}{\mathrm{dx}}=\mathrm{\θ}+\mathrm{\γ}$

式中：γ为梁截面的剪应变，即挠度的导数等于截面转角与剪应变之和。

显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若对本发明的这些修改和变型属于本发明权利要求及其同等技术的范围之内，则本发明也意图包含这些改动和变型在内。

Claims (1)

1.一种深水顶张式立管弯曲振动分析方法，其特征在于：该方法同时考虑大位移引起的截面转动和剪切变形，提供的深水顶张式立管弯曲振动分析模型的方程如下：

$\stackrel{\‾}{m}\frac{{\∂}^{2}y}{{\∂t}^2}+c\frac{\∂y}{\∂t}+(1+2.6{\mathrm{k\ϵ}}_T)\mathrm{EI}\frac{{\∂}^{4}y}{{\∂x}^4}-2.6{\mathrm{ks}}_{T}T\frac{{\∂}^{2}y}{{\∂x}^2}-(1+2.6{\mathrm{k\ϵ}}_T)\frac{\∂T}{\∂x}\frac{\∂y}{\∂x}=q(x,t)$

式中：y为立管横向弯曲位移；

x为立管的轴向坐标；

t为时间；

EI为立管横截面抗弯刚度；

T为立管张力，是时间和立管轴向坐标的函数，即：T＝T(x，t)；

为立管单位长度的质量；

c为阻尼系数；

k为剪切不均匀系数；

ε_T为张力引起的立管轴向应变；

q(x，t)为作用在立管上的流体荷载；

利用上述方程，计算立管弯曲振动的加速度、速度、位移、应力和应变随时间的变化，具体过程如下：

(1)将立管划分为若干个单元；

(2)将单元的位移函数表示为插值函数的形式：y＝[N]{a}，

其中，[N]是单元插值函数矩阵，{a}是单元节点位移向量；

(3)采用伽辽金方法，将深水顶张式立管弯曲振动分析模型的方程转换为矩阵方程如下：

$\left[M\right]\left\{\stackrel{\·\·}{z}\right\}+\left[C\right]\left\{\stackrel{\·}{a}\right\}+\left[K\right]\left\{a\right\}=\left\{F\right\}$

式中：

质量矩阵 $\left[M\right]=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}\stackrel{\‾}{m}\underset{0}{\overset{l}{\∫}}{\left[N\right]}^T\left[N\right]\mathrm{dx};$

刚度矩阵 $\left[K\right]=\stackrel{n}{\mathrm{\Σ}}\left[(1+2.6{\mathrm{k\ϵ}}_T)\underset{0}{\overset{l}{\∫}}{\left[N^{\′\′}\right]}^T\right[N^{\′\′}]\mathrm{dx}+2.6{\mathrm{k\ϵ}}_{T}T\underset{0}{\overset{l}{\∫}}{\left[N^{\′}\right]}^T[N^{\′}]\mathrm{dx}$ ；

$+(1+2.6{\mathrm{k\ϵ}}_T)\underset{0}{\overset{l}{\∫}}\frac{\∂T}{\∂x}{\left[N^{\′}\right]}^T\left[N\right]\mathrm{dx}]$

阻尼矩阵[C]＝α[M]+β[K]，α，β为瑞雷阻尼系数；

荷载向量 $\left\{F\right\}={\∫}_0^l{\left[N\right]}^{T}q(x,t)\mathrm{dx};$

为加速度矢量；

为速度矢量；

{a}为位移矢量；

l为单元长度；

n为单元数量；

(4)基于立管的初始张力求出立管的初始轴向应变ε_T，并将其带入步骤(3)中刚度矩阵[K]的表达式，计算出立管的初始刚度矩阵；

(5)由步骤(3)中质量矩阵[M]、阻尼矩阵[C]和荷载向量{F}的表达式分别计算出立管的质量矩阵、阻尼矩阵和荷载向量；

(6)将立管的质量矩阵、阻尼矩阵、刚度矩阵和荷载向量代入步骤(3)中的矩阵方程，采用逐步积分法按下式求出立管第一个时间增量Δt后的加速度增量

速度增量

和位移增量{Δa_i}，此时i＝0：

$\left(\frac{6}{{\mathrm{\Δt}}^2}\right[M]+\frac{3}{\mathrm{\Δt}}[C]+[K\left]\right)\left\{{\mathrm{\Δa}}_i\right\}=\left(\frac{6}{\mathrm{\Δt}}\right[M]+3[C\left]\right)\left\{{\stackrel{\·}{a}}_i\right\}+\left(3\right[M]+\frac{\mathrm{\Δt}}{2}[C\left]\right)\left\{{\stackrel{\·\·}{a}}_i\right\}+\left\{{\mathrm{\ΔF}}_i\right\},$

此处，{ΔF_i}是与时间增量Δt对应的荷载增量；

(7)由加速度增量、速度增量和位移增量按下式计算第一个时间增量后的加速度、速度和位移，此时i＝0：

a_i+1＝a_i+Δa_i

${\stackrel{\·}{a}}_{i-1}={\stackrel{\·}{a}}_i+\mathrm{\Δ}{\stackrel{\·}{a}}_i;$

${\stackrel{\·\·}{a}}_{i+1}={\stackrel{\·\·}{a}}_i+\mathrm{\Δ}{\stackrel{\·\·}{a}}_i$

(8)根据步骤(7)中得到的位移计算出立管的应力和轴向应变；

(9)将步骤(8)中得到的轴向应变代入步骤(3)中刚度矩阵[K]的表达式，计算出第一个时间增量后的立管刚度矩阵；

(10)重复步骤(5)～(9)，直至时间t达到要求的时长，即可计算出立管弯曲振动的加速度、速度、位移、应力和轴向应变随时间的变化。

Images